diff --git a/.github/workflows/python-test.yml b/.github/workflows/python-test.yml index d12735354..e892f30bd 100644 --- a/.github/workflows/python-test.yml +++ b/.github/workflows/python-test.yml @@ -23,6 +23,9 @@ jobs: pip install z3-solver pip install exact pip install pysdd + pip install pychoco + sudo snap install minizinc --classic + pip install minizinc - name: Test with pytest run: | python -m pytest tests/ diff --git a/.gitignore b/.gitignore index 286a1e803..400278d4a 100644 --- a/.gitignore +++ b/.gitignore @@ -140,4 +140,9 @@ cython_debug/ bugs/ # Intellij -.idea/ \ No newline at end of file +.idea/ + +# VeriPB proof log files (from Glasgow Constraint solver) +*.veripb +*.opb +*.pbp diff --git a/cpmpy/__init__.py b/cpmpy/__init__.py index f0b40f188..792919a34 100644 --- a/cpmpy/__init__.py +++ b/cpmpy/__init__.py @@ -12,9 +12,9 @@ - `solvers`: CPMpy classes that translate a model into approriate calls of a solver's API - `transformations`: common methods for transforming expressions into other expressions, used by `solvers` modules to simplify/rewrite expressions """ -# Tias Guns, 2019-2023 +# Tias Guns, 2019-2024 -__version__ = "0.9.18" +__version__ = "0.9.21" from .expressions import * diff --git a/cpmpy/exceptions.py b/cpmpy/exceptions.py index ced4d1be0..7ed0e2143 100644 --- a/cpmpy/exceptions.py +++ b/cpmpy/exceptions.py @@ -17,6 +17,12 @@ class MinizincNameException(CPMpyException): class MinizincBoundsException(CPMpyException): pass +class ChocoTypeException(CPMpyException): + pass + +class ChocoBoundsException(CPMpyException): + pass + class NotSupportedError(CPMpyException): pass diff --git a/cpmpy/expressions/__init__.py b/cpmpy/expressions/__init__.py index 9dc3b1f60..e86de9032 100644 --- a/cpmpy/expressions/__init__.py +++ b/cpmpy/expressions/__init__.py @@ -21,8 +21,11 @@ # others need to be imported by the developer explicitely from .variables import boolvar, intvar, cpm_array from .variables import BoolVar, IntVar, cparray # Old, to be deprecated -from .globalconstraints import AllDifferent, AllDifferentExcept0, AllEqual, Circuit, Inverse, Table, Xor, Cumulative, IfThenElse, GlobalCardinalityCount, DirectConstraint, InDomain +from .globalconstraints import AllDifferent, AllDifferentExcept0, AllDifferentExceptN, AllEqual, AllEqualExceptN, Circuit, Inverse, Table, Xor, Cumulative, \ + IfThenElse, GlobalCardinalityCount, DirectConstraint, InDomain, Increasing, Decreasing, IncreasingStrict, DecreasingStrict, \ + LexLess, LexLessEq, LexChainLess, LexChainLessEq, Precedence, NoOverlap, \ + NegativeTable from .globalconstraints import alldifferent, allequal, circuit # Old, to be deprecated -from .globalfunctions import Maximum, Minimum, Abs, Element, Count +from .globalfunctions import Maximum, Minimum, Abs, Element, Count, NValue, NValueExcept, Among from .core import BoolVal from .python_builtins import all, any, max, min, sum diff --git a/cpmpy/expressions/core.py b/cpmpy/expressions/core.py index 1b210e680..47d03be41 100644 --- a/cpmpy/expressions/core.py +++ b/cpmpy/expressions/core.py @@ -11,7 +11,7 @@ Here is a list of standard python operators and what object (with what expr.name) it creates: - Comparisons: + **Comparisons**: - x == y Comparison("==", x, y) - x != y Comparison("!=", x, y) @@ -20,9 +20,9 @@ - x > y Comparison(">", x, y) - x >= y Comparison(">=", x, y) - Mathematical operators: + **Mathematical operators**: - - -x Operator("-", [x]) + - −x Operator("-", [x]) - x + y Operator("sum", [x,y]) - sum([x,y,z]) Operator("sum", [x,y,z]) - sum([c0*x, c1*y, c2*z]) Operator("wsum", [[c0,c1,c2],[x,y,z]]) @@ -32,7 +32,7 @@ - x % y Operator("mod", [x,y]) - x ** y Operator("pow", [x,y]) - Logical operators: + **Logical operators**: - x & y Operator("and", [x,y]) - x | y Operator("or", [x,y]) @@ -47,14 +47,16 @@ Apart from operator overloading, expressions implement two important functions: - - `is_bool()` which returns whether the __return type__ of the expression is Boolean. - If it does, the expression can be used as top-level constraint - or in logical operators. + - `is_bool()` + which returns whether the return type of the expression is Boolean. + If it does, the expression can be used as top-level constraint + or in logical operators. - - `value()` computes the value of this expression, by calling .value() on its - subexpressions and doing the appropriate computation - this is used to conveniently print variable values, objective values - and any other expression value (e.g. during debugging). + - `value()` + computes the value of this expression, by calling .value() on its + subexpressions and doing the appropriate computation + this is used to conveniently print variable values, objective values + and any other expression value (e.g. during debugging). =============== List of classes @@ -72,7 +74,7 @@ import numpy as np -from .utils import is_num, is_any_list, flatlist, argval, get_bounds, is_boolexpr, is_true_cst, is_false_cst +from .utils import is_num, is_any_list, flatlist, argval, get_bounds, is_boolexpr, is_true_cst, is_false_cst, argvals from ..exceptions import IncompleteFunctionError, TypeError @@ -144,6 +146,7 @@ def is_bool(self): def value(self): return None # default + def get_bounds(self): if self.is_bool(): return 0, 1 #default for boolean expressions @@ -400,7 +403,8 @@ def __repr__(self): # return the value of the expression # optional, default: None def value(self): - arg_vals = [argval(a) for a in self.args] + arg_vals = argvals(self.args) + if any(a is None for a in arg_vals): return None if self.name == "==": return arg_vals[0] == arg_vals[1] elif self.name == "!=": return arg_vals[0] != arg_vals[1] @@ -526,11 +530,12 @@ def wrap_bracket(arg): return "{}({})".format(self.name, self.args) def value(self): + if self.name == "wsum": # wsum: arg0 is list of constants, no .value() use as is - arg_vals = [self.args[0], [argval(arg) for arg in self.args[1]]] + arg_vals = [self.args[0], argvals(self.args[1])] else: - arg_vals = [argval(arg) for arg in self.args] + arg_vals = argvals(self.args) if any(a is None for a in arg_vals): return None @@ -546,7 +551,8 @@ def value(self): try: return arg_vals[0] // arg_vals[1] except ZeroDivisionError: - raise IncompleteFunctionError(f"Division by zero during value computation for expression {self}") + raise IncompleteFunctionError(f"Division by zero during value computation for expression {self}" + + "\n Use argval(expr) to get the value of expr with relational semantics.") # boolean elif self.name == "and": return all(arg_vals) diff --git a/cpmpy/expressions/globalconstraints.py b/cpmpy/expressions/globalconstraints.py index 5ed2be859..756d79228 100644 --- a/cpmpy/expressions/globalconstraints.py +++ b/cpmpy/expressions/globalconstraints.py @@ -98,13 +98,23 @@ def my_circuit_decomp(self): AllDifferent AllDifferentExcept0 + AllDifferentExceptN AllEqual + AllEqualExceptN Circuit Inverse Table + NegativeTable Xor Cumulative + IfThenElse GlobalCardinalityCount + DirectConstraint + InDomain + Increasing + Decreasing + IncreasingStrict + DecreasingStrict """ import copy @@ -113,7 +123,7 @@ def my_circuit_decomp(self): from ..exceptions import CPMpyException, IncompleteFunctionError, TypeError from .core import Expression, Operator, Comparison from .variables import boolvar, intvar, cpm_array, _NumVarImpl, _IntVarImpl -from .utils import flatlist, all_pairs, argval, is_num, eval_comparison, is_any_list, is_boolexpr, get_bounds +from .utils import flatlist, all_pairs, argval, is_num, eval_comparison, is_any_list, is_boolexpr, get_bounds, argvals from .globalfunctions import * # XXX make this file backwards compatible @@ -171,25 +181,37 @@ def decompose(self): return [var1 != var2 for var1, var2 in all_pairs(self.args)], [] def value(self): - return len(set(a.value() for a in self.args)) == len(self.args) + return len(set(argvals(self.args))) == len(self.args) - -class AllDifferentExcept0(GlobalConstraint): +class AllDifferentExceptN(GlobalConstraint): """ - All nonzero arguments have a distinct value + All arguments except those equal to a value in n have a distinct value. """ - def __init__(self, *args): - super().__init__("alldifferent_except0", flatlist(args)) + def __init__(self, arr, n): + flatarr = flatlist(arr) + if not is_any_list(n): + n = [n] + super().__init__("alldifferent_except_n", [flatarr, n]) def decompose(self): - # equivalent to (var1 == 0) | (var2 == 0) | (var1 != var2) - return [(var1 == var2).implies(var1 == 0) for var1, var2 in all_pairs(self.args)], [] + from .python_builtins import any as cpm_any + # equivalent to (var1 == n) | (var2 == n) | (var1 != var2) + return [(var1 == var2).implies(cpm_any(var1 == a for a in self.args[1])) for var1, var2 in all_pairs(self.args[0])], [] def value(self): - vals = [a.value() for a in self.args if a.value() != 0] + vals = [argval(a) for a in self.args[0] if argval(a) not in argvals(self.args[1])] return len(set(vals)) == len(vals) +class AllDifferentExcept0(AllDifferentExceptN): + """ + All nonzero arguments have a distinct value + """ + def __init__(self, *arr): + flatarr = flatlist(arr) + super().__init__(arr, 0) + + def allequal(args): warnings.warn("Deprecated, use AllEqual(v1,v2,...,vn) instead, will be removed in stable version", DeprecationWarning) return AllEqual(*args) # unfold list as individual arguments @@ -208,7 +230,27 @@ def decompose(self): return [var1 == var2 for var1, var2 in zip(self.args[:-1], self.args[1:])], [] def value(self): - return len(set(a.value() for a in self.args)) == 1 + return len(set(argvals(self.args))) == 1 + + +class AllEqualExceptN(GlobalConstraint): + """ + All arguments except those equal to a value in n have the same value. + """ + + def __init__(self, arr, n): + flatarr = flatlist(arr) + if not is_any_list(n): + n = [n] + super().__init__("allequal_except_n", [flatarr, n]) + + def decompose(self): + from .python_builtins import any as cpm_any + return [(cpm_any(var1 == a for a in self.args[1]) | (var1 == var2) | cpm_any(var2 == a for a in self.args[1])) for var1, var2 in all_pairs(self.args[0])], [] + + def value(self): + vals = [argval(a) for a in self.args[0] if argval(a) not in argvals(self.args[1])] + return len(set(vals)) == 1 or len(set(vals)) == 0 def circuit(args): @@ -252,7 +294,8 @@ def value(self): pathlen = 0 idx = 0 visited = set() - arr = [argval(a) for a in self.args] + arr = argvals(self.args) + while idx not in visited: if idx is None: return False @@ -287,15 +330,22 @@ def decompose(self): return [all(rev[x] == i for i, x in enumerate(fwd))], [] def value(self): - fwd = [argval(a) for a in self.args[0]] - rev = [argval(a) for a in self.args[1]] - return all(rev[x] == i for i, x in enumerate(fwd)) + fwd = argvals(self.args[0]) + rev = argvals(self.args[1]) + # args are fine, now evaluate actual inverse cons + try: + return all(rev[x] == i for i, x in enumerate(fwd)) + except IndexError: # partiality of Element constraint + return False class Table(GlobalConstraint): """The values of the variables in 'array' correspond to a row in 'table' """ def __init__(self, array, table): + array = flatlist(array) + if not all(isinstance(x, Expression) for x in array): + raise TypeError("the first argument of a Table constraint should only contain variables/expressions") super().__init__("table", [array, table]) def decompose(self): @@ -305,10 +355,32 @@ def decompose(self): def value(self): arr, tab = self.args - arrval = [argval(a) for a in arr] + arrval = argvals(arr) return arrval in tab +class NegativeTable(GlobalConstraint): + """The values of the variables in 'array' do not correspond to any row in 'table' + """ + def __init__(self, array, table): + array = flatlist(array) + if not all(isinstance(x, Expression) for x in array): + raise TypeError("the first argument of a Table constraint should only contain variables/expressions") + super().__init__("negative_table", [array, table]) + + def decompose(self): + from .python_builtins import all as cpm_all + from .python_builtins import any as cpm_any + arr, tab = self.args + return [cpm_all(cpm_any(ai != ri for ai, ri in zip(arr, row)) for row in tab)], [] + + def value(self): + arr, tab = self.args + arrval = argvals(arr) + tabval = argvals(tab) + return arrval not in tabval + + # syntax of the form 'if b then x == 9 else x == 0' is not supported (no override possible) # same semantic as CPLEX IfThenElse constraint # https://www.ibm.com/docs/en/icos/12.9.0?topic=methods-ifthenelse-method @@ -320,10 +392,13 @@ def __init__(self, condition, if_true, if_false): def value(self): condition, if_true, if_false = self.args - if argval(condition): - return argval(if_true) - else: - return argval(if_false) + try: + if argval(condition): + return argval(if_true) + else: + return argval(if_false) + except IncompleteFunctionError: + return False def decompose(self): condition, if_true, if_false = self.args @@ -341,8 +416,6 @@ class InDomain(GlobalConstraint): """ def __init__(self, expr, arr): - assert not (is_boolexpr(expr) or any(is_boolexpr(a) for a in arr)), \ - "The expressions in the InDomain constraint should not be boolean" super().__init__("InDomain", [expr, arr]) def decompose(self): @@ -371,7 +444,7 @@ def decompose(self): def value(self): - return argval(self.args[0]) in argval(self.args[1]) + return argval(self.args[0]) in argvals(self.args[1]) def __repr__(self): return "{} in {}".format(self.args[0], self.args[1]) @@ -401,7 +474,7 @@ def decompose(self): return decomp, [] def value(self): - return sum(argval(a) for a in self.args) % 2 == 1 + return sum(argvals(self.args)) % 2 == 1 def __repr__(self): if len(self.args) == 2: @@ -412,7 +485,8 @@ def __repr__(self): class Cumulative(GlobalConstraint): """ Global cumulative constraint. Used for resource aware scheduling. - Ensures no overlap between tasks and never exceeding the capacity of the resource + Ensures that the capacity of the resource is never exceeded + Equivalent to noOverlap when demand and capacity are equal to 1 Supports both varying demand across tasks or equal demand for all jobs """ def __init__(self, start, duration, end, demand, capacity): @@ -470,14 +544,14 @@ def decompose(self): return cons, [] def value(self): - argvals = [np.array([argval(a) for a in arg]) if is_any_list(arg) + arg_vals = [np.array(argvals(arg)) if is_any_list(arg) else argval(arg) for arg in self.args] - if any(a is None for a in argvals): + if any(a is None for a in arg_vals): return None # start, dur, end are np arrays - start, dur, end, demand, capacity = argvals + start, dur, end, demand, capacity = arg_vals # start and end seperated by duration if not (start + dur == end).all(): return False @@ -491,22 +565,96 @@ def value(self): return True +class Precedence(GlobalConstraint): + """ + Constraint enforcing some values have precedence over others. + Given an array of variables X and a list of precedences P: + Then in order to satisfy the constraint, if X[i] = P[j+1], then there exists a X[i'] = P[j] with i' < i + """ + def __init__(self, vars, precedence): + if not is_any_list(vars): + raise TypeError("Precedence expects a list of variables, but got", vars) + if not is_any_list(precedence) or any(isinstance(x, Expression) for x in precedence): + raise TypeError("Precedence expects a list of values as precedence, but got", precedence) + super().__init__("precedence", [vars, precedence]) + + def decompose(self): + """ + Decomposition based on: + Law, Yat Chiu, and Jimmy HM Lee. "Global constraints for integer and set value precedence." + Principles and Practice of Constraint Programming–CP 2004: 10th International Conference, CP 2004 + """ + from .python_builtins import any as cpm_any + + args, precedence = self.args + constraints = [] + for s,t in zip(precedence[:-1], precedence[1:]): + for j in range(len(args)): + constraints += [(args[j] == t).implies(cpm_any(args[:j] == s))] + return constraints, [] + + def value(self): + + args, precedence = self.args + vals = np.array(argvals(args)) + for s,t in zip(precedence[:-1], precedence[1:]): + if vals[0] == t: return False + for j in range(len(args)): + if vals[j] == t and sum(vals[:j] == s) == 0: + return False + return True + + +class NoOverlap(GlobalConstraint): + + def __init__(self, start, dur, end): + assert is_any_list(start), "start should be a list" + assert is_any_list(dur), "duration should be a list" + assert is_any_list(end), "end should be a list" + + start = flatlist(start) + dur = flatlist(dur) + end = flatlist(end) + assert len(start) == len(dur) == len(end), "Start, duration and end should have equal length in NoOverlap constraint" + + super().__init__("no_overlap", [start, dur, end]) + + def decompose(self): + start, dur, end = self.args + cons = [s + d == e for s,d,e in zip(start, dur, end)] + for (s1, e1), (s2, e2) in all_pairs(zip(start, end)): + cons += [(e1 <= s2) | (e2 <= s1)] + return cons, [] + def value(self): + start, dur, end = argvals(self.args) + if any(s + d != e for s,d,e in zip(start, dur, end)): + return False + for (s1,d1, e1), (s2,d2, e2) in all_pairs(zip(start,dur, end)): + if e1 > s2 and e2 > s1: + return False + return True + + class GlobalCardinalityCount(GlobalConstraint): """ GlobalCardinalityCount(vars,vals,occ): The number of occurrences of each value vals[i] in the list of variables vars must be equal to occ[i]. """ - def __init__(self, vars, vals, occ): + def __init__(self, vars, vals, occ, closed=False): flatargs = flatlist([vars, vals, occ]) if any(is_boolexpr(arg) for arg in flatargs): raise TypeError("Only numerical arguments allowed for gcc global constraint: {}".format(flatargs)) super().__init__("gcc", [vars,vals,occ]) + self.closed = closed def decompose(self): from .globalfunctions import Count vars, vals, occ = self.args - return [Count(vars, i) == v for i, v in zip(vals, occ)], [] + constraints = [Count(vars, i) == v for i, v in zip(vals, occ)] + if self.closed: + constraints += [InDomain(v, vals) for v in vars] + return constraints, [] def value(self): from .python_builtins import all @@ -514,6 +662,220 @@ def value(self): return all(decomposed).value() +class Increasing(GlobalConstraint): + """ + The "Increasing" constraint, the expressions will have increasing (not strictly) values + """ + + def __init__(self, *args): + super().__init__("increasing", flatlist(args)) + + def decompose(self): + """ + Returns two lists of constraints: + 1) the decomposition of the Increasing constraint + 2) empty list of defining constraints + """ + args = self.args + return [args[i] <= args[i+1] for i in range(len(args)-1)], [] + + def value(self): + args = argvals(self.args) + return all(args[i] <= args[i+1] for i in range(len(args)-1)) + + +class Decreasing(GlobalConstraint): + """ + The "Decreasing" constraint, the expressions will have decreasing (not strictly) values + """ + + def __init__(self, *args): + super().__init__("decreasing", flatlist(args)) + + def decompose(self): + """ + Returns two lists of constraints: + 1) the decomposition of the Decreasing constraint + 2) empty list of defining constraints + """ + args = self.args + return [args[i] >= args[i+1] for i in range(len(args)-1)], [] + + def value(self): + args = argvals(self.args) + return all(args[i] >= args[i+1] for i in range(len(args)-1)) + + +class IncreasingStrict(GlobalConstraint): + """ + The "IncreasingStrict" constraint, the expressions will have increasing (strictly) values + """ + + def __init__(self, *args): + super().__init__("strictly_increasing", flatlist(args)) + + def decompose(self): + """ + Returns two lists of constraints: + 1) the decomposition of the IncreasingStrict constraint + 2) empty list of defining constraints + """ + args = self.args + return [args[i] < args[i+1] for i in range(len(args)-1)], [] + + def value(self): + args = argvals(self.args) + return all(args[i] < args[i+1] for i in range(len(args)-1)) + + +class DecreasingStrict(GlobalConstraint): + """ + The "DecreasingStrict" constraint, the expressions will have decreasing (strictly) values + """ + + def __init__(self, *args): + super().__init__("strictly_decreasing", flatlist(args)) + + def decompose(self): + """ + Returns two lists of constraints: + 1) the decomposition of the DecreasingStrict constraint + 2) empty list of defining constraints + """ + args = self.args + return [(args[i] > args[i+1]) for i in range(len(args)-1)], [] + + def value(self): + args = argvals(self.args) + return all(args[i] > args[i+1] for i in range(len(args)-1)) + + +class LexLess(GlobalConstraint): + """ Given lists X,Y, enforcing that X is lexicographically less than Y. + """ + def __init__(self, list1, list2): + X = flatlist(list1) + Y = flatlist(list2) + if len(X) != len(Y): + raise CPMpyException(f"The 2 lists given in LexLess must have the same size: X length is {len(X)} and Y length is {len(Y)}") + super().__init__("lex_less", [X, Y]) + + def decompose(self): + """ + Implementation inspired by Hakan Kjellerstrand (http://hakank.org/cpmpy/cpmpy_hakank.py) + + The decomposition creates auxiliary Boolean variables and constraints that + collectively ensure X is lexicographically less than Y + The auxiliary boolean vars are defined to represent if the given lists are lexicographically ordered + (less or equal) up to the given index. + Decomposition enforces through the constraining part that the first boolean variable needs to be true, and thus + through the defining part it is enforced that if it is not strictly lexicographically less in a given index, + then next index must be lexicographically less or equal. It needs to be strictly less in at least one index. + + The use of auxiliary Boolean variables bvar ensures that the constraints propagate immediately, + maintaining arc-consistency. Each bvar[i] enforces the lexicographic ordering at each position, ensuring that + every value in the domain of X[i] can be extended to a consistent value in the domain of $Y_i$ for all + subsequent positions. + """ + X, Y = cpm_array(self.args) + + bvar = boolvar(shape=(len(X) + 1)) + + # Constraint ensuring that each element in X is less than or equal to the corresponding element in Y, + # until a strict inequality is encountered. + defining = [bvar == ((X <= Y) & ((X < Y) | bvar[1:]))] + # enforce the last element to be true iff (X[-1] < Y[-1]), enforcing strict lexicographic order + defining.append(bvar[-1] == (X[-1] < Y[-1])) + constraining = [bvar[0]] + + return constraining, defining + + def value(self): + X, Y = argvals(self.args) + return any((X[i] < Y[i]) & all(X[j] <= Y[j] for j in range(i)) for i in range(len(X))) + + +class LexLessEq(GlobalConstraint): + """ Given lists X,Y, enforcing that X is lexicographically less than Y (or equal). + """ + def __init__(self, list1, list2): + X = flatlist(list1) + Y = flatlist(list2) + if len(X) != len(Y): + raise CPMpyException(f"The 2 lists given in LexLessEq must have the same size: X length is {len(X)} and Y length is {len(Y)}") + super().__init__("lex_lesseq", [X, Y]) + + def decompose(self): + """ + Implementation inspired by Hakan Kjellerstrand (http://hakank.org/cpmpy/cpmpy_hakank.py) + + The decomposition creates auxiliary Boolean variables and constraints that + collectively ensure X is lexicographically less than Y + The auxiliary boolean vars are defined to represent if the given lists are lexicographically ordered + (less or equal) up to the given index. + Decomposition enforces through the constraining part that the first boolean variable needs to be true, and thus + through the defining part it is enforced that if it is not strictly lexicographically less in a given index, + then next index must be lexicographically less or equal. + + The use of auxiliary Boolean variables bvar ensures that the constraints propagate immediately, + maintaining arc-consistency. Each bvar[i] enforces the lexicographic ordering at each position, ensuring that + every value in the domain of X[i] can be extended to a consistent value in the domain of $Y_i$ for all + subsequent positions. + """ + X, Y = cpm_array(self.args) + + bvar = boolvar(shape=(len(X) + 1)) + defining = [bvar == ((X <= Y) & ((X < Y) | bvar[1:]))] + defining.append(bvar[-1] == (X[-1] <= Y[-1])) + constraining = [bvar[0]] + + return constraining, defining + + def value(self): + X, Y = argvals(self.args) + return any((X[i] < Y[i]) & all(X[j] <= Y[j] for j in range(i)) for i in range(len(X))) | all(X[i] == Y[i] for i in range(len(X))) + + +class LexChainLess(GlobalConstraint): + """ Given a matrix X, LexChainLess enforces that all rows are lexicographically ordered. + """ + def __init__(self, X): + # Ensure the numpy array is 2D + X = cpm_array(X) + assert X.ndim == 2, "Input must be a 2D array or a list of lists" + super().__init__("lex_chain_less", X.tolist()) + + def decompose(self): + """ Decompose to a series of LexLess constraints between subsequent rows + """ + X = self.args + return [LexLess(prev_row, curr_row) for prev_row, curr_row in zip(X, X[1:])], [] + + def value(self): + X = argvals(self.args) + return all(LexLess(prev_row, curr_row).value() for prev_row, curr_row in zip(X, X[1:])) + + +class LexChainLessEq(GlobalConstraint): + """ Given a matrix X, LexChainLessEq enforces that all rows are lexicographically ordered. + """ + def __init__(self, X): + # Ensure the numpy array is 2D + X = cpm_array(X) + assert X.ndim == 2, "Input must be a 2D array or a list of lists" + super().__init__("lex_chain_lesseq", X.tolist()) + + def decompose(self): + """ Decompose to a series of LexLessEq constraints between subsequent rows + """ + X = self.args + return [LexLessEq(prev_row, curr_row) for prev_row, curr_row in zip(X, X[1:])], [] + + def value(self): + X = argvals(self.args) + return all(LexLessEq(prev_row, curr_row).value() for prev_row, curr_row in zip(X, X[1:])) + + class DirectConstraint(Expression): """ A DirectConstraint will directly call a function of the underlying solver when added to a CPMpy solver diff --git a/cpmpy/expressions/globalfunctions.py b/cpmpy/expressions/globalfunctions.py index 16c050ce8..8110a261d 100644 --- a/cpmpy/expressions/globalfunctions.py +++ b/cpmpy/expressions/globalfunctions.py @@ -1,9 +1,10 @@ #!/usr/bin/env python -# -*- coding:utf-8 -*- +#-*- coding:utf-8 -*- ## ## globalfunctions.py ## """ + Global functions conveniently express numerical global constraints. Using global functions ------------------------ @@ -57,6 +58,8 @@ def decompose_comparison(self): Maximum Element Count + Among + NValue Abs """ @@ -66,7 +69,7 @@ def decompose_comparison(self): from ..exceptions import CPMpyException, IncompleteFunctionError, TypeError from .core import Expression, Operator, Comparison from .variables import boolvar, intvar, cpm_array, _NumVarImpl -from .utils import flatlist, all_pairs, argval, is_num, eval_comparison, is_any_list, is_boolexpr, get_bounds +from .utils import flatlist, all_pairs, argval, is_num, eval_comparison, is_any_list, is_boolexpr, get_bounds, argvals class GlobalFunction(Expression): @@ -130,10 +133,6 @@ def decompose_comparison(self, cpm_op, cpm_rhs): they should be enforced toplevel. """ from .python_builtins import any, all - if cpm_op == "==": # can avoid creating aux var - return [any(x <= cpm_rhs for x in self.args), - all(x >= cpm_rhs for x in self.args)], [] - lb, ub = self.get_bounds() _min = intvar(lb, ub) return [eval_comparison(cpm_op, _min, cpm_rhs)], \ @@ -171,10 +170,6 @@ def decompose_comparison(self, cpm_op, cpm_rhs): they should be enforced toplevel. """ from .python_builtins import any, all - if cpm_op == "==": # can avoid creating aux var here - return [any(x >= cpm_rhs for x in self.args), - all(x <= cpm_rhs for x in self.args)], [] - lb, ub = self.get_bounds() _max = intvar(lb, ub) return [eval_comparison(cpm_op, _max, cpm_rhs)], \ @@ -256,7 +251,8 @@ def value(self): if idxval is not None: if idxval >= 0 and idxval < len(arr): return argval(arr[idxval]) - raise IncompleteFunctionError(f"Index {idxval} out of range for array of length {len(arr)} while calculating value for expression {self}") + raise IncompleteFunctionError(f"Index {idxval} out of range for array of length {len(arr)} while calculating value for expression {self}" + + "\n Use argval(expr) to get the value of expr with relational semantics.") return None # default def decompose_comparison(self, cpm_op, cpm_rhs): @@ -316,3 +312,135 @@ def get_bounds(self): """ arr, val = self.args return 0, len(arr) + + + +class Among(GlobalFunction): + """ + The Among (numerical) global constraint represents the number of variable that take values among the values in arr + """ + + def __init__(self,arr,vals): + if not is_any_list(arr) or not is_any_list(vals): + raise TypeError("Among takes as input two arrays, not: {} and {}".format(arr,vals)) + if any(isinstance(val, Expression) for val in vals): + raise TypeError(f"Among takes a set of values as input, not {vals}") + super().__init__("among", [arr,vals]) + + def decompose_comparison(self, cmp_op, cmp_rhs): + """ + Among(arr, vals) can only be decomposed if it's part of a comparison' + """ + from .python_builtins import sum, any + arr, values = self.args + count_for_each_val = [Count(arr, val) for val in values] + return [eval_comparison(cmp_op, sum(count_for_each_val), cmp_rhs)], [] + + def value(self): + return int(sum(np.isin(argvals(self.args[0]), self.args[1]))) + + def get_bounds(self): + return 0, len(self.args[0]) + + +class NValue(GlobalFunction): + + """ + The NValue constraint counts the number of distinct values in a given set of variables. + """ + + def __init__(self, arr): + if not is_any_list(arr): + raise ValueError("NValue takes an array as input") + super().__init__("nvalue", arr) + + def decompose_comparison(self, cmp_op, cpm_rhs): + """ + NValue(arr) can only be decomposed if it's part of a comparison + + Based on "simple decomposition" from: + Bessiere, Christian, et al. "Decomposition of the NValue constraint." + International Conference on Principles and Practice of Constraint Programming. + Berlin, Heidelberg: Springer Berlin Heidelberg, 2010. + """ + from .python_builtins import sum, any + + lbs, ubs = get_bounds(self.args) + lb, ub = min(lbs), max(ubs) + + constraints = [] + + # introduce boolvar for each possible value + bvars = boolvar(shape=(ub+1-lb)) + + args = cpm_array(self.args) + # bvar is true if the value is taken by any variable + for bv, val in zip(bvars, range(lb, ub+1)): + constraints += [any(args == val) == bv] + + return [eval_comparison(cmp_op, sum(bvars), cpm_rhs)], constraints + + def value(self): + return len(set(argval(a) for a in self.args)) + + def get_bounds(self): + """ + Returns the bounds of the (numerical) global constraint + """ + return 1, len(self.args) + + +class NValueExcept(GlobalFunction): + + """ + The NValueExceptN constraint counts the number of distinct values, + not including value N, if any argument is assigned to it. + """ + + def __init__(self, arr, n): + if not is_any_list(arr): + raise ValueError("NValueExcept takes an array as input") + if not is_num(n): + raise ValueError(f"NValueExcept takes an integer as second argument, but got {n} of type {type(n)}") + super().__init__("nvalue_except",[arr, n]) + + def decompose_comparison(self, cmp_op, cpm_rhs): + """ + NValue(arr) can only be decomposed if it's part of a comparison + + Based on "simple decomposition" from: + Bessiere, Christian, et al. "Decomposition of the NValue constraint." + International Conference on Principles and Practice of Constraint Programming. + Berlin, Heidelberg: Springer Berlin Heidelberg, 2010. + """ + from .python_builtins import sum, any + + arr, n = self.args + arr = cpm_array(arr) + lbs, ubs = get_bounds(arr) + lb, ub = min(lbs), max(ubs) + + constraints = [] + + # introduce boolvar for each possible value + bvars = boolvar(shape=(ub + 1 - lb)) + idx_of_n = n - lb + if 0 <= idx_of_n < len(bvars): + count_of_vals = sum(bvars[:idx_of_n]) + sum(bvars[idx_of_n+1:]) + else: + count_of_vals = sum(bvars) + + # bvar is true if the value is taken by any variable + for bv, val in zip(bvars, range(lb, ub + 1)): + constraints += [any(arr == val) == bv] + + return [eval_comparison(cmp_op, count_of_vals, cpm_rhs)], constraints + + def value(self): + return len(set(argval(a) for a in self.args[0]) - {self.args[1]}) + + def get_bounds(self): + """ + Returns the bounds of the (numerical) global constraint + """ + return 0, len(self.args) \ No newline at end of file diff --git a/cpmpy/expressions/utils.py b/cpmpy/expressions/utils.py index 9b0310773..ecd865fe2 100644 --- a/cpmpy/expressions/utils.py +++ b/cpmpy/expressions/utils.py @@ -12,27 +12,35 @@ .. autosummary:: :nosignatures: + is_bool is_int is_num + is_false_cst + is_true_cst + is_boolexpr is_pure_list is_any_list + is_transition flatlist all_pairs argval + argvals eval_comparison + get_bounds """ import numpy as np import math -from collections.abc import Iterable # for _flatten -from itertools import chain, combinations +from collections.abc import Iterable # for flatten +from itertools import combinations from cpmpy.exceptions import IncompleteFunctionError def is_bool(arg): """ is it a boolean (incl numpy variants) """ - return isinstance(arg, (bool, np.bool_)) + from cpmpy import BoolVal + return isinstance(arg, (bool, np.bool_, BoolVal)) def is_int(arg): @@ -120,11 +128,21 @@ def argval(a): We check with hasattr instead of isinstance to avoid circular dependency """ - try: - return a.value() if hasattr(a, "value") else a - except IncompleteFunctionError as e: - if a.is_bool(): return False - raise e + if hasattr(a, "value"): + try: + return a.value() + except IncompleteFunctionError as e: + if a.is_bool(): + return False + else: + raise e + return a + + +def argvals(arr): + if is_any_list(arr): + return [argvals(arg) for arg in arr] + return argval(arr) def eval_comparison(str_op, lhs, rhs): diff --git a/cpmpy/expressions/variables.py b/cpmpy/expressions/variables.py index 31c1ab33b..d39d70103 100644 --- a/cpmpy/expressions/variables.py +++ b/cpmpy/expressions/variables.py @@ -24,9 +24,18 @@ Boolean and Integer decision variables are the key elements of a CP model. - All variables in CPMpy are n-dimensional array objects and have defined dimensions. Following the numpy library, the dimension sizes of an n-dimenionsal array is called its __shape__. In CPMpy all variables are considered an array with a given shape. For 'single' variables the shape is '1'. For an array of length `n` the shape is 'n'. An `n*m` matrix has shape (n,m), and tensors with more than 2 dimensions are all supported too. For the implementation of this, CPMpy builts on numpy's n-dimensional ndarray and inherits many of its benefits (vectorized operators and advanced indexing). + All variables in CPMpy are n-dimensional array objects and have defined dimensions. + Following the numpy library, the dimension sizes of an n-dimenionsal array is called its __shape__. + In CPMpy all variables are considered an array with a given shape. For 'single' variables the shape + is '1'. For an array of length `n` the shape is 'n'. An `n*m` matrix has shape (n,m), and tensors + with more than 2 dimensions are all supported too. For the implementation of this, + CPMpy builts on numpy's n-dimensional ndarray and inherits many of its benefits + (vectorized operators and advanced indexing). - This module contains the cornerstone `boolvar()` and `intvar()` functions, which create (numpy arrays of) variables. There is also a helper function `cpm_array` for wrapping standard numpy arrays so they can be indexed by a variable. Apart from these 3 functions, none of the classes in this module should be directly created; they are created by these 3 helper functions. + This module contains the cornerstone `boolvar()` and `intvar()` functions, which create (numpy arrays of) + variables. There is also a helper function `cpm_array()` for wrapping standard numpy arrays so they can be + indexed by a variable. Apart from these 3 functions, none of the classes in this module should be directly + instantiated; they are created by these 3 helper functions. =============== @@ -126,7 +135,7 @@ def IntVar(lb, ub, shape=1, name=None): def intvar(lb, ub, shape=1, name=None): """ - Integer decision variables are constructed by specifying the lowest (lb) + Integer decision variables are constructed by specifying the lowest (lb) value the decision variable can take, as well as the highest value (ub). Arguments: diff --git a/cpmpy/model.py b/cpmpy/model.py index 147ac4f2b..a0249c4a7 100644 --- a/cpmpy/model.py +++ b/cpmpy/model.py @@ -15,9 +15,9 @@ See the examples for basic usage, which involves: - - creation, e.g. `m = Model(cons, minimize=obj)` - - solving, e.g. `m.solve()` - - optionally, checking status/runtime, e.g. `m.status()` + - creation, e.g. `m = Model(cons, minimize=obj)` + - solving, e.g. `m.solve()` + - optionally, checking status/runtime, e.g. `m.status()` =============== List of classes @@ -266,7 +266,7 @@ def from_file(fname): def copy(self): """ Makes a shallow copy of the model. - Constraints and variables are shared among the original and copied model. + Constraints and variables are shared among the original and copied model (references to the same Expression objects). The /list/ of constraints itself is different, so adding or removing constraints from one model does not affect the other. """ if self.objective_is_min: return Model(self.constraints, minimize=self.objective_) @@ -278,4 +278,4 @@ def copy(self): # keep for backwards compatibility def deepcopy(self, memodict={}): warnings.warn("Deprecated, use copy.deepcopy() instead, will be removed in stable version", DeprecationWarning) - return copy.deepcopy(self, memodict) \ No newline at end of file + return copy.deepcopy(self, memodict) diff --git a/cpmpy/solvers/TEMPLATE.py b/cpmpy/solvers/TEMPLATE.py index f7fd0458b..6ea5d5561 100644 --- a/cpmpy/solvers/TEMPLATE.py +++ b/cpmpy/solvers/TEMPLATE.py @@ -6,12 +6,23 @@ The functions are ordered in a way that could be convenient to start from the top and continue in that order + After you are done filling in the template, remove all comments starting with [GUIDELINE] + WARNING: do not include the python package at the top of the file, as CPMpy should also work without this solver installed. To ensure that, include it inside supported() and other functions that need it... """ -from ..transformations.decompose_global import decompose_in_tree + +from .solver_interface import SolverInterface, SolverStatus, ExitStatus +from ..expressions.core import Expression, Comparison, Operator +from ..expressions.variables import _BoolVarImpl, NegBoolView, _IntVarImpl, _NumVarImpl +from ..expressions.utils import is_num, is_any_list, is_boolexpr +from ..transformations.get_variables import get_variables from ..transformations.normalize import toplevel_list +from ..transformations.decompose_global import decompose_in_tree +from ..transformations.flatten_model import flatten_constraint +from ..transformations.comparison import only_numexpr_equality +from ..transformations.reification import reify_rewrite, only_bv_reifies """ Interface to TEMPLATE's API @@ -30,12 +41,6 @@ CPM_template """ -from .solver_interface import SolverInterface, SolverStatus, ExitStatus -from ..expressions.core import Expression, Comparison, Operator -from ..expressions.variables import _BoolVarImpl, NegBoolView, _IntVarImpl -from ..expressions.utils import is_num, is_any_list -from ..transformations.get_variables import get_variables -from ..transformations.flatten_model import flatten_constraint class CPM_template(SolverInterface): """ @@ -74,12 +79,16 @@ def __init__(self, cpm_model=None, subsolver=None): import TEMPLATEpy - assert(subsolver is None) # unless you support subsolvers, see pysat or minizinc + assert subsolver is None # unless you support subsolvers, see pysat or minizinc # initialise the native solver object - self.tpl_model = TEMPLATEpy.Model("cpmpy") + # [GUIDELINE] we commonly use 3-letter abbrivations to refer to native objects: + # OR-tools uses ort_solver, Gurobi grb_solver, Exact xct_solver... + self.TPL_solver = TEMPLATEpy.Solver("cpmpy") # initialise everything else and post the constraints/objective + # [GUIDELINE] this superclass call should happen AFTER all solver-native objects are created. + # internally, the constructor relies on __add__ which uses the above solver native object(s) super().__init__(name="TEMPLATE", cpm_model=cpm_model) @@ -92,22 +101,31 @@ def solve(self, time_limit=None, **kwargs): - kwargs: any keyword argument, sets parameters of solver object Arguments that correspond to solver parameters: - - + # [GUIDELINE] Please document key solver arguments that the user might wish to change + # for example: assumptions=[x,y,z], log_output=True, var_ordering=3, num_cores=8, ... + # [GUIDELINE] Add link to documentation of all solver parameters """ + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + if time_limit is not None: - raise NotImplementedError("TEMPLATE: TODO, implement time_limit") + self.TPL_solver.set_timelimit_seconds(time_limit) + + # [GUIDELINE] if your solver supports solving under assumptions, add `assumptions` as argument in header + # e.g., def solve(self, time_limit=None, assumptions=None, **kwargs): + # then translate assumptions here; assumptions are a list of Boolean variables or NegBoolViews # call the solver, with parameters - my_status = self.TEMPLATE_solver.solve(**kwargs) + my_status = self.TPL_solver.solve(**kwargs) + # [GUIDELINE] consider saving the status as self.TPL_status so that advanced CPMpy users can access the status object. + # This is mainly useful when more elaborate information about the solve-call is saved into the status # new status, translate runtime self.cpm_status = SolverStatus(self.name) - self.cpm_status.runtime = self.TEMPLATE_solver.time() + self.cpm_status.runtime = self.TPL_solver.time() # wallclock time in (float) seconds - # translate exit status + # translate solver exit status to CPMpy exit status if my_status is True: self.cpm_status.exitstatus = ExitStatus.FEASIBLE elif my_status is False: @@ -127,13 +145,12 @@ def solve(self, time_limit=None, **kwargs): # fill in variable values for cpm_var in self.user_vars: sol_var = self.solver_var(cpm_var) - cpm_var._value = None # if not in solution - #cpm_var._value = self.TEMPLATEpy.value(sol_var) + cpm_var._value = self.TPL_solver.value(sol_var) raise NotImplementedError("TEMPLATE: back-translating the solution values") # translate objective, for optimisation problems only if self.has_objective(): - self.objective_value_ = self.TEMPLATE_solver.ObjectiveValue() + self.objective_value_ = self.TPL_solver.ObjectiveValue() return has_sol @@ -146,6 +163,9 @@ def solver_var(self, cpm_var): if is_num(cpm_var): # shortcut, eases posting constraints return cpm_var + # [GUIDELINE] some solver interfaces explicitely create variables on a solver object + # then use self.TPL_solver.NewBoolVar(...) instead of TEMPLATEpy.NewBoolVar(...) + # special case, negative-bool-view # work directly on var inside the view if isinstance(cpm_var, NegBoolView): @@ -165,7 +185,7 @@ def solver_var(self, cpm_var): return self._varmap[cpm_var] - # if TEMPLATE does not support objective functions, you can delete objective()/_make_numexpr() + # [GUIDELINE] if TEMPLATE does not support objective functions, you can delete this function definition def objective(self, expr, minimize=True): """ Post the given expression to the solver as objective to minimize/maximize @@ -183,21 +203,30 @@ def objective(self, expr, minimize=True): # make objective function or variable and post obj = self._make_numexpr(flat_obj) + # [GUIDELINE] if the solver interface does not provide a solver native "numeric expression" object, + # _make_numexpr may be removed and an objective can be posted as: + # self.TPL_solver.MinimizeWeightedSum(obj.args[0], self.solver_vars(obj.args[1]) or similar + if minimize: - TEMPLATEpy.Minimize(obj) + self.TPL_solver.Minimize(obj) else: - TEMPLATEpy.Maximize(obj) + self.TPL_solver.Maximize(obj) def has_objective(self): - return TEMPLATEpy.hasObjective() + return self.TPL_solver.hasObjective() def _make_numexpr(self, cpm_expr): """ - Turns a numeric CPMpy 'flat' expression into a solver-specific - numeric expression + Converts a numeric CPMpy 'flat' expression into a solver-specific numeric expression - Used especially to post an expression as objective function + Primarily used for setting objective functions, and optionally in constraint posting """ + + # [GUIDELINE] not all solver interfaces have a native "numerical expression" object. + # in that case, this function may be removed and a case-by-case analysis of the numerical expression + # used in the constraint at hand is required in __add__ + # For an example of such solver interface, check out solvers/choco.py or solvers/exact.py + if is_num(cpm_expr): return cpm_expr @@ -205,22 +234,18 @@ def _make_numexpr(self, cpm_expr): if isinstance(cpm_expr, _NumVarImpl): # _BoolVarImpl is subclass of _NumVarImpl return self.solver_var(cpm_expr) - # if the solver only supports a decision variable as argument to its minimize/maximize then - # well, smth generic like - #(obj_var, obj_cons) = get_or_make_var(cpm_expr) - #self += obj_cons - #return self.solver_var(obj_var) - - # else if the solver support e.g. a linear expression as objective, built it here - # something like - #if isinstance(cpm_expr, Operator): - # if cpm_expr.name == 'sum': - # return sum(self.solver_vars(cpm_expr.args)) # if TEMPLATEpy supports this - # elif cpm_expr.name == 'wsum': - # w = cpm_expr.args[0] - # x = self.solver_vars(cpm_expr.args[1]) - # return sum(wi*xi for wi,xi in zip(w,x)) # if TEMPLATEpy supports this - + # any solver-native numerical expression + if isinstance(cpm_expr, Operator): + if cpm_expr.name == 'sum': + return self.TPL_solver.sum(self.solver_vars(cpm_expr.args)) + elif cpm_expr.name == 'wsum': + weights, vars = cpm_expr.args + return self.TPL_solver.weighted_sum(weights, self.solver_vars(vars)) + # [GUIDELINE] or more fancy ones such as max + # be aware this is not the Maximum CONSTRAINT, but rather the Maximum NUMERICAL EXPRESSION + elif cpm_expr.name == "max": + return self.TPL_solver.maximum_of_vars(self.solver_vars(cpm_expr.args)) + # ... raise NotImplementedError("TEMPLATE: Not a known supported numexpr {}".format(cpm_expr)) @@ -242,10 +267,11 @@ def transform(self, cpm_expr): # apply transformations # XXX chose the transformations your solver needs, see cpmpy/transformations/ cpm_cons = toplevel_list(cpm_expr) - cpm_cons = decompose_in_tree(cpm_cons, supported={"AllDifferent"}) - cpm_cons = flatten_constraint(cpm_expr) # flat normal form - #cpm_cons = reify_rewrite(cpm_cons, supported=frozenset(['sum', 'wsum'])) # constraints that support reification - #cpm_cons = only_numexpr_equality(cpm_cons, supported=frozenset(["sum", "wsum", "sub"])) # supports >, <, != + cpm_cons = decompose_in_tree(cpm_cons, supported={"alldifferent"}) + cpm_cons = flatten_constraint(cpm_cons) # flat normal form + cpm_cons = reify_rewrite(cpm_cons, supported=frozenset(['sum', 'wsum'])) # constraints that support reification + cpm_cons = only_bv_reifies(cpm_cons) + cpm_cons = only_numexpr_equality(cpm_cons, supported=frozenset(["sum", "wsum", "sub"])) # supports >, <, != # ... return cpm_cons @@ -274,18 +300,72 @@ def __add__(self, cpm_expr_orig): # transform and post the constraints for cpm_expr in self.transform(cpm_expr_orig): - if isinstance(cpm_con, _BoolVarImpl): + if isinstance(cpm_expr, _BoolVarImpl): # base case, just var or ~var - self.TEMPLATE_solver.add_clause([ self.solver_var(cpm_con) ]) - elif cpm_con.name == 'or': - self.TEMPLATE_solver.add_clause(self.solver_vars(cpm_con.args)) - elif hasattr(cpm_expr, 'decompose'): - # global constraint not known, try posting generic decomposition - # side-step `__add__()` as the decomposition can contain non-user (auxiliary) variables - for con in self.transform(cpm_expr.decompose()): - self._post_constraint(con) + self.TPL_solver.add_clause([ self.solver_var(cpm_expr) ]) + + elif isinstance(cpm_expr, Operator): + if cpm_expr.name == "or": + self.TPL_solver.add_clause(self.solver_vars(cpm_expr.args)) + elif cpm_expr.name == "->": # half-reification + bv, subexpr = cpm_expr.args + # [GUIDELINE] example code for a half-reified sum/wsum comparison e.g. BV -> sum(IVs) >= 5 + if isinstance(subexpr, Comparison): + lhs, rhs = subexpr.args + if isinstance(lhs, _NumVarImpl) or (isinstance(lhs, Operator) and lhs.name in {"sum", "wsum"}): + TPL_lhs = self._make_numexpr(lhs) + self.TPL_solver.add_half_reified_comparison(self.solver_var(bv), + TPL_lhs, subexpr.name, self.solver_var(rhs)) + else: + raise NotImplementedError("TEMPLATE: no support for half-reified comparison:", subexpr) + else: + raise NotImplementedError("TEMPLATE: no support for half-reified constraint:", subexpr) + + elif isinstance(cpm_expr, Comparison): + lhs, rhs = cpm_expr.args + + # [GUIDELINE] == is used for both double reification and numerical comparisons + # need case by case analysis here. Note that if your solver does not support full-reification, + # you can rely on the transformation only_implies to convert all reifications to half-reification + # for more information, please reach out on github! + if cpm_expr.name == "==" and is_boolexpr(lhs) and is_boolexpr(rhs): # reification + bv, subexpr = lhs, rhs + assert isinstance(lhs, _BoolVarImpl), "lhs of reification should be var because of only_bv_reifies" + + if isinstance(subexpr, Comparison): + lhs, rhs = subexpr.args + if isinstance(lhs, _NumVarImpl) or (isinstance(lhs, Operator) and lhs.name in {"sum", "wsum"}): + TPL_lhs = self._make_numexpr(lhs) + self.TPL_solver.add_reified_comparison(self.solver_var(bv), + TPL_lhs, subexpr.name, self.solver_var(rhs)) + else: + raise NotImplementedError("TEMPLATE: no support for reified comparison:", subexpr) + else: + raise NotImplementedError("TEMPLATE: no support for reified constraint:", subexpr) + + # otherwise, numerical comparisons + if isinstance(lhs, _NumVarImpl) or (isinstance(lhs, Operator) and lhs.name in {"sum", "wsum"}): + TPL_lhs = self._make_numexpr(lhs) + self.TPL_solver.add_comparison(TPL_lhs, cpm_expr.name, self.solver_var(rhs)) + # global functions + elif cpm_expr.name == "==": + TPL_rhs = self.solver_var(rhs) + if lhs.name == "max": + self.TPL_solver.add_max_constraint(self.solver_vars(lhs), TPL_rhs) + elif lhs.name == "element": + TPL_arr, TPL_idx = self.solver_vars(lhs.args) + self.TPL_solver.add_element_constraint(TPL_arr, TPL_idx, TPL_rhs) + # elif... + else: + raise NotImplementedError("TEMPLATE: unknown equality constraint:", cpm_expr) + else: + raise NotImplementedError("TEMPLATE: unknown comparison constraint", cpm_expr) - raise NotImplementedError("TEMPLATE: constraint not (yet) supported", cpm_con) + # global constraints + elif cpm_expr.name == "alldifferent": + self.TPL_solver.add_alldifferent(self.solver_vars(cpm_expr.args)) + else: + raise NotImplementedError("TEMPLATE: constraint not (yet) supported", cpm_expr) return self @@ -311,7 +391,7 @@ def solveAll(self, display=None, time_limit=None, solution_limit=None, call_from # check if objective function if self.has_objective(): - raise Exception("TEMPLATE does not support finding all optimal solutions") + raise NotSupportedError("TEMPLATE does not support finding all optimal solutions") # A. Example code if solver supports callbacks if is_any_list(display): @@ -320,12 +400,12 @@ def solveAll(self, display=None, time_limit=None, solution_limit=None, call_from callback = display self.solve(time_limit, callback=callback, enumerate_all_solutions=True, **kwargs) - return self.TEMPLATE_solver.SolutionCount() + return self.TPL_solver.SolutionCount() # B. Example code if solver does not support callbacks self.solve(time_limit, enumerate_all_solutions=True, **kwargs) solution_count = 0 - for solution in self.TEMPLATE_solver.GetSolutions(): + for solution in self.TPL_solver.GetAllSolutions(): solution_count += 1 # Translate solution to variables for cpm_var in self.user_vars: diff --git a/cpmpy/solvers/__init__.py b/cpmpy/solvers/__init__.py index 376416057..06fe6dc41 100644 --- a/cpmpy/solvers/__init__.py +++ b/cpmpy/solvers/__init__.py @@ -2,8 +2,8 @@ CPMpy interfaces to (the Python API interface of) solvers Solvers typically use some of the generic transformations in - `transformations` as well as specific reformulations to map the - CPMpy expression to the solver's Python API + :mod:`cpmpy.transformations` as well as specific reformulations to map the + CPMpy expression to the solver's Python API. ================== List of submodules @@ -18,6 +18,7 @@ pysdd z3 exact + choco utils =============== @@ -33,6 +34,7 @@ CPM_pysdd CPM_z3 CPM_exact + CPM_choco ================= List of functions @@ -51,3 +53,5 @@ from .pysdd import CPM_pysdd from .z3 import CPM_z3 from .exact import CPM_exact +from .choco import CPM_choco + diff --git a/cpmpy/solvers/choco.py b/cpmpy/solvers/choco.py new file mode 100644 index 000000000..1a10f878b --- /dev/null +++ b/cpmpy/solvers/choco.py @@ -0,0 +1,594 @@ +#!/usr/bin/env python +#-*- coding:utf-8 -*- +## +## choco.py +## +""" + Interface to Choco solver's Python API + + Requires that the 'pychoco' python package is installed: + + $ pip install pychoco + + + Documentation of the solver's own Python API: + https://pypi.org/project/pychoco/ + https://pychoco.readthedocs.io/en/latest/ + + =============== + List of classes + =============== + + .. autosummary:: + :nosignatures: + + CPM_choco + + ============== + Module details + ============== +""" +import time + +import numpy as np + +from ..transformations.normalize import toplevel_list +from .solver_interface import SolverInterface, SolverStatus, ExitStatus +from ..expressions.core import Expression, Comparison, Operator, BoolVal +from ..expressions.globalconstraints import DirectConstraint +from ..expressions.variables import _NumVarImpl, _IntVarImpl, _BoolVarImpl, NegBoolView, intvar +from ..expressions.globalconstraints import GlobalConstraint +from ..expressions.utils import is_num, is_int, is_boolexpr, is_any_list, get_bounds, argval, argvals +from ..transformations.decompose_global import decompose_in_tree +from ..transformations.get_variables import get_variables +from ..transformations.flatten_model import flatten_constraint, flatten_objective +from ..transformations.comparison import only_numexpr_equality +from ..transformations.linearize import canonical_comparison +from ..transformations.reification import only_bv_reifies, reify_rewrite +from ..exceptions import ChocoBoundsException, ChocoTypeException, NotSupportedError + + +class CPM_choco(SolverInterface): + """ + Interface to the Choco solver python API + + Requires that the 'pychoco' python package is installed: + $ pip install pychoco + + See detailed installation instructions at: + https://pypi.org/project/pychoco/ + https://pychoco.readthedocs.io/en/latest/ + + Creates the following attributes (see parent constructor for more): + chc_model: the pychoco.Model() created by _model() + chc_solver: the choco Model().get_solver() instance used in solve() + + """ + + @staticmethod + def supported(): + # try to import the package + try: + import pychoco as chc + return True + except ImportError: + return False + + def __init__(self, cpm_model=None, subsolver=None): + """ + Constructor of the native solver object + + Requires a CPMpy model as input, and will create the corresponding + choco model and solver object (chc_model and chc_solver) + + chc_model and chc_solver can both be modified externally before + calling solve(), a prime way to use more advanced solver features + + Arguments: + - cpm_model: Model(), a CPMpy Model() (optional) + - subsolver: None + """ + if not self.supported(): + raise Exception("Install the python 'pychoco' package to use this solver interface") + + import pychoco as chc + + assert (subsolver is None), "Choco does not support any subsolver" + + # initialise the native solver objects + self.chc_model = chc.Model() + + # for the objective + self.obj = None + self.minimize_obj = None + self.helper_var = None + # for solving with assumption variables, TO-CHECK + + # initialise everything else and post the constraints/objective + super().__init__(name="choco", cpm_model=cpm_model) + + def solve(self, time_limit=None, **kwargs): + """ + Call the Choco solver + + Arguments: + - time_limit: maximum solve time in seconds (float, optional) + - kwargs: any keyword argument, sets parameters of solver object + + """ + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + + # call the solver, with parameters + self.chc_solver = self.chc_model.get_solver() + + start = time.time() + + if time_limit is not None: + self.chc_solver.limit_time(str(time_limit) + "s") + + if self.has_objective(): + sol = self.chc_solver.find_optimal_solution(maximize= not self.minimize_obj, + objective=self.solver_var(self.obj), + **kwargs) + else: + sol = self.chc_solver.find_solution() + end = time.time() + + # new status, get runtime + self.cpm_status = SolverStatus(self.name) + self.cpm_status.runtime = end - start + + # translate exit status + if sol is not None: + if time_limit is None or self.cpm_status.runtime < time_limit: # solved to optimality + self.cpm_status.exitstatus = ExitStatus.OPTIMAL + else: # solved, but optimality not proven + self.cpm_status.exitstatus = ExitStatus.FEASIBLE + elif time_limit is None or self.cpm_status.runtime < time_limit: # proven unsat + self.cpm_status.exitstatus = ExitStatus.UNSATISFIABLE + else: + self.cpm_status.exitstatus = ExitStatus.UNKNOWN # can happen when timeout is reached... + + + # True/False depending on self.chc_status + has_sol = sol is not None + + # translate solution values (of user specified variables only) + self.objective_value_ = None + if has_sol: + # fill in variable values + for cpm_var in self.user_vars: + value = sol.get_int_val(self.solver_var(cpm_var)) + if isinstance(cpm_var, _BoolVarImpl): + cpm_var._value = bool(value) + else: + cpm_var._value = value + + # translate objective + if self.has_objective(): + self.objective_value_ = sol.get_int_val(self.solver_var(self.obj)) + + return has_sol + + def solveAll(self, display=None, time_limit=None, solution_limit=None, call_from_model=False, **kwargs): + """ + Compute all (optimal) solutions, map them to CPMpy and optionally display the solutions. + + Arguments: + - display: either a list of CPMpy expressions, OR a callback function, called with the variables after value-mapping + default/None: nothing displayed + - solution_limit: stop after this many solutions (default: None) + - time_limit: maximum solve time in seconds (float, default: None) + + Returns: number of solutions found + """ + + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + + if time_limit is not None: + self.chc_solver.limit_time(str(time_limit) + "s") + + self.chc_solver = self.chc_model.get_solver() + start = time.time() + if self.has_objective(): + sols = self.chc_solver.find_all_optimal_solutions(maximize=not self.minimize_obj, + solution_limit=solution_limit, + objective=self.solver_var(self.obj), + **kwargs) + else: + sols = self.chc_solver.find_all_solutions(solution_limit=solution_limit) + end = time.time() + + # new status, get runtime + self.cpm_status = SolverStatus(self.name) + self.cpm_status.runtime = end - start + + # display if needed + if display is not None: + for sol in sols: + # map the solution to user vars + for cpm_var in self.user_vars: + value = sol.get_int_val(self.solver_var(cpm_var)) + if isinstance(cpm_var, _BoolVarImpl): + cpm_var._value = bool(value) + else: + cpm_var._value = value + # print the desired display + if isinstance(display, Expression): + print(argval(display)) + elif isinstance(display, list): + print(argvals(display)) + else: + display() # callback + + return len(sols) + + def solver_var(self, cpm_var): + """ + Creates solver variable for cpmpy variable + or returns from cache if previously created + """ + if is_num(cpm_var): # shortcut, eases posting constraints + if not is_int(cpm_var): + raise ValueError(f"Choco only accepts integer constants, got {cpm_var} of type {type(cpm_var)}") + if cpm_var < -2147483646 or cpm_var > 2147483646: + raise ChocoBoundsException( + "Choco does not accept integer literals with bounds outside of range (-2147483646..2147483646)") + return int(cpm_var) + + # special case, negative-bool-view + # work directly on var inside the view + if isinstance(cpm_var, NegBoolView): + return self.chc_model.bool_not_view(self.solver_var(cpm_var._bv)) + + # create if it does not exist + if cpm_var not in self._varmap: + if isinstance(cpm_var, _BoolVarImpl): + revar = self.chc_model.boolvar(name=str(cpm_var.name)) + elif isinstance(cpm_var, _IntVarImpl): + if cpm_var.lb < -2147483646 or cpm_var.ub > 2147483646: + raise ChocoBoundsException( + "Choco does not accept variables with bounds outside of range (-2147483646..2147483646)") + revar = self.chc_model.intvar(cpm_var.lb, cpm_var.ub, name=str(cpm_var.name)) + else: + raise NotImplementedError("Not a known var {}".format(cpm_var)) + self._varmap[cpm_var] = revar + + return self._varmap[cpm_var] + + def objective(self, expr, minimize): + """ + Post the given expression to the solver as objective to minimize/maximize + + - expr: Expression, the CPMpy expression that represents the objective function + - minimize: Bool, whether it is a minimization problem (True) or maximization problem (False) + + 'objective()' can be called multiple times, only the last one is stored + + (technical side note: constraints created during conversion of the objective + are premanently posted to the solver. Choco accepts variables to maximize or minimize + so it is needed to post constraints and create auxiliary variables) + """ + + # make objective function non-nested + obj_var = intvar(*get_bounds(expr)) + self += obj_var == expr + + self.obj = obj_var + self.minimize_obj = minimize # Choco has as default to maximize + + def has_objective(self): + return self.obj is not None + + + def _to_var(self, val): + from pychoco.variables.intvar import IntVar + if is_int(val): + # Choco accepts only int32, not int64 + if val < -2147483646 or val > 2147483646: + raise ChocoBoundsException( + "Choco does not accept integer literals with bounds outside of range (-2147483646..2147483646)") + return self.chc_model.intvar(int(val), int(val)) # convert to "variable" + elif isinstance(val, _NumVarImpl): + return self.solver_var(val) # use variable + else: + raise ValueError(f"Cannot convert {val} of type {type(val)} to Choco variable, expected int or NumVarImpl") + + # elif isinstance(val, IntVar): + # return val + # return None + + def _to_vars(self, vals): + if is_any_list(vals): + return [self._to_vars(v) for v in vals] + return self._to_var(vals) + + def transform(self, cpm_expr): + """ + Transform arbitrary CPMpy expressions to constraints the solver supports + + Implemented through chaining multiple solver-independent **transformation functions** from + the `cpmpy/transformations/` directory. + + See the 'Adding a new solver' docs on readthedocs for more information. + + :param cpm_expr: CPMpy expression, or list thereof + :type cpm_expr: Expression or list of Expression + + :return: list of Expression + """ + + cpm_cons = toplevel_list(cpm_expr) + supported = {"min", "max", "abs", "count", "element", "alldifferent", "alldifferent_except0", "allequal", + "table", 'negative_table', "InDomain", "cumulative", "circuit", "gcc", "inverse", "nvalue", "increasing", + "decreasing","strictly_increasing","strictly_decreasing","lex_lesseq", "lex_less", "among", "precedence"} + + # choco supports reification of any constraint, but has a bug in increasing and decreasing + supported_reified = {"min", "max", "abs", "count", "element", "alldifferent", "alldifferent_except0", + "allequal", "table", 'negative_table', "InDomain", "cumulative", "circuit", "gcc", "inverse", "nvalue", + "lex_lesseq", "lex_less", "among"} + + # for when choco new release comes, fixing the bug on increasing and decreasing + #supported_reified = supported + cpm_cons = decompose_in_tree(cpm_cons, supported, supported_reified) + cpm_cons = flatten_constraint(cpm_cons) # flat normal form + cpm_cons = canonical_comparison(cpm_cons) + cpm_cons = reify_rewrite(cpm_cons, supported = supported_reified | {"sum", "wsum"}) # constraints that support reification + cpm_cons = only_numexpr_equality(cpm_cons, supported=frozenset(["sum", "wsum", "sub"])) # support >, <, != + + return cpm_cons + + def __add__(self, cpm_expr): + """ + Eagerly add a constraint to the underlying solver. + + Any CPMpy expression given is immediately transformed (through `transform()`) + and then posted to the solver in this function. + + This can raise 'NotImplementedError' for any constraint not supported after transformation + + The variables used in expressions given to add are stored as 'user variables'. Those are the only ones + the user knows and cares about (and will be populated with a value after solve). All other variables + are auxiliary variables created by transformations. + + :param cpm_expr: CPMpy expression, or list thereof + :type cpm_expr: Expression or list of Expression + + :return: self + """ + # add new user vars to the set + get_variables(cpm_expr, collect=self.user_vars) + # ensure all vars are known to solver + + # transform and post the constraints + for con in self.transform(cpm_expr): + c = self._get_constraint(con) + if c is not None: # Reification constraints are not posted + c.post() + + return self + + def _get_constraint(self, cpm_expr): + """ + Get a solver's constraint by a supported CPMpy constraint + + :param cpm_expr: CPMpy expression + :type cpm_expr: Expression + + """ + + # Operators: base (bool), lhs=numexpr, lhs|rhs=boolexpr (reified ->) + if isinstance(cpm_expr, Operator): + # 'and'/n, 'or'/n, '->'/2 + if cpm_expr.name == 'and': + return self.chc_model.and_(self.solver_vars(cpm_expr.args)) + elif cpm_expr.name == 'or': + return self.chc_model.or_(self.solver_vars(cpm_expr.args)) + + # elif cpm_expr.name == "->": # prepared for if pychoco releases if_then addition + # cond, subexpr = cpm_expr.args + # if isinstance(cond, _BoolVarImpl) and isinstance(subexpr, _BoolVarImpl): + # return self.chc_model.or_(self.solver_vars([~cond, subexpr])) + # elif isinstance(cond, _BoolVarImpl): + # return self._get_constraint(subexpr).implied_by(self.solver_var(cond)) + # elif isinstance(subexpr, _BoolVarImpl): + # return self._get_constraint(cond).implies(self.solver_var(subexpr)) + # else: + # ValueError(f"Unexpected implication: {cpm_expr}") + + elif cpm_expr.name == '->': + cond, subexpr = cpm_expr.args + if isinstance(cond, _BoolVarImpl) and isinstance(subexpr, _BoolVarImpl): # bv -> bv + chc_cond, chc_subexpr = self.solver_vars([cond, subexpr]) + elif isinstance(cond, _BoolVarImpl): # bv -> expr + chc_cond = self.solver_var(cond) + chc_subexpr = self._get_constraint(subexpr).reify() + elif isinstance(subexpr, _BoolVarImpl): # expr -> bv + chc_cond = self._get_constraint(cond).reify() + chc_subexpr = self.solver_var(subexpr) + else: + raise ValueError(f"Unexpected implication {cpm_expr}") + + return self.chc_model.or_([self.chc_model.bool_not_view(chc_cond), chc_subexpr]) + + else: + raise NotImplementedError("Not a known supported Choco Operator '{}' {}".format( + cpm_expr.name, cpm_expr)) + + # Comparisons: both numeric and boolean ones + # numexpr `comp` bvar|const + elif isinstance(cpm_expr, Comparison): + lhs, rhs = cpm_expr.args + op = cpm_expr.name if cpm_expr.name != "==" else "=" + + if is_boolexpr(lhs) and is_boolexpr(rhs): #boolean equality -- Reification + # # prepared for if pychoco releases reify_with addition + # if isinstance(lhs, _BoolVarImpl) and isinstance(lhs, _BoolVarImpl): + # return self.chc_model.all_equal(self.solver_vars([lhs, rhs])) + # elif isinstance(lhs, _BoolVarImpl): + # return self._get_constraint(rhs).reify_with(self.solver_var(lhs)) + # elif isinstance(rhs, _BoolVarImpl): + # return self._get_constraint(lhs).reify_with(self.solver_var(rhs)) + # else: + # raise ValueError(f"Unexpected reification {cpm_expr}") + + if isinstance(lhs, _BoolVarImpl) and isinstance(lhs, _BoolVarImpl): + chc_var, bv = self.solver_vars([lhs, rhs]) + elif isinstance(lhs, _BoolVarImpl): + bv = self._get_constraint(rhs).reify() + chc_var = self.solver_var(lhs) + elif isinstance(rhs, _BoolVarImpl): + bv = self._get_constraint(lhs).reify() + chc_var = self.solver_var(rhs) + else: + raise ValueError(f"Unexpected reification {cpm_expr}") + return self.chc_model.all_equal([chc_var, bv]) + + elif isinstance(lhs, _NumVarImpl): + return self.chc_model.arithm(self.solver_var(lhs), op, self.solver_var(rhs)) + elif isinstance(lhs, Operator) and lhs.name in {'sum','wsum','sub'}: + if lhs.name == 'sum': + return self.chc_model.sum(self.solver_vars(lhs.args), op, self.solver_var(rhs)) + elif lhs.name == "sub": + a, b = self.solver_vars(lhs.args) + return self.chc_model.arithm(a, "-", b, op, self.solver_var(rhs)) + elif lhs.name == 'wsum': + wgt, x = lhs.args + w = np.array(wgt).tolist() + x = self.solver_vars(lhs.args[1]) + return self.chc_model.scalar(x, w, op, self.solver_var(rhs)) + + elif cpm_expr.name == '==': + + chc_rhs = self._to_var(rhs) # result is always var + all_vars = {"min", "max", "abs", "div", "mod", "element", "nvalue"} + if lhs.name in all_vars: + + chc_args = self._to_vars(lhs.args) + + if lhs.name == 'min': # min(vars) = var + return self.chc_model.min(chc_rhs, chc_args) + elif lhs.name == 'max': # max(vars) = var + return self.chc_model.max(chc_rhs, chc_args) + elif lhs.name == 'abs': # abs(var) = var + assert len(chc_args) == 1, f"Expected one argument of abs constraint, but got {chc_args}" + return self.chc_model.absolute(chc_rhs, chc_args[0]) + elif lhs.name == "div": # var / var = var + dividend, divisor = chc_args + return self.chc_model.div(dividend, divisor, chc_rhs) + elif lhs.name == 'mod': # var % var = var + dividend, divisor = chc_args + return self.chc_model.mod(dividend, divisor, chc_rhs) + elif lhs.name == "element": # varsvar[var] = var + # TODO: actually, Choco also supports ints[var] = var, but no mix of var and int in array + arr, idx = chc_args + return self.chc_model.element(chc_rhs, arr, idx) + elif lhs.name == "nvalue": # nvalue(vars) = var + # TODO: should look into leaving nvalue <= arg so can post atmost_nvalues here + return self.chc_model.n_values(chc_args, chc_rhs) + + elif lhs.name == 'count': # count(vars, var/int) = var + arr, val = lhs.args + return self.chc_model.count(self.solver_var(val), self._to_vars(arr), chc_rhs) + elif lhs.name == "among": + arr, vals = lhs.args + return self.chc_model.among(chc_rhs, self._to_vars(arr), vals) + elif lhs.name == 'mul': # var * var/int = var/int + a,b = self.solver_vars(lhs.args) + if isinstance(a, int): + a,b = b,a # int arg should always be second + return self.chc_model.times(a,b, self.solver_var(rhs)) + elif lhs.name == 'pow': # var ^ int = var + chc_rhs = self._to_var(rhs) + return self.chc_model.pow(*self.solver_vars(lhs.args),chc_rhs) + + + + raise NotImplementedError( + "Not a known supported Choco left-hand-side '{}' {}".format(lhs.name, cpm_expr)) + + # base (Boolean) global constraints + elif isinstance(cpm_expr, GlobalConstraint): + + # many globals require all variables as arguments + if cpm_expr.name in {"alldifferent", "alldifferent_except0", "allequal", "circuit", "inverse","increasing","decreasing","strictly_increasing","strictly_decreasing","lex_lesseq","lex_less"}: + chc_args = self._to_vars(cpm_expr.args) + if cpm_expr.name == 'alldifferent': + return self.chc_model.all_different(chc_args) + elif cpm_expr.name == 'alldifferent_except0': + return self.chc_model.all_different_except_0(chc_args) + elif cpm_expr.name == 'allequal': + return self.chc_model.all_equal(chc_args) + elif cpm_expr.name == "circuit": + return self.chc_model.circuit(chc_args) + elif cpm_expr.name == "inverse": + return self.chc_model.inverse_channeling(*chc_args) + elif cpm_expr.name == "increasing": + return self.chc_model.increasing(chc_args,0) + elif cpm_expr.name == "decreasing": + return self.chc_model.decreasing(chc_args,0) + elif cpm_expr.name == "strictly_increasing": + return self.chc_model.increasing(chc_args,1) + elif cpm_expr.name == "strictly_decreasing": + return self.chc_model.decreasing(chc_args,1) + elif cpm_expr.name in ["lex_lesseq", "lex_less"]: + if cpm_expr.name == "lex_lesseq": + return self.chc_model.lex_less_eq(*chc_args) + return self.chc_model.lex_less(*chc_args) +# Ready for when it is fixed in pychoco (https://github.com/chocoteam/pychoco/issues/30) +# elif cpm_expr.name == "lex_chain_less": +# return self.chc_model.lex_chain_less(chc_args) + + # but not all + elif cpm_expr.name == 'table': + assert (len(cpm_expr.args) == 2) # args = [array, table] + array, table = self.solver_vars(cpm_expr.args) + return self.chc_model.table(array, table) + elif cpm_expr.name == 'negative_table': + assert (len(cpm_expr.args) == 2) # args = [array, table] + array, table = self.solver_vars(cpm_expr.args) + return self.chc_model.table(array, table, False) + elif cpm_expr.name == 'InDomain': + assert len(cpm_expr.args) == 2 # args = [array, list of vals] + expr, table = self.solver_vars(cpm_expr.args) + return self.chc_model.member(expr, table) + elif cpm_expr.name == "cumulative": + start, dur, end, demand, cap = cpm_expr.args + # start, end, demand and cap should be var + start, end, demand, cap = self._to_vars([start, end, demand, cap]) + # duration can be var or int + dur = self.solver_vars(dur) + # Create task variables. Choco can create them only one by one + tasks = [self.chc_model.task(s, d, e) for s, d, e in zip(start, dur, end)] + return self.chc_model.cumulative(tasks, demand, cap) + elif cpm_expr.name == "precedence": + return self.chc_model.int_value_precede_chain(self._to_vars(cpm_expr.args[0]), cpm_expr.args[1]) + elif cpm_expr.name == "gcc": + vars, vals, occ = cpm_expr.args + return self.chc_model.global_cardinality(*self.solver_vars([vars, vals]), self._to_vars(occ), cpm_expr.closed) + else: + raise NotImplementedError(f"Unknown global constraint {cpm_expr}, should be decomposed! If you reach this, please report on github.") + + # unlikely base case: Boolean variable + elif isinstance(cpm_expr, _BoolVarImpl): + return self.chc_model.and_([self.solver_var(cpm_expr)]) + + # unlikely base case: True or False + elif isinstance(cpm_expr, BoolVal): + # Choco does not allow to post True or False. Post "certainly True or False" constraints instead + if cpm_expr.args[0] is True: + return None + else: + if self.helper_var is None: + self.helper_var = self.chc_model.intvar(0, 0) + return self.chc_model.arithm(self.helper_var, "<", 0) + + # a direct constraint, pass to solver + elif isinstance(cpm_expr, DirectConstraint): + c = cpm_expr.callSolver(self, self.chc_model) + return c + + # else + raise NotImplementedError(cpm_expr) # if you reach this... please report on github \ No newline at end of file diff --git a/cpmpy/solvers/exact.py b/cpmpy/solvers/exact.py index f19b71377..2573c9fa6 100644 --- a/cpmpy/solvers/exact.py +++ b/cpmpy/solvers/exact.py @@ -6,8 +6,11 @@ """ Interface to Exact - Exact solves decision and optimization problems formulated as integer linear programs. Under the hood, it converts integer variables to binary (0-1) variables and applies highly efficient propagation routines and strong cutting-planes / pseudo-Boolean conflict analysis. + Exact solves decision and optimization problems formulated as integer linear programs. + Under the hood, it converts integer variables to binary (0-1) variables and applies highly efficient + propagation routines and strong cutting-planes / pseudo-Boolean conflict analysis. + The solver's git repository: https://gitlab.com/JoD/exact =============== @@ -18,6 +21,10 @@ :nosignatures: CPM_exact + + ============== + Module details + ============== """ import sys # for stdout checking import time @@ -37,7 +44,7 @@ from ..transformations.normalize import toplevel_list from ..expressions.globalconstraints import DirectConstraint from ..exceptions import NotSupportedError -from ..expressions.utils import flatlist +from ..expressions.utils import flatlist, argvals import numpy as np import numbers @@ -84,7 +91,7 @@ def __init__(self, cpm_model=None, subsolver=None): Arguments: - cpm_model: Model(), a CPMpy Model() (optional) - - subsolver: None + - subsolver: None, not used """ if not self.supported(): raise Exception("Install 'exact' as a Python package to use this solver interface") @@ -155,6 +162,9 @@ def solve(self, time_limit=None, assumptions=None, **kwargs): """ from exact import Exact as xct + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + if not self.solver_is_initialized: assert not self.objective_given # NOTE: initialization of exact is also how it fixes the objective function. @@ -223,6 +233,9 @@ def solveAll(self, display=None, time_limit=None, solution_limit=None, call_from Returns: number of solutions found """ + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + if self.objective_given: raise NotSupportedError("Exact does not support finding all optimal solutions.") @@ -257,9 +270,9 @@ def solveAll(self, display=None, time_limit=None, solution_limit=None, call_from if display is not None: self._fillObjAndVars() if isinstance(display, Expression): - print(display.value()) + print(argval(display)) elif isinstance(display, list): - print([v.value() for v in display]) + print(argvals(display)) else: display() # callback elif my_status == 2: # found inconsistency diff --git a/cpmpy/solvers/gcs.py b/cpmpy/solvers/gcs.py new file mode 100644 index 000000000..c54ccccd7 --- /dev/null +++ b/cpmpy/solvers/gcs.py @@ -0,0 +1,593 @@ +#!/usr/bin/env python +#-*- coding:utf-8 -*- +## +## ortools.py +## +""" + Interface to the Glasgow Constraint Solver's API for the cpmpy library. + The key feature of this solver is the ability to produce proof logs. + + See: + https://github.com/ciaranm/glasgow-constraint-solver + =============== + List of classes + =============== + .. autosummary:: + :nosignatures: + GlasgowConstraintSolver +""" +from cpmpy.transformations.comparison import only_numexpr_equality +from cpmpy.transformations.reification import reify_rewrite, only_bv_reifies, only_implies +from ..exceptions import NotSupportedError +from .solver_interface import SolverInterface, SolverStatus, ExitStatus +from ..expressions.core import Expression, Comparison, Operator, BoolVal +from ..expressions.variables import _BoolVarImpl, _IntVarImpl, _NumVarImpl, NegBoolView, boolvar +from ..expressions.globalconstraints import GlobalConstraint +from ..expressions.utils import is_num, is_any_list, argval, argvals +from ..transformations.decompose_global import decompose_in_tree +from ..transformations.get_variables import get_variables +from ..transformations.flatten_model import flatten_constraint, flatten_objective, get_or_make_var + +from ..transformations.normalize import toplevel_list + +# For proof file handling and verifying +import sys +from os import path +from shutil import which +import subprocess + +class CPM_gcs(SolverInterface): + """ + Interface to Glasgow Constraint Solver's API. + + Requires that the 'gcspy' python package is installed: + Current installation instructions: + - Ensure you have C++20 compiler such as GCC 10.3 / clang 15 + - (on Debian-based systems, see https://apt.llvm.org for easy installation) + - If necessary `export CXX=` + - Ensure you have Boost installed + - `git clone https://github.com/ciaranm/glasgow-constraint-solver.git` + - `cd glasgow-constraint-solver/python` + - `pip install .` + NB: if for any reason you need to retry the build, ensure you remove glasgow-constraints-solver/generator before rebuilding. + + For the verifier functionality, the 'veripb' tool is also required. + See `https://gitlab.com/MIAOresearch/software/VeriPB#installation` for installation instructions. + + Creates the following attributes (see parent constructor for more): + - gcs: the gcspy solver object + - objective_var: optional: the variable used as objective + - proof_location: location of the last proof produced by the solver + - proof_name: name of the last proof (means .opb and .pbp will be present at the proof location) + - proof_verification_failed: whether the last VeriPB check failed to verify the proof. + - proof_check_timeout: whether the last VeriPB check timed out. + """ + + @staticmethod + def supported(): + # try to import the package + try: + import gcspy + return True + except ImportError as e: + return False + + def __init__(self, cpm_model=None, subsolver=None): + """ + Constructor of the native solver object + Arguments: + - cpm_model: Model(), a CPMpy Model() (optional) + - subsolver: None (not supported) + """ + if not self.supported(): + raise Exception("Glasgow Constraint Solver: Install the python package 'gcspy'") + + import gcspy + + assert(subsolver is None) # unless you support subsolvers, see pysat or minizinc + + # initialise the native solver object + self.gcs = gcspy.GCS() + self.objective_var = None + self.proof_location = None + self.proof_name = None + self.proof_check_timeout = True + self.proof_verification_failed = False + + # initialise everything else and post the constraints/objective + super().__init__(name="Glasgow Constraint Solver", cpm_model=cpm_model) + + def has_objective(self): + return self.objective_var is not None + + def solve(self, time_limit=None, prove=False, proof_name=None, proof_location=".", + verify=False, verify_time_limit=None, veripb_args = [], display_verifier_output=True, **kwargs): + """ + Run the Glasgow Constraint Solver, get just one (optimal) solution. + Arguments: + - time_limit: maximum solve time in seconds (float, optional). + - prove: whether to produce a VeriPB proof (.opb model file and .pbp proof file). + - proof_name: name for the the proof files. + - proof_location: location for the proof files (default to current working directory). + - verify: whether to verify the result of the solve run (overrides prove if prove is False) + - verify_time_limit: time limit for verification (ignored if verify=False) + - veripb_args: list of command line arguments to pass to veripb e.g. `--trace --useColor` (run `veripb --help` for a full list) + - display_verifier_output: whether to print the output from VeriPB + - kwargs: currently GCS does not support any additional keyword arguments. + + Returns: whether a solution was found. + """ + # ensure all user vars are known to solver + self.solver_vars(list(self.user_vars)) + + # If we're verifying we must be proving + prove |= verify + # Set default proof name to name of file containing __main__ + if prove and proof_name is None: + if hasattr(sys.modules['__main__'], "__file__"): + self.proof_name = path.splitext(path.basename(sys.modules['__main__'].__file__))[0] + else: + self.proof_name = "gcs_proof" + else: + self.proof_name = proof_name + self.proof_location = proof_location + + # call the solver, with parameters + gcs_stats = self.gcs.solve( + all_solutions=self.has_objective(), + timeout=time_limit, + callback=None, + prove=prove, + proof_name=self.proof_name, + proof_location=proof_location, + **kwargs) + + # new status, translate runtime + self.cpm_status = SolverStatus(self.name) + self.cpm_status.runtime = gcs_stats["solve_time"] + + # translate exit status + if gcs_stats['solutions'] != 0: + self.cpm_status.exitstatus = ExitStatus.FEASIBLE + elif not gcs_stats['completed']: + self.cpm_status.exitstatus = ExitStatus.UNKNOWN + else: + self.cpm_status.exitstatus = ExitStatus.UNSATISFIABLE + + # True/False depending on self.cpm_status + has_sol = self._solve_return(self.cpm_status) + + # translate solution values (of user specified variables only) + self.objective_value_ = None + if has_sol: + # fill in variable values + for cpm_var in self.user_vars: + sol_var = self.solver_var(cpm_var) + if isinstance(cpm_var, _BoolVarImpl): + # Convert back to bool + cpm_var._value = bool(self.gcs.get_solution_value(sol_var)) + else: + cpm_var._value = self.gcs.get_solution_value(sol_var) + + # translate objective, for optimisation problems only + if self.has_objective(): + self.objective_value_ = self.gcs.get_solution_value(self.solver_var(self.objective_var)) + + # Verify proof, if requested + if verify: + self.verify(name=self.proof_name, location=proof_location, time_limit=verify_time_limit, + veripb_args=veripb_args, display_output=display_verifier_output) + return has_sol + + def solveAll(self, time_limit=None, display=None, solution_limit=None, call_from_model=False, + prove=False, proof_name=None, proof_location=".", verify=False, verify_time_limit=None, veripb_args = [], + display_verifier_output=True, **kwargs): + """ + Run the Glasgow Constraint Solver, and get a number of solutions, with optional solution callbacks. + + Arguments: + - display: either a list of CPMpy expressions, OR a callback function, called with the variables after value-mapping + default/None: nothing displayed + - solution_limit: stop after this many solutions (default: None) + - time_limit: maximum solve time in seconds (float, default: None) + - call_from_model: whether the method is called from a CPMpy Model instance or not + - prove: whether to produce a VeriPB proof (.opb model file and .pbp proof file). + - proof_name: name for the the proof files. + - proof_location: location for the proof files (default to current working directory). + - verify: whether to verify the result of the solve run (overrides prove if prove is False) + - verify_time_limit: time limit for verification (ignored if verify=False) + - veripb_args: list of command line arguments to pass to veripb e.g. `--trace --useColor` (run `veripb --help` for a full list) + - display_verifier_output: whether to print the output from VeriPB + - kwargs: currently GCS does not support any additional keyword arguments. + Returns: number of solutions found + """ + if self.has_objective(): + raise NotSupportedError("Glasgow Constraint Solver: does not support finding all optimal solutions.") + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + + # If we're verifying we must be proving + prove |= verify + # Set default proof name to name of file containing __main__ + if prove and proof_name is None: + if hasattr(sys.modules['__main__'], "__file__"): + self.proof_name = path.splitext(path.basename(sys.modules['__main__'].__file__))[0] + else: + self.proof_name = "gcs_proof" + self.proof_location = proof_location + + # Set display callback + def display_callback(solution_map): + for cpm_var in self.user_vars: + sol_var = self.solver_var(cpm_var) + if isinstance(cpm_var, _BoolVarImpl): + # Convert back to bool + cpm_var._value = bool(solution_map[sol_var]) + else: + cpm_var._value = solution_map[sol_var] + + if isinstance(display, Expression): + print(argval(display)) + elif isinstance(display, list): + # explicit list of expressions to display + print(argvals(display)) + elif callable(display): + display() + else: + raise NotImplementedError("Glasgow Constraint Solver: Unknown display type.".format(cpm_var)) + return + sol_callback = None + if display: + sol_callback=display_callback + + gcs_stats = self.gcs.solve( + all_solutions=True, + timeout=time_limit, + solution_limit=solution_limit, + callback=sol_callback, + prove=prove, + proof_name=proof_name, + proof_location=proof_location, **kwargs) + + # new status, get runtime + self.cpm_status = SolverStatus(self.name) + self.cpm_status.runtime = gcs_stats["solve_time"] + + # Verify proof, if requested + if verify: + self.verify(name=self.proof_name, location=proof_location, time_limit=verify_time_limit, + veripb_args=veripb_args, display_output=display_verifier_output) + + return gcs_stats["solutions"] + + def solver_var(self, cpm_var): + """ + Creates solver variable for cpmpy variable + or returns from cache if previously created + """ + if is_num(cpm_var): # shortcut, eases posting constraints + return self.gcs.create_integer_constant(cpm_var) + + # special case, negative-bool-view + # work directly on var inside the view + if isinstance(cpm_var, NegBoolView): + # gcs only works with integer variables, so not(x) = -x + 1 + return self.gcs.add_constant(self.gcs.negate(self.solver_var(cpm_var._bv)), 1) + + # create if it does not exist + if cpm_var not in self._varmap: + if isinstance(cpm_var, _BoolVarImpl): + # Bool vars are just int vars with [0, 1] domain + revar = self.gcs.create_integer_variable(0, 1, str(cpm_var)) + elif isinstance(cpm_var, _IntVarImpl): + revar = self.gcs.create_integer_variable(cpm_var.lb, cpm_var.ub, str(cpm_var)) + else: + raise NotImplementedError("Not a known var {}".format(cpm_var)) + self._varmap[cpm_var] = revar + + return self._varmap[cpm_var] + + def objective(self, expr, minimize=True): + """ + Post the given expression to the solver as objective to minimize/maximize + 'objective()' can be called multiple times, only the last one is stored + (technical side note: any constraints created during conversion of the objective + are permanently posted to the solver) + """ + # make objective function non-nested + (flat_obj, flat_cons) = flatten_objective(expr) + self += flat_cons # 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) + self += obj_cons + + self.objective_var = obj + + if minimize: + self.gcs.minimise(self.solver_var(obj)) + else: + self.gcs.maximise(self.solver_var(obj)) + + def transform(self, cpm_expr): + """ + Transform arbitrary CPMpy expressions to constraints the solver supports + + Implemented through chaining multiple solver-independent **transformation functions** from + the `cpmpy/transformations/` directory. + + See the 'Adding a new solver' docs on readthedocs for more information. + + :param cpm_expr: CPMpy expression, or list thereof + :type cpm_expr: Expression or list of Expression + + :return: list of Expression + """ + cpm_cons = toplevel_list(cpm_expr) + supported = { + "min", + "max", + "abs", + "alldifferent", + "element", + 'table', + 'negative_table', + 'count', + 'nvalue', + 'inverse', + 'circuit', + 'xor'} + cpm_cons = decompose_in_tree(cpm_cons, supported) + cpm_cons = flatten_constraint(cpm_cons) # flat normal form + + # NB: GCS supports full reification for linear equality and linear inequaltiy constraints + # but no reification for linear not equals and not half reification for linear equality. + # Maybe a future transformation (or future work on the GCS solver). + cpm_cons = reify_rewrite(cpm_cons, supported=frozenset(['=='])) + cpm_cons = only_numexpr_equality(cpm_cons, supported=frozenset(["sum", "wsum"])) # supports >, <, != + + # NB: GCS supports a small number of simple expressions as the reifying term + # e.g. (x > 3) -> constraint could in principle be supported in future. + cpm_cons = only_bv_reifies(cpm_cons) + str_rep = "" + for c in cpm_cons: + str_rep += str(c) + '\n' + return cpm_cons + + def verify(self, name=None, location=".", time_limit=None, display_output=False, veripb_args=[]): + """ + Verify a solver-produced proof using VeriPB. + + Arguments: + - name: name for the the proof files (default to self.proof_name) + - location: location for the proof files (default to current working directory). + - time_limit: time limit for verification (ignored if verify=False) + - veripb_args: list of command line arguments to pass to veripb e.g. `--trace --useColor` (run `veripb --help` for a full list) + - display_output: whether to print the output from VeriPB + """ + if not which("veripb"): + raise Exception("Unable to run VeriPB: check it is installed and on system path") + + if name is None: + name = self.proof_name + if name is None: # Still None? + raise ValueError("No proof to verify") + + try: + result = subprocess.run(["veripb", path.join(location, name +".opb"), path.join(location, name +".pbp")] + veripb_args, + capture_output=True, text=True, timeout=time_limit) + self.proof_check_timeout = False + self.proof_verification_failed = result.returncode + except subprocess.TimeoutExpired: + self.proof_check_timeout = True + self.proof_verification_failed = False + + if display_output: + print(result.stdout) + print(result.stderr) + + return self.proof_verification_failed + + def __add__(self, cpm_cons): + """ + Post a (list of) CPMpy constraints(=expressions) to the solver + Note that we don't store the constraints in a cpm_model, + we first transform the constraints into primitive constraints, + then post those primitive constraints directly to the native solver + :param cpm_con CPMpy constraint, or list thereof + :type cpm_con (list of) Expression(s) + """ + # add new user vars to the set + # add new user vars to the set + get_variables(cpm_cons, collect=self.user_vars) + + for con in self.transform(cpm_cons): + cpm_expr = con + if isinstance(cpm_expr, _BoolVarImpl): + # base case, just var or ~var + self.gcs.post_or([self.solver_var(cpm_expr)]) + elif isinstance(cpm_expr, BoolVal): + if not cpm_expr: + self.gcs.post_or([]) + else: + return + elif isinstance(cpm_expr, Operator) or \ + (cpm_expr.name == '==' and isinstance(cpm_expr.args[0], _BoolVarImpl) \ + and not isinstance(cpm_expr.args[1], _NumVarImpl)): + # ^ Somewhat awkward, but want to deal with full and half reifications + # in one go here, and then deal with regular == comparisons later.l + + # 'and'/n, 'or'/n, '->'/2 + if cpm_expr.name == 'and': + self.gcs.post_and(self.solver_vars(cpm_expr.args)) + elif cpm_expr.name == 'or': + self.gcs.post_or(self.solver_vars(cpm_expr.args)) + + # Reified constraint: BoolVar -> Boolexpr or BoolVar == Boolexpr + # LHS must be boolvar due to only_bv_reifies + elif cpm_expr.name == '->' or cpm_expr.name == '==': + fully_reify = (cpm_expr.name == '==') + assert(isinstance(cpm_expr.args[0], _BoolVarImpl)) + bool_lhs = cpm_expr.args[0] + reif_var = self.solver_var(bool_lhs) + bool_expr = cpm_expr.args[1] + + # Just a plain implies or equals between two boolvars + if isinstance(bool_expr, _BoolVarImpl): # bv1 -> bv2 + (self.gcs.post_equals if fully_reify else self.gcs.post_implies)\ + (*self.solver_vars([bool_lhs, bool_expr])) + elif isinstance(bool_expr, Operator): # bv -> and(...), bv -> or(...) bv -> [->] + if bool_expr.name == 'and': + self.gcs.post_and_reif(self.solver_vars(bool_expr.args), reif_var, fully_reify) + elif bool_expr.name == 'or': + self.gcs.post_or_reif(self.solver_vars(bool_expr.args), reif_var, fully_reify) + elif bool_expr.name == '->': + self.gcs.post_implies_reif(self.solver_vars(bool_expr.args), reif_var, fully_reify) + else: + # Shouldn't happen if reify_rewrite worked? + raise NotImplementedError("Not currently supported by Glasgow Constraint Solver API '{}' {}".format) + + # Reified Comparison + elif isinstance(bool_expr, Comparison): + lhs = bool_expr.args[0] + rhs = bool_expr.args[1] + if bool_expr.name == '==': + self.gcs.post_equals_reif(*self.solver_vars([lhs, rhs]), reif_var, fully_reify) + elif bool_expr.name == '<=': + self.gcs.post_less_than_equal_reif(*self.solver_vars([lhs, rhs]), reif_var, fully_reify) + elif bool_expr.name == '<': + self.gcs.post_less_than_reif(*self.solver_vars([lhs, rhs]), reif_var, fully_reify) + elif bool_expr.name == '>=': + self.gcs.post_greater_than_equal_reif(*self.solver_vars([lhs, rhs]), reif_var, fully_reify) + elif bool_expr.name == '>': + self.gcs.post_greater_than_reif(*self.solver_vars([lhs, rhs]), reif_var, fully_reify) + elif bool_expr.name == '!=': + # Note: GCS doesn't currently support reified NotEquals, so we need this ugly workaround for now: + # bv -> x != y can be written as + # bv -> OR(lt, gt) with lt, gt being BoolVars and the additional constraints + # lt == x < y + # gt == x > y + lt_bool, gt_bool = boolvar(shape=2) + self += (lhs < rhs) == lt_bool + self += (lhs > rhs) == gt_bool + if fully_reify: + self += (~bool_lhs).implies(lhs == rhs) + 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) + else: + # Shouldn't happen if reify_rewrite worked + raise NotImplementedError("Not currently supported by Glasgow Constraint Solver API '{}' {}".format) + + # Normal comparison + elif isinstance(cpm_expr, Comparison): + lhs = cpm_expr.args[0] + rhs = cpm_expr.args[1] + + # Due to only_numexpr_equality we can only have '!=', "<=", etc. + # when the lhs is a variable, sum or wsum + if isinstance(lhs, _NumVarImpl) or lhs.name == 'sum' or lhs.name == 'wsum': + if lhs.name == 'sum' or lhs.name == 'wsum': + if lhs.name == 'sum': + summands = self.solver_vars(lhs.args) + summands.append(self.solver_var(rhs)) + coeffs = [1]*len(lhs.args) + [-1] + else: + summands = self.solver_vars(lhs.args[1]) + summands.append(self.solver_var(rhs)) + coeffs = list(lhs.args[0]) + [-1] + + if cpm_expr.name == '==': + self.gcs.post_linear_equality(summands, coeffs, 0) + elif cpm_expr.name == '!=': + self.gcs.post_linear_not_equal(summands, coeffs, 0) + elif cpm_expr.name == '<=': + self.gcs.post_linear_less_equal(summands, coeffs, 0) + elif cpm_expr.name == '<': + self.gcs.post_linear_less_equal(summands, coeffs, -1) + elif cpm_expr.name == '>=': + self.gcs.post_linear_greater_equal(summands, coeffs, 0) + elif cpm_expr.name == '>': + self.gcs.post_linear_greater_equal(summands, coeffs, 1) + else: + raise NotImplementedError("Not currently supported by Glasgow Constraint Solver API '{}'".format(cpm_expr)) + else: + if cpm_expr.name == '==': + self.gcs.post_equals(*self.solver_vars([lhs, rhs])) + elif cpm_expr.name == '!=': + self.gcs.post_not_equals(*self.solver_vars([lhs, rhs])) + elif cpm_expr.name == '<=': + self.gcs.post_less_than_equal(*self.solver_vars([lhs, rhs])) + elif cpm_expr.name == '<': + self.gcs.post_less_than(*self.solver_vars([lhs, rhs])) + elif cpm_expr.name == '>=': + self.gcs.post_greater_than_equal(*self.solver_vars([lhs, rhs])) + elif cpm_expr.name == '>': + self.gcs.post_greater_than(*self.solver_vars([lhs, rhs])) + else: + raise NotImplementedError("Not currently supported by Glasgow Constraint Solver API '{}'".format(cpm_expr)) + + # If the comparison is '==' we can have a NumExpr on the lhs + elif cpm_expr.name == '==': + if lhs.name == 'abs': + assert(len(lhs.args) == 1) # Should not have a nested expression inside abs + self.gcs.post_abs(*self.solver_vars(list(lhs.args) + [rhs])) + elif lhs.name in ['mul', 'div', 'pow', 'mod']: + self.gcs.post_arithmetic(*self.solver_vars(list(lhs.args) + [rhs]), lhs.name) + elif lhs.name == 'sub': + var1 = self.solver_var(lhs.args[0]) + nVar2 = self.gcs.negate(self.solver_var(lhs.args[1])) + self.gcs.post_arithmetic(var1, nVar2, self.solver_var(rhs), 'sum') + elif lhs.name == 'sum' and len(lhs.args) == 2: + var1 = self.solver_var(lhs.args[0]) + var2 = self.solver_var(lhs.args[1]) + self.gcs.post_arithmetic(var1, var2, self.solver_var(rhs), 'sum') + elif lhs.name == 'sum' and len(lhs.args) > 2: + summands = self.solver_vars(lhs.args) + summands.append(self.gcs.negate(self.solver_var(rhs))) + self.gcs.post_linear_equality(summands, [1]*len(summands), 0) + elif lhs.name == 'wsum': + summands = self.solver_vars(lhs.args[1]) + summands.append(self.gcs.negate(self.solver_var(rhs))) + self.gcs.post_linear_equality(summands, list(lhs.args[0]) + [1], 0) + elif lhs.name == 'max': + self.gcs.post_max(self.solver_vars(lhs.args), self.solver_var(rhs)) + elif lhs.name == 'min': + self.gcs.post_min(self.solver_vars(lhs.args), self.solver_var(rhs)) + elif lhs.name == 'element': + self.gcs.post_element(self.solver_var(rhs), self.solver_vars(lhs.args[1]), self.solver_vars(lhs.args[0])) + elif lhs.name == 'count': + self.gcs.post_count(self.solver_vars(lhs.args[0]), self.solver_var(lhs.args[1]), self.solver_var(rhs)) + elif lhs.name == 'nvalue': + self.gcs.post_nvalue(self.solver_var(rhs), self.solver_vars(lhs.args)) + else: + # Think that's all the possible NumExprs? + raise NotImplementedError("Not currently supported by Glasgow Constraint Solver API '{}'".format(cpm_expr)) + else: + raise NotImplementedError("Not currently supported by Glasgow Constraint Solver API '{}'".format(cpm_expr)) + + # rest: base (Boolean) global constraints + elif cpm_expr.name == 'xor': + self.gcs.post_xor(self.solver_vars(cpm_expr.args)) + elif cpm_expr.name == 'circuit': + self.gcs.post_circuit(self.solver_vars(cpm_expr.args)) + elif cpm_expr.name == 'inverse': + self.gcs.post_inverse(self.solver_vars(cpm_expr.args[0]), self.solver_vars(cpm_expr.args[1])) + elif cpm_expr.name == 'alldifferent': + self.gcs.post_alldifferent(self.solver_vars(cpm_expr.args)) + elif cpm_expr.name == 'table': + self.gcs.post_table(self.solver_vars(cpm_expr.args[0]), cpm_expr.args[1]) + elif cpm_expr.name == 'negative_table': + self.gcs.post_negative_table(self.solver_vars(cpm_expr.args[0]), cpm_expr.args[1]) + 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... + else: + # Hopefully we don't end up here. + raise NotImplementedError(cpm_expr) + + return self + + + diff --git a/cpmpy/solvers/gurobi.py b/cpmpy/solvers/gurobi.py index ccf557769..de6a58499 100644 --- a/cpmpy/solvers/gurobi.py +++ b/cpmpy/solvers/gurobi.py @@ -1,4 +1,8 @@ #!/usr/bin/env python +#-*- coding:utf-8 -*- +## +## gurobi.py +## """ Interface to the python 'gurobi' package @@ -6,12 +10,13 @@ $ pip install gurobipy - as well as the Gurobi bundled binary packages, downloadable from: - https://www.gurobi.com/ In contrast to other solvers in this package, Gurobi is not free to use and requires an active licence You can read more about available licences at https://www.gurobi.com/downloads/ + Documentation of the solver's own Python API: + https://www.gurobi.com/documentation/current/refman/py_python_api_details.html + =============== List of classes =============== @@ -28,6 +33,7 @@ from .solver_interface import SolverInterface, SolverStatus, ExitStatus from ..expressions.core import * +from ..expressions.utils import argvals from ..expressions.variables import _BoolVarImpl, NegBoolView, _IntVarImpl, _NumVarImpl, intvar from ..expressions.globalconstraints import DirectConstraint from ..transformations.comparison import only_numexpr_equality @@ -82,6 +88,7 @@ def __init__(self, cpm_model=None, subsolver=None): Arguments: - cpm_model: a CPMpy Model() + - subsolver: None, not used """ if not self.supported(): raise Exception( @@ -115,6 +122,9 @@ def solve(self, time_limit=None, solution_callback=None, **kwargs): """ from gurobipy import GRB + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + if time_limit is not None: self.grb_model.setParam("TimeLimit", time_limit) @@ -458,9 +468,9 @@ def solveAll(self, display=None, time_limit=None, solution_limit=None, call_from if display is not None: if isinstance(display, Expression): - print(display.value()) + print(argval(display)) elif isinstance(display, list): - print([v.value() for v in display]) + print(argvals(display)) else: display() # callback diff --git a/cpmpy/solvers/minizinc.py b/cpmpy/solvers/minizinc.py index e601a7268..5026c7b1f 100644 --- a/cpmpy/solvers/minizinc.py +++ b/cpmpy/solvers/minizinc.py @@ -1,8 +1,18 @@ #!/usr/bin/env python +#-*- coding:utf-8 -*- +## +## minizinc.py +## """ Interface to MiniZinc's Python API - CPMpy can translate CPMpy models to the (text-based) MiniZinc language. + Requires that the 'minizinc' python package is installed: + + $ pip install minizinc + + as well as the Minizinc bundled binary packages, downloadable from: + https://github.com/MiniZinc/MiniZincIDE/releases + MiniZinc is a free and open-source constraint modeling language. MiniZinc is used to model constraint satisfaction and optimization problems in @@ -14,6 +24,8 @@ Documentation of the solver's own Python API: https://minizinc-python.readthedocs.io/ + CPMpy can translate CPMpy models to the (text-based) MiniZinc language. + =============== List of classes =============== @@ -22,25 +34,29 @@ :nosignatures: CPM_minizinc + + ============== + Module details + ============== """ import re import warnings import sys import os -from datetime import timedelta # for mzn's timeout +from datetime import timedelta # for mzn's timeout import numpy as np from .solver_interface import SolverInterface, SolverStatus, ExitStatus from ..exceptions import MinizincNameException, MinizincBoundsException from ..expressions.core import Expression, Comparison, Operator, BoolVal -from ..expressions.variables import _NumVarImpl, _IntVarImpl, _BoolVarImpl, NegBoolView +from ..expressions.variables import _NumVarImpl, _IntVarImpl, _BoolVarImpl, NegBoolView, cpm_array from ..expressions.globalconstraints import DirectConstraint -from ..expressions.utils import is_num, is_any_list +from ..expressions.utils import is_num, is_any_list, argvals, argval from ..transformations.decompose_global import decompose_in_tree -from ..transformations.get_variables import get_variables from ..exceptions import MinizincPathException, NotSupportedError -from ..transformations.normalize import toplevel_list, simplify_boolean +from ..transformations.get_variables import get_variables +from ..transformations.normalize import toplevel_list class CPM_minizinc(SolverInterface): @@ -71,15 +87,40 @@ class CPM_minizinc(SolverInterface): The `DirectConstraint`, when used, adds a constraint with that name and the given args to the MiniZinc model. """ + required_version = (2, 8, 2) + @staticmethod def supported(): + return CPM_minizinc.installed() and CPM_minizinc.executable_installed() and not CPM_minizinc.outdated() + + @staticmethod + def installed(): # try to import the package try: + # check if MiniZinc Python is installed import minizinc return True except ImportError as e: return False + @staticmethod + def executable_installed(): + # check if MiniZinc executable is installed + from minizinc import default_driver + if default_driver is None: + warnings.warn("MiniZinc Python is installed, but the MiniZinc executable is missing in path.") + return False + return True + + @staticmethod + def outdated(): + from minizinc import default_driver + if default_driver.parsed_version >= CPM_minizinc.required_version: + return False + else: + # outdated + return True + @staticmethod def solvernames(): """ @@ -107,6 +148,7 @@ def solvernames(): 'symdiff', 'test', 'then', 'true', 'tuple', 'type', 'union', 'var', 'where', 'xor']) # variable names must have this pattern mzn_name_pattern = re.compile('^[A-Za-z][A-Za-z0-9_]*$') + def __init__(self, cpm_model=None, subsolver=None): """ Constructor of the native solver object @@ -116,8 +158,15 @@ def __init__(self, cpm_model=None, subsolver=None): - subsolver: str, name of a subsolver (optional) has to be one of solvernames() """ - if not self.supported(): + if not self.installed(): raise Exception("CPM_minizinc: Install the python package 'minizinc'") + elif not self.executable_installed(): + raise Exception("CPM_minizinc: Install the MiniZinc executable and make it available in path.") + elif self.outdated(): + version = str(self.required_version[0]) + for x in self.required_version[1:]: + version = version + "." + str(x) + raise ImportError("Your Minizinc compiler is outdated, please upgrade to a version >= " + version) import minizinc @@ -126,7 +175,7 @@ def __init__(self, cpm_model=None, subsolver=None): # default solver subsolver = "gecode" elif subsolver.startswith('minizinc:'): - subsolver = subsolver[9:] # strip 'minizinc:' + subsolver = subsolver[9:] # strip 'minizinc:' # initialise the native solver object # (so its params can still be changed before calling solve) @@ -138,11 +187,9 @@ def __init__(self, cpm_model=None, subsolver=None): self.mzn_txt_solve = "solve satisfy;" self.mzn_result = None - # initialise everything else and post the constraints/objective super().__init__(name="minizinc:"+subsolver, cpm_model=cpm_model) - def _pre_solve(self, time_limit=None, **kwargs): """ shared by solve() and solveAll() """ import minizinc @@ -152,12 +199,12 @@ def _pre_solve(self, time_limit=None, **kwargs): # hack, we need to add the objective in a way that it can be changed # later, so make copy of the mzn_model - copy_model = self.mzn_model.__copy__() # it is implemented + copy_model = self.mzn_model.__copy__() # it is implemented copy_model.add_string(self.mzn_txt_solve) # Transform Model into an instance mzn_inst = minizinc.Instance(self.mzn_solver, copy_model) - kwargs['output-time'] = True # required for time getting + kwargs['output-time'] = True # required for time getting return (kwargs, mzn_inst) def solve(self, time_limit=None, **kwargs): @@ -179,6 +226,10 @@ def solve(self, time_limit=None, **kwargs): Does not store the minizinc.Instance() or minizinc.Result() """ + + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + # make mzn_inst (mzn_kwargs, mzn_inst) = self._pre_solve(time_limit=time_limit, **kwargs) @@ -187,7 +238,7 @@ def solve(self, time_limit=None, **kwargs): try: self.mzn_result = mzn_inst.solve(**mzn_kwargs) except minizinc.error.MiniZincError as e: - if sys.platform == "win32" or sys.platform == "cygwin": #path error can occur in windows + if sys.platform == "win32" or sys.platform == "cygwin": # path error can occur in windows path = os.environ.get("path") if "MiniZinc" in str(path): warnings.warn('You might have the wrong minizinc PATH set (windows user Environment Variables') @@ -204,7 +255,7 @@ def solve(self, time_limit=None, **kwargs): # translate solution values (of user specified variables only) self.objective_value_ = None - if has_sol: #mzn_result.status.has_solution(): + if has_sol: # mzn_result.status.has_solution(): mznsol = self.mzn_result.solution if is_any_list(mznsol): print("Warning: multiple solutions found, only returning last one") @@ -216,7 +267,7 @@ def solve(self, time_limit=None, **kwargs): if hasattr(mznsol, sol_var): cpm_var._value = getattr(mznsol, sol_var) else: - print("Warning, no value for ", sol_var) + raise ValueError(f"Var {cpm_var} is unknown to the Minizinc solver, this is unexpected - please report on github...") # translate objective, for optimisation problems only (otherwise None) self.objective_value_ = self.mzn_result.objective @@ -239,7 +290,7 @@ def _post_solve(self, mzn_result): if runtime != 0: self.cpm_status.runtime = runtime else: - raise NotImplementedError #Please report on github, minizinc probably changed their time names/types + raise NotImplementedError # Please report on github, minizinc probably changed their time names/types # translate exit status mzn_status = mzn_result.status @@ -268,7 +319,7 @@ def mzn_time_to_seconds(self, time): elif isinstance(time, timedelta): return time.total_seconds() # --output-time else: - raise NotImplementedError #unexpected type for time + raise NotImplementedError # unexpected type for time async def _solveAll(self, display=None, time_limit=None, solution_limit=None, **kwargs): """ Special 'async' function because mzn.solutions() is async """ @@ -297,11 +348,11 @@ async def _solveAll(self, display=None, time_limit=None, solution_limit=None, ** # and the actual displaying if isinstance(display, Expression): - print(display.value()) + print(argval(display)) elif isinstance(display, list): - print([v.value() for v in display]) + print(argvals(display)) else: - display() # callback + display() # callback # count and stop solution_count += 1 @@ -316,7 +367,6 @@ async def _solveAll(self, display=None, time_limit=None, solution_limit=None, ** return solution_count - def solver_var(self, cpm_var) -> str: """ Creates solver variable for cpmpy variable @@ -334,13 +384,11 @@ def solver_var(self, cpm_var) -> str: return str(cpm_var) if cpm_var not in self._varmap: - # we assume all variables are user variables (because no transforms) - self.user_vars.add(cpm_var) # clean the varname varname = cpm_var.name mzn_var = varname.replace(',', '_').replace('.', '_').replace(' ', '_').replace('[', '_').replace(']', '') - #test if the name is a valid minizinc identifier + # test if the name is a valid minizinc identifier if not self.mzn_name_pattern.search(mzn_var): raise MinizincNameException("Minizinc only accept names with alphabetic characters, digits and underscores. " "First character must be an alphabetic character") @@ -357,7 +405,6 @@ def solver_var(self, cpm_var) -> str: return self._varmap[cpm_var] - def objective(self, expr, minimize): """ Post the given expression to the solver as objective to minimize/maximize @@ -367,7 +414,7 @@ def objective(self, expr, minimize): 'objective()' can be called multiple times, only the last one is stored """ - #get_variables(expr, collect=self.user_vars) # add objvars to vars # all are user vars + # get_variables(expr, collect=self.user_vars) # add objvars to vars # all are user vars # make objective function or variable and post obj = self._convert_expression(expr) @@ -391,10 +438,12 @@ def transform(self, cpm_expr): :return: list of Expression """ cpm_cons = toplevel_list(cpm_expr) - supported = {"min", "max", "abs", "element", "count", "alldifferent", "alldifferent_except0", "allequal", - "inverse", "ite" "xor", "table", "cumulative", "circuit", "gcc"} - return decompose_in_tree(cpm_cons, supported) - + supported = {"min", "max", "abs", "element", "count", "nvalue", "alldifferent", "alldifferent_except0", "allequal", + "inverse", "ite" "xor", "table", "cumulative", "circuit", "gcc", "increasing", "decreasing", + "precedence", "no_overlap", + "strictly_increasing", "strictly_decreasing", "lex_lesseq", "lex_less", "lex_chain_less", + "lex_chain_lesseq", "among"} + return decompose_in_tree(cpm_cons, supported, supported_reified=supported - {"circuit", "precedence"}) def __add__(self, cpm_expr): """ @@ -414,8 +463,7 @@ def __add__(self, cpm_expr): :return: self """ - # all variables are user variables, handled in `solver_var()` - + get_variables(cpm_expr, collect=self.user_vars) # transform and post the constraints for cpm_con in self.transform(cpm_expr): # Get text expression, add to the solver @@ -432,23 +480,18 @@ def _convert_expression(self, expr) -> str: function that returns strings """ if is_any_list(expr): - if len(expr) == 1: - # unary special case, don't put in list - # continue with later code - expr = expr[0] + if isinstance(expr, np.ndarray): + # must flatten + expr_str = [self._convert_expression(e) for e in expr.flat] else: - if isinstance(expr, np.ndarray): - # must flatten - expr_str = [self._convert_expression(e) for e in expr.flat] - else: - expr_str = [self._convert_expression(e) for e in expr] - return "[{}]".format(",".join(expr_str)) + expr_str = [self._convert_expression(e) for e in expr] + return "[{}]".format(",".join(expr_str)) - if isinstance(expr,(bool,np.bool_)): + if isinstance(expr, (bool, np.bool_)): expr = BoolVal(expr) if not isinstance(expr, Expression): - return self.solver_var(expr) # constants + return self.solver_var(expr) # constants if isinstance(expr, BoolVal): return str(expr.args[0]).lower() @@ -479,29 +522,27 @@ def zero_based(array): if expr.name == "alldifferent_except0": args_str = [self._convert_expression(e) for e in expr.args] - return "alldifferent_except_0({})".format(args_str) - - # count: we need the lhs and rhs together - if isinstance(expr, Comparison) and expr.args[0].name == 'count': - name = expr.name - lhs, rhs = expr.args - c = self._convert_expression(rhs) # count - x = [self._convert_expression(countable) for countable in lhs.args[0]] # array - y = self._convert_expression(lhs.args[1]) # value to count in array - functionmap = {'==': 'count_eq', '!=': 'count_neq', - '<=': 'count_geq', '>=': 'count_leq', - '>': 'count_lt', '<': 'count_gt'} - if name in functionmap: - name = functionmap[name] - return "{}({},{},{})".format(name, x, y, c) + return "alldifferent_except_0([{}])".format(",".join(args_str)) - args_str = [self._convert_expression(e) for e in expr.args] + if expr.name in ["lex_lesseq", "lex_less"]: + X = [self._convert_expression(e) for e in expr.args[0]] + Y = [self._convert_expression(e) for e in expr.args[1]] + return f"{expr.name}({{}}, {{}})".format(X, Y) + + if expr.name in ["lex_chain_less", "lex_chain_lesseq"]: + X = cpm_array([[self._convert_expression(e) for e in row] for row in expr.args]) + str_X = "[|\n" # opening + for row in X.T: # Minizinc enforces lexicographic order on columns + str_X += ",".join(map(str, row)) + " |" # rows + str_X += "\n|]" # closing + return f"{expr.name}({{}})".format(str_X) + args_str = [self._convert_expression(e) for e in expr.args] # standard expressions: comparison, operator, element if isinstance(expr, Comparison): # wrap args that are a subexpression in () for i, arg_str in enumerate(args_str): - if isinstance(expr.args[i], Expression): #(Comparison, Operator) + if isinstance(expr.args[i], Expression): # (Comparison, Operator) args_str[i] = "(" + args_str[i] + ")" # infix notation return "{} {} {}".format(args_str[0], expr.name, args_str[1]) @@ -529,7 +570,7 @@ def zero_based(array): # I don't think there is a more direct way unfortunately w = [self._convert_expression(wi) for wi in expr.args[0]] x = [self._convert_expression(xi) for xi in expr.args[1]] - args_str = [f"{wi}*({xi})" for wi,xi in zip(w,x)] + args_str = [f"{wi}*({xi})" for wi, xi in zip(w, x)] return "{}([{}])".format("sum", ",".join(args_str)) # special case, infix: two args @@ -571,15 +612,21 @@ def zero_based(array): elif expr.name == "cumulative": start, dur, end, _, _ = expr.args - self += [s + d == e for s,d,e in zip(start,dur,end)] - if len(start) == 1: - assert len(start) == 1 - format_str = "cumulative([{}],[{}],[{}],{})" - else: - format_str = "cumulative({},{},{},{})" + + durstr = self._convert_expression([s + d == e for s, d, e in zip(start, dur, end)]) + format_str = "forall(" + durstr + " ++ [cumulative({},{},{},{})])" return format_str.format(args_str[0], args_str[1], args_str[3], args_str[4]) + elif expr.name == "precedence": + return "value_precede_chain({},{})".format(args_str[1], args_str[0]) + + elif expr.name == "no_overlap": + start, dur, end = expr.args + durstr = self._convert_expression([s + d == e for s, d, e in zip(start, dur, end)]) + format_str = "forall(" + durstr + " ++ [disjunctive({},{})])" + return format_str.format(args_str[0], args_str[1]) + elif expr.name == 'ite': cond, tr, fal = expr.args return "if {} then {} else {} endif".format(self._convert_expression(cond), self._convert_expression(tr), @@ -590,17 +637,33 @@ def zero_based(array): vars = self._convert_expression(vars) vals = self._convert_expression(vals) occ = self._convert_expression(occ) - return "global_cardinality({},{},{})".format(vars,vals,occ) + if expr.closed is False: + name = "global_cardinality" + else: + name = "global_cardinality_closed" + return "{}({},{},{})".format(name, vars, vals, occ) elif expr.name == "abs": return "abs({})".format(args_str[0]) + elif expr.name == "count": + vars, val = expr.args + vars = self._convert_expression(vars) + val = self._convert_expression(val) + return "count({},{})".format(vars, val) + + elif expr.name == "among": + vars, vals = expr.args + vars = self._convert_expression(vars) + vals = self._convert_expression(vals).replace("[", "{").replace("]", "}") # convert to set + return "among({},{})".format(vars, vals) + # a direct constraint, treat differently for MiniZinc, a text-based language # use the name as, unpack the arguments from the argument tuple elif isinstance(expr, DirectConstraint): return "{}({})".format(expr.name, ",".join(args_str)) - print_map = {"allequal":"all_equal", "xor":"xorall"} + print_map = {"allequal": "all_equal", "xor": "xorall"} if expr.name in print_map: return "{}([{}])".format(print_map[expr.name], ",".join(args_str)) @@ -619,7 +682,7 @@ def solveAll(self, display=None, time_limit=None, solution_limit=None, call_from - time_limit: stop after this many seconds (default: None) - solution_limit: stop after this many solutions (default: None) - call_from_model: whether the method is called from a CPMpy Model instance or not - - any other keyword argument + - kwargs: any keyword argument, sets parameters of solver object, overwrites construction-time kwargs Returns: number of solutions found """ diff --git a/cpmpy/solvers/ortools.py b/cpmpy/solvers/ortools.py index 1e47fb14b..cb58e6598 100644 --- a/cpmpy/solvers/ortools.py +++ b/cpmpy/solvers/ortools.py @@ -4,7 +4,12 @@ ## ortools.py ## """ - Interface to ortools' CP-SAT Python API + Interface to OR-Tools' CP-SAT Python API. + + The 'ortools' python package is bundled by default with CPMpy. + It can be installed through `pip`: + + $ pip install ortools Google OR-Tools is open source software for combinatorial optimization, which seeks to find the best solution to a problem out of a very large set of possible solutions. @@ -22,6 +27,10 @@ :nosignatures: CPM_ortools + + ============== + Module details + ============== """ import sys # for stdout checking import numpy as np @@ -32,7 +41,7 @@ from ..expressions.globalconstraints import DirectConstraint from ..expressions.variables import _NumVarImpl, _IntVarImpl, _BoolVarImpl, NegBoolView, boolvar from ..expressions.globalconstraints import GlobalConstraint -from ..expressions.utils import is_num, is_any_list, eval_comparison, flatlist +from ..expressions.utils import is_num, is_any_list, eval_comparison, flatlist, argval, argvals from ..transformations.decompose_global import decompose_in_tree from ..transformations.get_variables import get_variables from ..transformations.flatten_model import flatten_constraint, flatten_objective @@ -42,7 +51,7 @@ class CPM_ortools(SolverInterface): """ - Interface to the python 'ortools' CP-SAT API + Interface to the Python 'ortools' CP-SAT API Requires that the 'ortools' python package is installed: $ pip install ortools @@ -79,7 +88,7 @@ def __init__(self, cpm_model=None, subsolver=None): Arguments: - cpm_model: Model(), a CPMpy Model() (optional) - - subsolver: None + - subsolver: None, not used """ if not self.supported(): raise Exception("Install the python 'ortools' package to use this solver interface") @@ -133,6 +142,8 @@ def solve(self, time_limit=None, assumptions=None, solution_callback=None, **kwa """ from ortools.sat.python import cp_model as ort + # ensure all user vars are known to solver + self.solver_vars(list(self.user_vars)) # set time limit? if time_limit is not None: @@ -196,7 +207,7 @@ def solve(self, time_limit=None, assumptions=None, solution_callback=None, **kwa if isinstance(cpm_var, _BoolVarImpl): cpm_var._value = bool(cpm_var._value) # ort value is always an int except IndexError: - cpm_var._value = None # probably got optimized away by our transformations + raise ValueError(f"Var {cpm_var} is unknown to the OR-Tools solver, this is unexpected - please report on github...") # translate objective if self.has_objective(): @@ -328,7 +339,7 @@ def transform(self, cpm_expr): :return: list of Expression """ cpm_cons = toplevel_list(cpm_expr) - supported = {"min", "max", "abs", "element", "alldifferent", "xor", "table", "cumulative", "circuit", "inverse"} + supported = {"min", "max", "abs", "element", "alldifferent", "xor", "table", "negative_table", "cumulative", "circuit", "inverse", "no_overlap"} cpm_cons = decompose_in_tree(cpm_cons, supported) cpm_cons = flatten_constraint(cpm_cons) # flat normal form cpm_cons = reify_rewrite(cpm_cons, supported=frozenset(['sum', 'wsum'])) # constraints that support reification @@ -464,10 +475,18 @@ def _post_constraint(self, cpm_expr, reifiable=False): assert (len(cpm_expr.args) == 2) # args = [array, table] array, table = self.solver_vars(cpm_expr.args) return self.ort_model.AddAllowedAssignments(array, table) + elif cpm_expr.name == 'negative_table': + assert (len(cpm_expr.args) == 2) # args = [array, table] + array, table = self.solver_vars(cpm_expr.args) + return self.ort_model.AddForbiddenAssignments(array, table) elif cpm_expr.name == "cumulative": start, dur, end, demand, cap = self.solver_vars(cpm_expr.args) intervals = [self.ort_model.NewIntervalVar(s,d,e,f"interval_{s}-{d}-{e}") for s,d,e in zip(start,dur,end)] return self.ort_model.AddCumulative(intervals, demand, cap) + elif cpm_expr.name == "no_overlap": + start, dur, end = self.solver_vars(cpm_expr.args) + intervals = [self.ort_model.NewIntervalVar(s,d,e, f"interval_{s}-{d}-{d}") for s,d,e in zip(start,dur,end)] + return self.ort_model.add_no_overlap(intervals) elif cpm_expr.name == "circuit": # ortools has a constraint over the arcs, so we need to create these # when using an objective over arcs, using these vars direclty is recommended @@ -555,12 +574,12 @@ def tunable_params(cls): List compiled based on a conversation with OR-tools' Laurent Perron (issue #138). """ return { - 'use_branching_in_lp': [False, True], + # 'use_branching_in_lp': [False, True], # removed in OR-tools >= v9.9 'optimize_with_core' : [False, True], 'search_branching': [0,1,2,3,4,5,6], 'boolean_encoding_level' : [0,1,2,3], 'linearization_level': [0, 1, 2], - 'core_minimization_level' : [0,1,2], # new in OR-tools>=v9.8 + 'core_minimization_level' : [0,1,2], # new in OR-tools>= v9.8 'cp_model_probing_level': [0, 1, 2, 3], 'cp_model_presolve' : [False, True], 'clause_cleanup_ordering' : [0,1], @@ -572,7 +591,7 @@ def tunable_params(cls): @classmethod def default_params(cls): return { - 'use_branching_in_lp': False, + # 'use_branching_in_lp': False, # removed in OR-tools >= v9.9 'optimize_with_core': False, 'search_branching': 0, 'boolean_encoding_level': 1, @@ -687,10 +706,10 @@ def on_solution_callback(self): cpm_var._value = self.Value(self._varmap[cpm_var]) if isinstance(self._display, Expression): - print(self._display.value()) + print(argval(self._display)) elif isinstance(self._display, list): # explicit list of expressions to display - print([v.value() for v in self._display]) + print(argvals(self._display)) else: # callable self._display() diff --git a/cpmpy/solvers/pysat.py b/cpmpy/solvers/pysat.py index 690cbf21a..dba10e431 100644 --- a/cpmpy/solvers/pysat.py +++ b/cpmpy/solvers/pysat.py @@ -6,6 +6,10 @@ """ Interface to PySAT's API + Requires that the 'python-sat' python package is installed: + + $ pip install python-sat[aiger,approxmc,cryptosat,pblib] + PySAT is a Python (2.7, 3.4+) toolkit, which aims at providing a simple and unified interface to a number of state-of-art Boolean satisfiability (SAT) solvers as well as to a variety of cardinality and pseudo-Boolean encodings. @@ -28,6 +32,10 @@ :nosignatures: CPM_pysat + + ============== + Module details + ============== """ from .solver_interface import SolverInterface, SolverStatus, ExitStatus from ..exceptions import NotSupportedError @@ -38,7 +46,7 @@ from ..transformations.decompose_global import decompose_in_tree from ..transformations.get_variables import get_variables from ..transformations.flatten_model import flatten_constraint -from ..transformations.normalize import toplevel_list +from ..transformations.normalize import toplevel_list, simplify_boolean from ..transformations.reification import only_implies, only_bv_reifies @@ -137,6 +145,10 @@ def solve(self, time_limit=None, assumptions=None): For use with s.get_core(): if the model is UNSAT, get_core() returns a small subset of assumption variables that are unsat together. Note: the PySAT interface is statefull, so you can incrementally call solve() with assumptions and it will reuse learned clauses """ + + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + if assumptions is None: pysat_assum_vars = [] # default if no assumptions else: @@ -185,9 +197,8 @@ def solve(self, time_limit=None, assumptions=None): cpm_var._value = True elif -lit in sol: cpm_var._value = False - else: - # not specified... - cpm_var._value = None + else: # not specified, dummy val + cpm_var._value = True return has_sol @@ -230,6 +241,7 @@ def transform(self, cpm_expr): """ cpm_cons = toplevel_list(cpm_expr) cpm_cons = decompose_in_tree(cpm_cons) + cpm_cons = simplify_boolean(cpm_cons) cpm_cons = flatten_constraint(cpm_cons) cpm_cons = only_bv_reifies(cpm_cons) cpm_cons = only_implies(cpm_cons) diff --git a/cpmpy/solvers/pysdd.py b/cpmpy/solvers/pysdd.py index 6fd6238e7..c30bf373b 100644 --- a/cpmpy/solvers/pysdd.py +++ b/cpmpy/solvers/pysdd.py @@ -1,6 +1,15 @@ +#!/usr/bin/env python +#-*- coding:utf-8 -*- +## +## pysdd.py +## """ Interface to PySDD's API + Requires that the 'PySDD' python package is installed: + + $ pip install PySDD + PySDD is a knowledge compilation package for Sentential Decision Diagrams (SDD) https://pysdd.readthedocs.io/en/latest/ @@ -19,6 +28,10 @@ :nosignatures: CPM_pysdd + + ============== + Module details + ============== """ from functools import reduce from .solver_interface import SolverInterface, SolverStatus, ExitStatus @@ -26,7 +39,7 @@ from ..expressions.core import Expression, Comparison, Operator, BoolVal from ..expressions.variables import _BoolVarImpl, NegBoolView, boolvar from ..expressions.globalconstraints import DirectConstraint -from ..expressions.utils import is_any_list, is_bool +from ..expressions.utils import is_any_list, is_bool, argval, argvals from ..transformations.decompose_global import decompose_in_tree from ..transformations.get_variables import get_variables from ..transformations.normalize import toplevel_list, simplify_boolean @@ -94,6 +107,10 @@ def solve(self, time_limit=None, assumptions=None): - building it is the (computationally) hard part - checking for a solution is trivial after that """ + + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + has_sol = True if self.pysdd_root is not None: # if root node is false (empty), no solutions @@ -117,7 +134,8 @@ def solve(self, time_limit=None, assumptions=None): if lit in sol: cpm_var._value = bool(sol[lit]) else: - cpm_var._value = None # not specified... + cpm_var._value = cpm_var.get_bounds()[0] # dummy value - TODO: ensure Pysdd assigns an actual value + # cpm_var._value = None # not specified... return has_sol @@ -135,6 +153,9 @@ def solveAll(self, display=None, time_limit=None, solution_limit=None, call_from Returns: number of solutions found """ + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + if time_limit is not None: raise NotImplementedError("PySDD.solveAll(), time_limit not (yet?) supported") if solution_limit is not None: @@ -143,27 +164,36 @@ def solveAll(self, display=None, time_limit=None, solution_limit=None, call_from if self.pysdd_root is None: return 0 + sddmodels = [x for x in self.pysdd_root.models()] + if len(sddmodels) != self.pysdd_root.model_count: + #pysdd doesn't always have correct solution count.. + projected_sols = set() + for sol in sddmodels: + projectedsol = [] + for cpm_var in self.user_vars: + lit = self.solver_var(cpm_var).literal + projectedsol.append(bool(sol[lit])) + projected_sols.add(tuple(projectedsol)) + else: + projected_sols = set(sddmodels) if display is None: # the desired, fast computation - return self.pysdd_root.model_count() + return len(projected_sols) + else: # manually walking over the tree, much slower... solution_count = 0 - for sol in self.pysdd_root.models(): + for sol in projected_sols: solution_count += 1 # fill in variable values - for cpm_var in self.user_vars: - lit = self.solver_var(cpm_var).literal - if lit in sol: - cpm_var._value = bool(sol[lit]) - else: - cpm_var._value = None + for i, cpm_var in enumerate(self.user_vars): + cpm_var._value = sol[i] # display is not None: if isinstance(display, Expression): - print(display.value()) + print(argval(display)) elif isinstance(display, list): - print([v.value() for v in display]) + print(argvals(display)) else: display() # callback return solution_count diff --git a/cpmpy/solvers/utils.py b/cpmpy/solvers/utils.py index 3d3b2465a..41a180581 100644 --- a/cpmpy/solvers/utils.py +++ b/cpmpy/solvers/utils.py @@ -6,9 +6,6 @@ """ Utilities for handling solvers - Contains a static variable `builtin_solvers` that lists - CPMpy solvers (first one is the default solver by default) - ================= List of functions ================= @@ -26,17 +23,20 @@ from .minizinc import CPM_minizinc from .pysat import CPM_pysat from .z3 import CPM_z3 +from .gcs import CPM_gcs from .pysdd import CPM_pysdd from .exact import CPM_exact +from .choco import CPM_choco def param_combinations(all_params, remaining_keys=None, cur_params=None): """ Recursively yield all combinations of param values For example usage, see `examples/advanced/hyperparameter_search.py` + https://github.com/CPMpy/cpmpy/blob/master/examples/advanced/hyperparameter_search.py - all_params is a dict of {key: list} items, e.g.: - {'val': [1,2], 'opt': [True,False]} + {'val': [1,2], 'opt': [True,False]} - output is an generator over all {key:value} combinations of the keys and values. For the example above: @@ -61,8 +61,8 @@ def param_combinations(all_params, remaining_keys=None, cur_params=None): cur_params=cur_params) class SolverLookup(): - @staticmethod - def base_solvers(): + @classmethod + def base_solvers(cls): """ Return ordered list of (name, class) of base CPMpy solvers @@ -72,16 +72,18 @@ def base_solvers(): return [("ortools", CPM_ortools), ("z3", CPM_z3), ("minizinc", CPM_minizinc), + ("gcs", CPM_gcs), ("gurobi", CPM_gurobi), ("pysat", CPM_pysat), ("pysdd", CPM_pysdd), ("exact", CPM_exact), + ("choco", CPM_choco), ] - @staticmethod - def solvernames(): + @classmethod + def solvernames(cls): names = [] - for (basename, CPM_slv) in SolverLookup.base_solvers(): + for (basename, CPM_slv) in cls.base_solvers(): if CPM_slv.supported(): names.append(basename) if hasattr(CPM_slv, "solvernames"): @@ -90,23 +92,23 @@ def solvernames(): names.append(basename+":"+subn) return names - @staticmethod - def get(name=None, model=None): + @classmethod + def get(cls, name=None, model=None): """ get a specific solver (by name), with 'model' passed to its constructor This is the preferred way to initialise a solver from its name """ - cls = SolverLookup.lookup(name=name) + solver_cls = cls.lookup(name=name) # check for a 'solver:subsolver' name subname = None if name is not None and ':' in name: _,subname = name.split(':',maxsplit=1) - return cls(model, subsolver=subname) + return solver_cls(model, subsolver=subname) - @staticmethod - def lookup(name=None): + @classmethod + def lookup(cls, name=None): """ lookup a solver _class_ by its name @@ -115,7 +117,7 @@ def lookup(name=None): """ if name is None: # first solver class - return SolverLookup.base_solvers()[0][1] + return cls.base_solvers()[0][1] # split name if relevant solvername = name @@ -123,18 +125,16 @@ def lookup(name=None): if ':' in solvername: solvername,_ = solvername.split(':',maxsplit=1) - - for (basename, CPM_slv) in SolverLookup.base_solvers(): + for (basename, CPM_slv) in cls.base_solvers(): if basename == solvername: # found the right solver return CPM_slv - - return None + raise ValueError(f"Unknown solver '{name}', chose from {cls.solvernames()}") # using `builtin_solvers` is DEPRECATED, use `SolverLookup` object instead # Order matters! first is default, then tries second, etc... -builtin_solvers = [CPM_ortools, CPM_gurobi, CPM_minizinc, CPM_pysat, CPM_exact] +builtin_solvers = [CPM_ortools, CPM_gurobi, CPM_minizinc, CPM_pysat, CPM_exact, CPM_choco] def get_supported_solvers(): """ Returns a list of solvers supported on this machine. diff --git a/cpmpy/solvers/z3.py b/cpmpy/solvers/z3.py index 4302c6677..907eb53af 100644 --- a/cpmpy/solvers/z3.py +++ b/cpmpy/solvers/z3.py @@ -1,7 +1,15 @@ #!/usr/bin/env python +#-*- coding:utf-8 -*- +## +## z3.py +## """ Interface to z3's API + Requires that the 'z3-solver' python package is installed: + + $ pip install z3-solver + Z3 is a highly versatile and effective theorem prover from Microsoft. Underneath, it is an SMT solver with a wide scala of theory solvers. We will interface to the finite-domain integer related parts of the API @@ -19,6 +27,10 @@ :nosignatures: CPM_z3 + + ============== + Module details + ============== """ from .solver_interface import SolverInterface, SolverStatus, ExitStatus from ..exceptions import NotSupportedError @@ -112,6 +124,9 @@ def solve(self, time_limit=None, assumptions=[], **kwargs): """ import z3 + # ensure all vars are known to solver + self.solver_vars(list(self.user_vars)) + if time_limit is not None: # z3 expects milliseconds in int self.z3_solver.set(timeout=int(time_limit*1000)) diff --git a/cpmpy/tools/__init__.py b/cpmpy/tools/__init__.py index bba0efc1e..15fc0fdf0 100644 --- a/cpmpy/tools/__init__.py +++ b/cpmpy/tools/__init__.py @@ -1,5 +1,17 @@ """ Set of independent tools that users might appreciate. + + ============= + List of tools + ============= + + .. autosummary:: + :nosignatures: + + explain + dimacs + maximal_propagate + tune_solver """ from .tune_solver import ParameterTuner, GridSearchTuner diff --git a/cpmpy/tools/dimacs.py b/cpmpy/tools/dimacs.py new file mode 100644 index 000000000..0c7e11a57 --- /dev/null +++ b/cpmpy/tools/dimacs.py @@ -0,0 +1,116 @@ +#!/usr/bin/env python +#-*- coding:utf-8 -*- +## +## dimacs.py +## +""" + This file implements helper functions for exporting CPMpy models from and to DIMACS format. + DIMACS is a textual format to represent CNF problems. + The header of the file should be formatted as `p cnf ` + If the number of variables and constraints are not given, it is inferred by the parser. + + Each remaining line of the file is formatted as a list of integers. + An integer represents a Boolean variable and a negative Boolean variable is represented using a `-` sign. +""" + +import cpmpy as cp + +from cpmpy.expressions.variables import _BoolVarImpl, NegBoolView +from cpmpy.expressions.core import Operator, Comparison + +from cpmpy.transformations.normalize import toplevel_list +from cpmpy.transformations.to_cnf import to_cnf +from cpmpy.transformations.get_variables import get_variables + +import re + + +def write_dimacs(model, fname=None): + """ + Writes CPMpy model to DIMACS format + Uses the "to_cnf" transformation from CPMpy + + # TODO: implement pseudoboolean constraints in to_cnf + :param model: a CPMpy model + :fname: optional, file name to write the DIMACS output to + """ + + constraints = toplevel_list(model.constraints) + constraints = to_cnf(constraints) + + vars = get_variables(constraints) + mapping = {v : i+1 for i, v in enumerate(vars)} + + out = f"p cnf {len(vars)} {len(constraints)}\n" + for cons in constraints: + + if isinstance(cons, _BoolVarImpl): + cons = Operator("or", [cons]) + + if not (isinstance(cons, Operator) and cons.name == "or"): + raise NotImplementedError(f"Unsupported constraint {cons}") + + # write clause to cnf format + ints = [] + for v in cons.args: + if isinstance(v, NegBoolView): + ints.append(str(-mapping[v._bv])) + elif isinstance(v, _BoolVarImpl): + ints.append(str(mapping[v])) + else: + raise ValueError(f"Expected Boolean variable in clause, but got {v} which is of type {type(v)}") + + out += " ".join(ints + ["0"]) + "\n" + + if fname is not None: + with open(fname, "w") as f: + f.write(out) + + return out + + +def read_dimacs(fname): + """ + Read a CPMpy model from a DIMACS formatted file + If the number of variables and constraints is not present in the header, they are inferred. + :param: fname: the name of the DIMACS file + :param: sep: optional, separator used in the DIMACS file, will try to infer if None + """ + + m = cp.Model() + + with open(fname, "r") as f: + + lines = f.readlines() + for i, line in enumerate(lines): # DIMACS allows for comments, skip comment lines + if line.startswith("p cnf"): + break + else: + assert line.startswith("c"), f"Expected comment on line {i}, but got {line}" + + cnf = "\n".join(lines[i+1:]) # part of file containing clauses + + bvs = [] + txt_clauses = re.split(r"\n* \n*0", cnf) # clauses end with ` 0` but can have arbitrary newlines + + for txt_ints in txt_clauses: + if txt_ints is None or len(txt_ints.strip()) == 0: + continue # empty clause or weird format + + clause = [] + ints = [int(idx.strip()) for idx in txt_ints.split(" ") if len(idx.strip())] + + for i in ints: + if abs(i) >= len(bvs): # var does not exist yet, create + bvs += [cp.boolvar() for _ in range(abs(i) - len(bvs))] + bv = bvs[abs(i) - 1] + clause.append(bv if i > 0 else ~bv) + + m += cp.any(clause) + + return m + + + + + diff --git a/cpmpy/tools/explain/__init__.py b/cpmpy/tools/explain/__init__.py index b17289035..bc0e2b10f 100644 --- a/cpmpy/tools/explain/__init__.py +++ b/cpmpy/tools/explain/__init__.py @@ -1,3 +1,24 @@ +#!/usr/bin/env python +#-*- coding:utf-8 -*- +## +## __init__.py +## +""" + Collection of tools for explanation techniques. + + ============= + List of tools + ============= + + .. autosummary:: + :nosignatures: + + mus + mss + mcs + utils +""" + from .mus import * from .mss import * from .mcs import * \ No newline at end of file diff --git a/cpmpy/tools/explain/mcs.py b/cpmpy/tools/explain/mcs.py index cde1323b2..28ffa0db9 100644 --- a/cpmpy/tools/explain/mcs.py +++ b/cpmpy/tools/explain/mcs.py @@ -5,7 +5,7 @@ def mcs(soft, hard=[], solver="ortools"): """ Compute Minimal Correction Subset of unsatisfiable model. - Remvoving these contraints will result in a satisfiable model. + Removing these contraints will result in a satisfiable model. Computes a subset of constraints which minimizes the total number of constraints to be removed """ return mcs_opt(soft, hard, 1, solver) diff --git a/cpmpy/tools/explain/utils.py b/cpmpy/tools/explain/utils.py index 9e11e2131..8db2fd691 100644 --- a/cpmpy/tools/explain/utils.py +++ b/cpmpy/tools/explain/utils.py @@ -1,3 +1,21 @@ +#!/usr/bin/env python +#-*- coding:utf-8 -*- +## +## utils.py +## +""" + Utilities for explanation techniques + + ================= + List of functions + ================= + + .. autosummary:: + :nosignatures: + + make_assump_model +""" + import cpmpy as cp from cpmpy.transformations.normalize import toplevel_list @@ -11,7 +29,7 @@ def make_assump_model(soft, hard=[], name=None): soft2 = toplevel_list(soft, merge_and=False) # make assumption variables - assump = cp.boolvar(shape=(len(soft),), name=name) + assump = cp.boolvar(shape=(len(soft2),), name=name) # hard + implied soft constraints model = cp.Model(hard + [assump.implies(soft2)]) # each assumption variable implies a candidate diff --git a/cpmpy/tools/tune_solver.py b/cpmpy/tools/tune_solver.py index e6e37b8e3..c56e1099a 100644 --- a/cpmpy/tools/tune_solver.py +++ b/cpmpy/tools/tune_solver.py @@ -6,9 +6,9 @@ Programming, Artificial Intelligence, and Operations Research. CPAIOR 2022. Lecture Notes in Computer Science, vol 13292. Springer, Cham. https://doi.org/10.1007/978-3-031-08011-1_6 - DOI: https://doi.org/10.1007/978-3-031-08011-1_6 - Link to paper: https://rdcu.be/cQyWR - Link to original code of paper: https://github.com/ML-KULeuven/DeCaprio + - DOI: https://doi.org/10.1007/978-3-031-08011-1_6 + - Link to paper: https://rdcu.be/cQyWR + - Link to original code of paper: https://github.com/ML-KULeuven/DeCaprio This code currently only implements the author's 'Hamming' surrogate function. The parameter tuner iteratively finds better hyperparameters close to the current best configuration during the search. @@ -23,13 +23,15 @@ from ..solvers.solver_interface import ExitStatus class ParameterTuner: + """ + Parameter tuner based on DeCaprio method [ref_to_decaprio] + """ def __init__(self, solvername, model, all_params=None, defaults=None): - """ - Parameter tuner based on DeCaprio method [ref_to_decaprio] - :param: solvername: Name of solver to tune - :param: model: CPMpy model to tune parameters on - :param: all_params: optional, dictionary with parameter names and values to tune. If None, use predefined parameter set. + """ + :param solvername: Name of solver to tune + :param model: CPMpy model to tune parameters on + :param all_params: optional, dictionary with parameter names and values to tune. If None, use predefined parameter set. """ self.solvername = solvername self.model = model @@ -44,9 +46,9 @@ def __init__(self, solvername, model, all_params=None, defaults=None): def tune(self, time_limit=None, max_tries=None, fix_params={}): """ - :param: time_limit: Time budget to run tuner in seconds. Solver will be interrupted when time budget is exceeded - :param: max_tries: Maximum number of configurations to test - :param: fix_params: Non-default parameters to run solvers with. + :param time_limit: Time budget to run tuner in seconds. Solver will be interrupted when time budget is exceeded + :param max_tries: Maximum number of configurations to test + :param fix_params: Non-default parameters to run solvers with. """ if time_limit is not None: start_time = time.time() diff --git a/cpmpy/transformations/__init__.py b/cpmpy/transformations/__init__.py index 516b75c65..7acb14e12 100644 --- a/cpmpy/transformations/__init__.py +++ b/cpmpy/transformations/__init__.py @@ -1,5 +1,5 @@ """ - Methods to transform CPMpy expressions in other CPMpy expressions + Methods to transform CPMpy expressions into other CPMpy expressions Input and output are always CPMpy expressions, so transformations can be chained and called multiple times, as needed. @@ -15,6 +15,13 @@ .. autosummary:: :nosignatures: + comparison + decompose_global flatten_model get_variables + linearize + negation + normalize + reification + to_cnf """ diff --git a/cpmpy/transformations/comparison.py b/cpmpy/transformations/comparison.py index 6f4f40b91..eee30461c 100644 --- a/cpmpy/transformations/comparison.py +++ b/cpmpy/transformations/comparison.py @@ -1,26 +1,30 @@ -import copy - -from .flatten_model import get_or_make_var -from ..expressions.core import Comparison -from ..expressions.variables import _NumVarImpl - """ - Transformations regarding Comparison constraints. + Transformations regarding Comparison constraints (originally). + Now, it is regarding numeric expressions in general, including nested ones. - Comparisons in Flat Normal Form are of the kind - - NumExpr == IV - - BoolExpr == BV + Let with one of == or !=,<,<=,>,>= + Numeric expressions in Flat Normal Form are of the kind + - NumExpr IV + - BoolVar == NumExpr IV + - BoolVar -> NumExpr IV + - NumExpr IV -> BoolVar + + The NumExpr can be a sum, wsum or global function with a non-bool return type. - The latter is a reified expression, not considered here. - (for handling of all reified expressions, see `reification.py` transformations) - This file implements: - - only_numexpr_equality(): transforms `NumExpr IV` to `(NumExpr == A) & (A IV)` if not supported + - only_numexpr_equality(): transforms `NumExpr IV` (also reified) to `(NumExpr == A) & (A IV)` if not supported """ +import copy +from .flatten_model import get_or_make_var +from ..expressions.core import Comparison, Operator +from ..expressions.utils import is_boolexpr +from ..expressions.variables import _NumVarImpl, _BoolVarImpl + def only_numexpr_equality(constraints, supported=frozenset()): """ transforms `NumExpr IV` to `(NumExpr == A) & (A IV)` if not supported + also for the reified uses of NumExpr :param supported a (frozen)set of expression names that supports all comparisons in the solver """ @@ -28,17 +32,54 @@ def only_numexpr_equality(constraints, supported=frozenset()): # shallow copy (could support inplace too this way...) newcons = copy.copy(constraints) - for i,con in enumerate(newcons): - if isinstance(con, Comparison) and con.name != '==': - # LHS IV with one of !=,<,<=,>,>= - lhs = con.args[0] - if not isinstance(lhs, _NumVarImpl) and not lhs.name in supported: - # LHS is unsupported for LHS IV, rewrite to `(LHS == A) & (A IV)` - (lhsvar, lhscons) = get_or_make_var(lhs) - # replace comparison by A IV - newcons[i] = Comparison(con.name, lhsvar, con.args[1]) - # add lhscon(s), which will be [(LHS == A)] - assert(len(lhscons) == 1), "only_numexpr_eq: lhs surprisingly non-flat" - newcons.insert(i, lhscons[0]) + for i,cpm_expr in enumerate(newcons): + + if isinstance(cpm_expr, Operator) and cpm_expr.name == "->": + cond, subexpr = cpm_expr.args + if not isinstance(cond, _BoolVarImpl): # expr -> bv + res = only_numexpr_equality([cond], supported) + if len(res) > 1: + newcons[i] = res[1].implies(subexpr) + newcons.insert(i, res[0]) + + elif not isinstance(subexpr, _BoolVarImpl): # bv -> expr + res = only_numexpr_equality([subexpr], supported) + if len(res) > 1: + newcons[i] = cond.implies(res[1]) + newcons.insert(i, res[0]) + else: #bv -> bv + pass + + + if isinstance(cpm_expr, Comparison): + lhs, rhs = cpm_expr.args + + if cpm_expr.name == "==" and is_boolexpr(lhs) and is_boolexpr(rhs): # reification, check recursively + + if not isinstance(lhs, _BoolVarImpl): # expr == bv + res = only_numexpr_equality([lhs], supported) + if len(res) > 1: + newcons[i] = res[1] == rhs + newcons.insert(i, res[0]) + + elif not isinstance(rhs, _BoolVarImpl): # bv == expr + res = only_numexpr_equality([rhs], supported) + if len(res) > 1: + newcons[i] = lhs == res[1] + newcons.insert(i, res[0]) + else: # bv == bv + pass + + elif cpm_expr.name != "==": + # LHS IV with one of !=,<,<=,>,>= + lhs = cpm_expr.args[0] + if not isinstance(lhs, _NumVarImpl) and lhs.name not in supported: + # LHS is unsupported for LHS IV, rewrite to `(LHS == A) & (A IV)` + (lhsvar, lhscons) = get_or_make_var(lhs) + # replace comparison by A IV + newcons[i] = Comparison(cpm_expr.name, lhsvar, cpm_expr.args[1]) + # add lhscon(s), which will be [(LHS == A)] + assert(len(lhscons) == 1), "only_numexpr_eq: lhs surprisingly non-flat" + newcons.insert(i, lhscons[0]) return newcons diff --git a/cpmpy/transformations/decompose_global.py b/cpmpy/transformations/decompose_global.py index 48e8028ac..79228d0a9 100644 --- a/cpmpy/transformations/decompose_global.py +++ b/cpmpy/transformations/decompose_global.py @@ -27,11 +27,12 @@ def decompose_in_tree(lst_of_expr, supported=set(), supported_reified=set(), _to Supported numerical global functions remain in the expression tree as is. They can be rewritten using `cpmpy.transformations.reification.reify_rewrite` The following `bv -> NumExpr Var/Const` can be rewritten as [bv -> IV0 Var/Const, NumExpr == IV0]. - So even if numerical constraints are not supported in reified context, we can rewrite them to non-reified versions. + So even if numerical constraints are not supported in reified context, we can rewrite them to non-reified versions if they are total. """ if _toplevel is None: _toplevel = [] + # swap the arguments of a comparison while maintaining its semantics flipmap = {"==": "==", "!=": "!=", "<": ">", "<=": ">=", ">": "<", ">=": "<="} newlist = [] # decomposed constraints will go here @@ -78,7 +79,7 @@ def decompose_in_tree(lst_of_expr, supported=set(), supported_reified=set(), _to _toplevel.extend(define) # definitions should be added toplevel # the `decomposed` expression might contain other global constraints, check it - decomposed = decompose_in_tree(decomposed, supported, supported_reified, [], nested=nested) + decomposed = decompose_in_tree(decomposed, supported, supported_reified, _toplevel, nested=nested) newlist.append(all(decomposed)) else: @@ -145,7 +146,7 @@ def decompose_in_tree(lst_of_expr, supported=set(), supported_reified=set(), _to return toplevel_list(newlist) else: # we are toplevel and some new constraints are introduced, decompose new constraints! - return toplevel_list(newlist) + decompose_in_tree(_toplevel, supported, supported_reified, nested=False) + return toplevel_list(newlist) + decompose_in_tree(toplevel_list(_toplevel), supported, supported_reified, nested=False) # DEPRECATED! diff --git a/cpmpy/transformations/flatten_model.py b/cpmpy/transformations/flatten_model.py index 401e7f2d5..29b93acd5 100644 --- a/cpmpy/transformations/flatten_model.py +++ b/cpmpy/transformations/flatten_model.py @@ -245,7 +245,12 @@ def flatten_constraint(expr): # Reification (double implication): Boolexpr == Var # normalize the lhs (does not have to be a var, hence we call normalize instead of get_or_make_var if exprname == '==' and lexpr.is_bool(): - (lhs, lcons) = normalized_boolexpr(lexpr) + if rvar.is_bool(): + # this is a reification + (lhs, lcons) = normalized_boolexpr(lexpr) + else: + # integer comparison + (lhs, lcons) = get_or_make_var(lexpr) else: (lhs, lcons) = normalized_numexpr(lexpr) @@ -319,10 +324,10 @@ def get_or_make_var(expr): (flatexpr, flatcons) = normalized_boolexpr(expr) if isinstance(flatexpr,_BoolVarImpl): - #avoids unnecessary bv == bv or bv == ~bv assignments + # avoids unnecessary bv == bv or bv == ~bv assignments return flatexpr,flatcons bvar = _BoolVarImpl() - return (bvar, [flatexpr == bvar]+flatcons) + return bvar, [flatexpr == bvar] + flatcons else: # normalize expr into a numexpr LHS, @@ -330,8 +335,12 @@ def get_or_make_var(expr): (flatexpr, flatcons) = normalized_numexpr(expr) lb, ub = flatexpr.get_bounds() + if not(isinstance(lb,int) and isinstance(ub,int)): + warnings.warn("CPMPy only uses integer variables, non-integer expression detected that will be reified " + "into an intvar with rounded bounds. \n Your constraints will stay the same.", UserWarning) + lb, ub = math.floor(lb), math.ceil(ub) ivar = _IntVarImpl(lb, ub) - return (ivar, [flatexpr == ivar]+flatcons) + return ivar, [flatexpr == ivar] + flatcons def get_or_make_var_or_list(expr): """ Like get_or_make_var() but also accepts and recursively transforms lists diff --git a/cpmpy/transformations/linearize.py b/cpmpy/transformations/linearize.py index 63ad3bb68..d7cc4dd53 100644 --- a/cpmpy/transformations/linearize.py +++ b/cpmpy/transformations/linearize.py @@ -222,7 +222,7 @@ def only_positive_bv(lst_of_expr): weights, args = lhs.args idxes = {i for i, a in enumerate(args) if isinstance(a, NegBoolView)} nw, na = zip(*[(-w,a._bv) if i in idxes else (w,a) for i, (w,a) in enumerate(zip(weights, args))]) - lhs = Operator("wsum", [nw, na]) # force making wsum, even for arity = 1 + lhs = Operator("wsum", [list(nw), list(na)]) # force making wsum, even for arity = 1 rhs -= sum(weights[i] for i in idxes) if isinstance(lhs, Operator) and lhs.name not in {"sum","wsum"}: @@ -269,8 +269,10 @@ def canonical_comparison(lst_of_expr): elif isinstance(lhs, Comparison): lhs = canonical_comparison(lhs)[0] newlist.append(lhs.implies(rhs)) + else: + newlist.append(cpm_expr) - if isinstance(cpm_expr, Comparison): + elif isinstance(cpm_expr, Comparison): lhs, rhs = cpm_expr.args if isinstance(lhs, Comparison) and cpm_expr.name == "==": # reification of comparison lhs = canonical_comparison(lhs)[0] @@ -299,24 +301,30 @@ def canonical_comparison(lst_of_expr): assert not is_num(lhs), "lhs cannot be an integer at this point!" # bring all const to rhs - if lhs.name == "sum": - new_args = [] - for i, arg in enumerate(lhs.args): - if is_num(arg): - rhs -= arg - else: - new_args.append(arg) - lhs = Operator("sum", new_args) - - elif lhs.name == "wsum": - new_weights, new_args = [], [] - for i, (w, arg) in enumerate(zip(*lhs.args)): - if is_num(arg): - rhs -= w * arg - else: - new_weights.append(w) - new_args.append(arg) - lhs = Operator("wsum", [new_weights, new_args]) + if isinstance(lhs, Operator): + if lhs.name == "sum": + new_args = [] + for i, arg in enumerate(lhs.args): + if is_num(arg): + rhs -= arg + else: + new_args.append(arg) + lhs = Operator("sum", new_args) + + elif lhs.name == "wsum": + new_weights, new_args = [], [] + for i, (w, arg) in enumerate(zip(*lhs.args)): + if is_num(arg): + rhs -= w * arg + else: + new_weights.append(w) + new_args.append(arg) + lhs = Operator("wsum", [new_weights, new_args]) + else: + raise ValueError(f"lhs should be sum or wsum, but got {lhs}") + else: + assert isinstance(lhs, _NumVarImpl) + lhs = Operator("sum", [lhs]) newlist.append(eval_comparison(cpm_expr.name, lhs, rhs)) else: # rest of expressions diff --git a/cpmpy/transformations/normalize.py b/cpmpy/transformations/normalize.py index ec0c5e756..b5cf35d17 100644 --- a/cpmpy/transformations/normalize.py +++ b/cpmpy/transformations/normalize.py @@ -124,6 +124,8 @@ def simplify_boolean(lst_of_expr, num_context=False): """ if is_boolexpr(lhs) and is_num(rhs): # direct simplification of boolean comparisons + if isinstance(rhs, BoolVal): + rhs = int(rhs.value()) # ensure proper comparisons below if rhs < 0: newlist.append(BoolVal(name in {"!=", ">", ">="})) # all other operators evaluate to False elif rhs == 0: @@ -132,15 +134,15 @@ def simplify_boolean(lst_of_expr, num_context=False): if name == "==" or name == "<=": newlist.append(recurse_negation(lhs)) if name == "<": - newlist.append(BoolVal(False)) + newlist.append(0 if num_context else BoolVal(False)) if name == ">=": - newlist.append(BoolVal(True)) + newlist.append(1 if num_context else BoolVal(True)) elif 0 < rhs < 1: # a floating point value if name == "==": - newlist.append(BoolVal(False)) + newlist.append(0 if num_context else BoolVal(False)) if name == "!=": - newlist.append(BoolVal(True)) + newlist.append(1 if num_context else BoolVal(True)) if name == "<" or name == "<=": newlist.append(recurse_negation(lhs)) if name == ">" or name == ">=": @@ -151,9 +153,9 @@ def simplify_boolean(lst_of_expr, num_context=False): if name == "!=" or name == "<": newlist.append(recurse_negation(lhs)) if name == ">": - newlist.append(BoolVal(False)) + newlist.append(0 if num_context else BoolVal(False)) if name == "<=": - newlist.append(BoolVal(True)) + newlist.append(1 if num_context else BoolVal(True)) elif rhs > 1: newlist.append(BoolVal(name in {"!=", "<", "<="})) # all other operators evaluate to False else: diff --git a/cpmpy/transformations/reification.py b/cpmpy/transformations/reification.py index 3861398c1..21b41b769 100644 --- a/cpmpy/transformations/reification.py +++ b/cpmpy/transformations/reification.py @@ -1,13 +1,3 @@ -import copy -from ..expressions.core import Operator, Comparison, Expression -from ..expressions.globalconstraints import GlobalConstraint -from ..expressions.globalfunctions import Element -from ..expressions.variables import _BoolVarImpl, _NumVarImpl -from ..expressions.python_builtins import all -from ..expressions.utils import is_any_list -from .flatten_model import flatten_constraint, get_or_make_var -from .negation import recurse_negation - """ Transformations regarding reification constraints. @@ -23,6 +13,15 @@ - only_implies(): transforms all reifications to BV -> BE form - reify_rewrite(): rewrites reifications not supported by a solver to ones that are """ +import copy +from ..expressions.core import Operator, Comparison, Expression +from ..expressions.globalconstraints import GlobalConstraint +from ..expressions.globalfunctions import Element +from ..expressions.variables import _BoolVarImpl, _NumVarImpl +from ..expressions.python_builtins import all +from ..expressions.utils import is_any_list +from .flatten_model import flatten_constraint, get_or_make_var +from .negation import recurse_negation def only_bv_reifies(constraints): newcons = [] diff --git a/cpmpy/transformations/to_cnf.py b/cpmpy/transformations/to_cnf.py index 918d9cc4e..c982b5276 100644 --- a/cpmpy/transformations/to_cnf.py +++ b/cpmpy/transformations/to_cnf.py @@ -1,7 +1,3 @@ -from ..expressions.core import Operator, Comparison -from ..expressions.variables import _BoolVarImpl, NegBoolView -from .reification import only_implies -from .flatten_model import flatten_constraint """ Converts the logical constraints into disjuctions using the tseitin transform, including flattening global constraints that are is_bool() and not in `supported`. @@ -22,6 +18,10 @@ - BE -> BV - BV -> BE """ +from ..expressions.core import Operator, Comparison +from ..expressions.variables import _BoolVarImpl, NegBoolView +from .reification import only_implies +from .flatten_model import flatten_constraint def to_cnf(constraints): """ diff --git a/docs/adding_solver.md b/docs/adding_solver.md index f579f3bea..2b584de60 100644 --- a/docs/adding_solver.md +++ b/docs/adding_solver.md @@ -12,7 +12,7 @@ Implementing the template consists of the following parts: * `solve()` where you call the solver, get the status and runtime, and reverse-map the variable values after solving * `objective()` if your solver supports optimisation * `transform()` where you call the necessary transformations in `cpmpy.transformations` to transform CPMpy expressions to those that the solver supports - * `__add__()` where you call transform and map the resulting CPMpy expressions that the solver supports, to API function calls on the underlying solver + * `__add__()` where you call transform and map the resulting CPMpy expressions, that the solver supports, to API function calls on the underlying solver * `solveAll()` optionally, if the solver natively supports solution enumeration ## Transformations and posting constraints @@ -29,13 +29,13 @@ So for any solver you wish to add, chances are that most of the transformations ## Stateless transformation functions -Because CPMpy solver interfaces transform and post constraints *eagerly*, they can be used *incremental*, meaning that you can add some constraints, call `solve()` add some more constraints and solve again. If the underlying solver is also incremental, it will reuse knowledge of the previous solve call to speed up this solve call. +Because CPMpy's solver-interfaces transform and post constraints *eagerly*, they can be used *incremental*, meaning that you can add some constraints, call `solve()` add some more constraints and solve again. If the underlying solver is also incremental, it will reuse knowledge of the previous solve call to speed up this solve call. The way that CPMpy succeeds to be an incremental modeling language, is by making all transformation functions *stateless*. Every transformation function is a python *function* that maps a (list of) CPMpy expressions to (a list of) equivalent CPMpy expressions. Transformations are not classes, they do not store state, they do not know (or care) what model a constraint belongs to. They take expressions as input and compute expressions as output. That means they can be called over and over again, and chained in any combination or order. That also makes them modular, and any solver can use any combination of transformations that it needs. We continue to add and improve the transformations, and we are happy to discuss transformations you are missing, or variants of existing transformations that can be refined. -Most transformations do not need any state, they just do a bit of rewriting. Some transformations do, for example in the case of common subexpression elimination. In that case, the solver interface (you who are reading this), should store a dictionary in your solver interface class, and pass that as (optional) argument to the transformation function. The transformation function will read and write to that dictionary as it needs, while still remaining stateless on its own. Each transformation function documents when it supports an optional state dictionary, see all available transformations in `cpmpy/transformations/`. +Most transformations do not need any state, they just do a bit of rewriting. Some transformations do, for example in the case of [Common Subexpression Elimination (CSE)](https://en.wikipedia.org/wiki/Common_subexpression_elimination). In that case, the solver interface (you who are reading this), should store a dictionary in your solver interface class, and pass that as (optional) argument to the transformation function. The transformation function will read and write to that dictionary as it needs, while still remaining stateless on its own. Each transformation function documents when it supports an optional state dictionary, see all available transformations in `cpmpy/transformations/`. ## What is a good Python interface for a solver? @@ -87,7 +87,7 @@ class SolverX { } ``` -If you have such a C++ API, then there exist automatic python packages that can make Python bindings, such as [CPPYY](https://cppyy.readthedocs.io/en/latest/). +If you have such a C++ API, then there exist automatic python packages that can make Python bindings, such as [CPPYY](https://cppyy.readthedocs.io/en/latest/) or [pybind11](https://pybind11.readthedocs.io/en/stable/). We have not done this ourselves yet, so get in touch to share your experience and advice! @@ -102,5 +102,6 @@ After posting the constraint, the answer of your solver is checked so you will b ## Tunable hyperparameters CPMpy offers a tool for searching the best hyperparameter configuration for a given model on a solver (see [corresponding documentation](solver_parameters.md)). -Solver wanting to support this tool should add the following attributes to their solver interface: `tunable_params` and `default_params` (see [ortools](https://github.com/CPMpy/cpmpy/blob/11ae35b22357ad9b8d6f47317df2c236c3ef5997/cpmpy/solvers/ortools.py#L473) for an example). +Solvers wanting to support this tool should add the following attributes to their interface: `tunable_params` and `default_params` (see [OR-Tools](https://github.com/CPMpy/cpmpy/blob/11ae35b22357ad9b8d6f47317df2c236c3ef5997/cpmpy/solvers/ortools.py#L473) for an example). + diff --git a/docs/api/expressions/globalfunctions.rst b/docs/api/expressions/globalfunctions.rst new file mode 100644 index 000000000..8240d85b0 --- /dev/null +++ b/docs/api/expressions/globalfunctions.rst @@ -0,0 +1,7 @@ +Global Functions (:mod:`cpmpy.expressions.globalfunctions`) +=================================================== + +.. automodule:: cpmpy.expressions.globalfunctions + :members: + :undoc-members: + :inherited-members: diff --git a/docs/api/solvers/choco.rst b/docs/api/solvers/choco.rst new file mode 100644 index 000000000..66f739ed8 --- /dev/null +++ b/docs/api/solvers/choco.rst @@ -0,0 +1,7 @@ +CPMpy exact interface (:mod:`cpmpy.solvers.choco`) +===================================================== + +.. automodule:: cpmpy.solvers.choco + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/tools.rst b/docs/api/tools.rst new file mode 100644 index 000000000..022761aa1 --- /dev/null +++ b/docs/api/tools.rst @@ -0,0 +1,6 @@ +Tools (:mod:`cpmpy.tools`) +======================================== + +.. automodule:: cpmpy.tools + :members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/tools/dimacs.rst b/docs/api/tools/dimacs.rst new file mode 100644 index 000000000..8b785778c --- /dev/null +++ b/docs/api/tools/dimacs.rst @@ -0,0 +1,7 @@ +DIMACS (:mod:`cpmpy.tools.dimacs`) +===================================================== + +.. automodule:: cpmpy.tools.dimacs + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/tools/explain.rst b/docs/api/tools/explain.rst new file mode 100644 index 000000000..22f15dd12 --- /dev/null +++ b/docs/api/tools/explain.rst @@ -0,0 +1,33 @@ +Explanation tools (:mod:`cpmpy.tools.explain`) +===================================================== + +.. automodule:: cpmpy.tools.explain + :members: + :undoc-members: + :inherited-members: + +.. include:: ./explain/mus.rst +.. include:: ./explain/mss.rst +.. include:: ./explain/mcs.rst +.. include:: ./explain/utils.rst + +.. .. automodule:: cpmpy.tools.explain.mcs +.. :title: +.. :members: +.. :undoc-members: +.. :inherited-members: + +.. .. automodule:: cpmpy.tools.explain.mss +.. :members: +.. :undoc-members: +.. :inherited-members: + +.. .. automodule:: cpmpy.tools.explain.mus +.. :members: +.. :undoc-members: +.. :inherited-members: + +.. .. automodule:: cpmpy.tools.explain.utils +.. :members: +.. :undoc-members: +.. :inherited-members: \ No newline at end of file diff --git a/docs/api/tools/explain/mcs.rst b/docs/api/tools/explain/mcs.rst new file mode 100644 index 000000000..f098d0dc1 --- /dev/null +++ b/docs/api/tools/explain/mcs.rst @@ -0,0 +1,7 @@ +MCS (:mod:`cpmpy.tools.explain.mcs`) +===================================================== + +.. automodule:: cpmpy.tools.explain.mcs + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/tools/explain/mss.rst b/docs/api/tools/explain/mss.rst new file mode 100644 index 000000000..cc77290e2 --- /dev/null +++ b/docs/api/tools/explain/mss.rst @@ -0,0 +1,7 @@ +MSS (:mod:`cpmpy.tools.explain.mss`) +===================================================== + +.. automodule:: cpmpy.tools.explain.mss + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/tools/explain/mus.rst b/docs/api/tools/explain/mus.rst new file mode 100644 index 000000000..8fdcd532c --- /dev/null +++ b/docs/api/tools/explain/mus.rst @@ -0,0 +1,7 @@ +MUS (:mod:`cpmpy.tools.explain.mus`) +===================================================== + +.. automodule:: cpmpy.tools.explain.mus + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/tools/explain/utils.rst b/docs/api/tools/explain/utils.rst new file mode 100644 index 000000000..a4a44657f --- /dev/null +++ b/docs/api/tools/explain/utils.rst @@ -0,0 +1,7 @@ +Utils (:mod:`cpmpy.tools.explain.utils`) +===================================================== + +.. automodule:: cpmpy.tools.explain.utils + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/tools/miximal_propagate.rst b/docs/api/tools/miximal_propagate.rst new file mode 100644 index 000000000..9fe5dad3a --- /dev/null +++ b/docs/api/tools/miximal_propagate.rst @@ -0,0 +1,7 @@ +Maximal Propagate (:mod:`cpmpy.tools.maximal_propagate`) +===================================================== + +.. automodule:: cpmpy.tools.maximal_propagate + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/tools/tune_solver.rst b/docs/api/tools/tune_solver.rst new file mode 100644 index 000000000..503b3eb36 --- /dev/null +++ b/docs/api/tools/tune_solver.rst @@ -0,0 +1,7 @@ +Solver Parameter Tuning (:mod:`cpmpy.tools.tune_solver`) +===================================================== + +.. automodule:: cpmpy.tools.tune_solver + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/transformations.rst b/docs/api/transformations.rst index 29174f126..266ea5802 100644 --- a/docs/api/transformations.rst +++ b/docs/api/transformations.rst @@ -3,4 +3,15 @@ Expression transformations (:mod:`cpmpy.transformations`) .. automodule:: cpmpy.transformations :members: + :undoc-members: :inherited-members: + +.. .. include:: ./transformations/comparison.rst +.. .. include:: ./transformations/decompose_global.rst +.. .. include:: ./transformations/flatten_model.rst +.. .. include:: ./transformations/get_variables.rst +.. .. include:: ./transformations/linearize.rst +.. .. include:: ./transformations/negation.rst +.. .. include:: ./transformations/normalize.rst +.. .. include:: ./transformations/reification.rst +.. .. include:: ./transformations/to_cnf.rst \ No newline at end of file diff --git a/docs/api/transformations/comparison.rst b/docs/api/transformations/comparison.rst new file mode 100644 index 000000000..a534380e7 --- /dev/null +++ b/docs/api/transformations/comparison.rst @@ -0,0 +1,9 @@ +.. orphan:: + +Comparison (:mod:`cpmpy.transformations.comparison`) +======================================================================== + +.. automodule:: cpmpy.transformations.comparison + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/transformations/decompose_global.rst b/docs/api/transformations/decompose_global.rst new file mode 100644 index 000000000..e8f6f9973 --- /dev/null +++ b/docs/api/transformations/decompose_global.rst @@ -0,0 +1,7 @@ +Decompose Global (:mod:`cpmpy.transformations.decompose_global`) +======================================================================== + +.. automodule:: cpmpy.transformations.decompose_global + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/transformations/negation.rst b/docs/api/transformations/negation.rst new file mode 100644 index 000000000..9496ad5ef --- /dev/null +++ b/docs/api/transformations/negation.rst @@ -0,0 +1,7 @@ +Negation (:mod:`cpmpy.transformations.negation`) +======================================================================== + +.. automodule:: cpmpy.transformations.negation + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/api/transformations/normalize.rst b/docs/api/transformations/normalize.rst new file mode 100644 index 000000000..1c5f0dfc1 --- /dev/null +++ b/docs/api/transformations/normalize.rst @@ -0,0 +1,7 @@ +Normalize (:mod:`cpmpy.transformations.normalize`) +======================================================================== + +.. automodule:: cpmpy.transformations.normalize + :members: + :undoc-members: + :inherited-members: \ No newline at end of file diff --git a/docs/conf.py b/docs/conf.py index 11328e77d..7f5b1daa6 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -24,7 +24,7 @@ author = 'Tias Guns' # The full version, including alpha/beta/rc tags -release = '0.9.18' +release = '0.9.21' # variables to be accessed from html html_context = { diff --git a/docs/developers.md b/docs/developers.md index f8d4db9af..c3ee62e25 100644 --- a/docs/developers.md +++ b/docs/developers.md @@ -34,27 +34,28 @@ python -m pytest tests/test_solvers.py ## Code structure - * _tests/_ contains the tests - * _docs/_ contains the docs. Any change there is automatically updated, with some delay, on [https://cpmpy.readthedocs.io/](https://cpmpy.readthedocs.io/) - * _examples/_ our examples, always happy to include more - * _cpmpy/_ the python module that you install when doing `pip install cpmpy` + * _tests/_ : contains the tests + * _docs/_ : contains the docs. Any change there is automatically updated, with some delay, on [https://cpmpy.readthedocs.io/](https://cpmpy.readthedocs.io/) + * _examples/_ : our examples, always happy to include more + * _cpmpy/_ : the python module that you install when doing `pip install cpmpy` The module is structured as such: - * _model.py_ contains the omnipresent `Model()` container - * _expressions/_ Classes and functions that represent and create expressions (constraints and objectives) - * _solvers/_ CPMpy interfaces to (the Python API interface of) solvers - * _transformations/_ Methods to transform CPMpy expressions in other CPMpy expressions - * _tools/_ Set of independent tools that users might appreciate. + * _model.py_ : contains the omnipresent `Model()` container + * _exceptions.py_ : contains a collection of CPMpy specific exceptions + * _expressions/_ : Classes and functions that represent and create expressions (constraints and objectives) + * _solvers/_ : CPMpy interfaces to (the Python API interface of) solvers + * _transformations/_ : Methods to transform CPMpy expressions into other CPMpy expressions + * _tools/_ : Set of independent tools that users might appreciate. -The typical flow in which these submodules are used when using CPMpy is: the user creates _expressions_ which they put into a _model_ object. This is then given to a _solver_ object to solve, which will first _transform_ the expressions into expressions that it supports, which it then posts to the Python API interface of that particular solver. +The typical flow in which these submodules are used when programming with CPMpy is: the user creates _expressions_ which they put into a _model_ object. This is then given to a _solver_ object to solve, which will first _transform_ the original expressions into expressions that it supports, which it then posts to the Python API interface of that particular solver. Tools are not part of the core of CPMpy. They are additional tools that _use_ CPMpy, e.g. for debugging, parameter tuning etc. -## Github practices +## GitHub practices -When filing a bug, please add a small case that allows us to reproduce it. If the testcase is confidential, mail Tias directly. +When filing a bug, please add a small case that allows us to reproduce it. If the testcase is confidential, mail [Tias](mailto:tias.guns@kuleuven.be) directly. Only documentation changes can be directly applied on the master branch. All other changes should be submitted as a pull request. diff --git a/docs/docs_todo.md b/docs/docs_todo.md new file mode 100644 index 000000000..dae494d6d --- /dev/null +++ b/docs/docs_todo.md @@ -0,0 +1,33 @@ +# TODO + +A list of things in CPMpy's docs which need fixed / which are optional improvements. + +A lot of proof-reading is still needed to catch all formatting mistakes (especially in the API section, as taken for the code's docstrings) + + +- [ ] Outdated copyright + - whilst updated on the repo, ReadTheDocs still gives 2021 as copyright year + +- [ ] Mention Exact + - In many places where solvers get listed, Exact is not inluded (especially for things that make Exact stand out; e.g. unsat core extraction) + +- [ ] Solver overview table + - available solvers and their supported features + e.g. incrementality, core exctraction, proof logging, parallelisation, ... + +- [ ] Overview of available CPMpy tools +- [ ] "frietkot" link in ocus is down +- [ ] docs already refer to cse in "how to add a solver" + - CSE has not been added to master yet + +- [ ] Decision variables shape + - decision variable docs say that all variables are arrays and thus have the shape attribute + -> only NDVarArrays + +- [ ] Sphinx's python autodoc annotations: + - see Exact's interface + - [The Python Domain — Sphinx documentation (sphinx-doc.org)](https://www.sphinx-doc.org/en/master/usage/domains/python.html#directive-py-method) + +- [ ] Gurobi GRB_ENV subsolver + document this + diff --git a/docs/how_to_debug.md b/docs/how_to_debug.md index 5c9f5d253..1def89b6f 100644 --- a/docs/how_to_debug.md +++ b/docs/how_to_debug.md @@ -7,7 +7,7 @@ The bug can be situated in one of three layers: - the CPMpy library - the solver -coincidentally, they are ordered from most likely to least likely. So let's start at the bottom. +Coincidentally, they are ordered from most likely to least likely. So let's start at the bottom. If you don't have a bug yet, but are curious, here is some general advise from expert modeller [Håkan Kjellerstrand](http://www.hakank.org/): - Test the model early and often. This makes it easier to detect problems in the model. @@ -21,22 +21,23 @@ If you get an error and have difficulty understanding it, try searching on the i If you don't find it, or if the solver runs fine and without error, but you don't get the answer you expect; then try swapping out the solver for another solver and see what gives... -Replace `model.solve()` by `model.solve(solver='minizinc')` for example. You do need to install MiniZinc and minizinc-python first though. +Replace `model.solve()` by `model.solve(solver='minizinc')` for example. You do need to install MiniZinc and minizinc-python first though. Take a look at [the solver API interface](api/solvers.html) for the install instructions. Either you have the same output, and it is not the solver's fault, or you have a different output and you actually found one of these rare solver bugs. Report on the bugtracker of the solver, or on the CPMpy github page where we will help you file a bug 'upstream' (or maybe even work around it in CPMpy). ## Debugging a modeling error -You get an error when you create an expression? Then you are probably writing it wrongly. Check the documentation and the running examples for similar examples of what you wish to express. +You get an error when you create an expression? Then you are probably writing it wrongly. Check the documentation and the running examples for similar instances of what you wish to express. Here are a few quirks in Python/CPMpy: - when using `&` and `|`, make sure to always put the subexpressions in brackets. E.g. `(x == 1) & (y == 0)` instead of `x == 1 & y == 0`. The latter wont work, because Python will unfortunately think you meant `x == (1 & (y == 0))`. - you can write `vars[other_var]` but you can't write `non_var_list[a_var]`. That is because the `vars` list knows CPMpy, and the `non_var_list` does not. Wrap it: `non_var_list = cpm_array(non_var_list)` first, or write `Element(non_var_list, a_var)` instead. - only write `sum(v)` on lists, don't write it if `v` is a matrix or tensor, as you will get a list in response. Instead, use NumPy's `v.sum()` instead. + - when providing names for decision variables, make shure that they are unique. Many solvers depend on this uniqueness and you will encounter strange (and hard to debug) behaviour if you don't enforce this. Try printing the expression `print(e)` or subexpressions, and check that the output matches what you wish to express. Decompose the expression and try printing the individual components and their piecewice composition to see what works and when it starts to break. -If you don't find it, report it on the CPMpy github Issues page and we'll help you (and maybe even extend the above list of quirks). +If you don't find it, report it on the CPMpy GitHub Issues page and we'll help you (and maybe even extend the above list of quirks). ## Debugging a `solve()` error @@ -52,9 +53,9 @@ for c in model.constraints: Model(c).solve() ``` -The last constraint printed before the exception is the curlpit... Please report on Github. We want to catch corner cases in CPMpy, even if it is a solver limitation, so please report on the CPMpy github Issues page. +The last constraint printed before the exception is the curlpit... Please report on GitHub. We want to catch corner cases in CPMpy, even if it is a solver limitation, so please report on the CPMpy GitHub Issues page. -Or maybe, you got one of CPMpy's NotImplementedErrors. Share your use case with us on Github and we will implement it. Or implemented it yourself first, that is also very welcome ; ) +Or maybe, you got one of CPMpy's NotImplementedErrors. Share your use case with us on GitHub and we will implement it. Or implemented it yourself first, that is also very welcome ; ) ## Debugging an UNSATisfiable model @@ -62,7 +63,7 @@ First, print the model: ```print(model)``` -and check that the output matches what you want to express. Do you see anything unusual? Start there, see why the expression is not what you intended to express, as described in 'Debugging a modeling error'. +and check that the output matches what you want to express. Do you see anything unusual? Start there, see why the expression is not what you intended to express, as described in [Debugging a modeling error](#debugging-a-modeling-error). If that does not help, try printing the 'transformed' **constraints**, the way that the solver actually sees them, including decompositions and rewrites: @@ -95,7 +96,7 @@ print(f"Optimizing {obj_var} subject to", s.transform(obj_expr)) ### Automatically minimising the UNSAT program If the above is unwieldy because your constraint problem is too large, then consider automatically reducing it to a 'Minimal Unsatisfiable Subset' (MUS). -This is now part of our standard tools, that you can use as follows: +This is now part of our [standard tools](api/tools.html), that you can use as follows: ```python from cpmpy.tools import mus @@ -113,7 +114,7 @@ unsat_cons = mus(model.constraints) With this smaller set of constraints, repeat the visual inspection steps above. -(Note that for an UNSAT problem there can be many MUSes, the `examples/advanced/` folder has the MARCO algorithm that can enumerate all MSS/MUSes.) +(Note that for an UNSAT problem there can be many MUSes, the `examples/advanced/` [folder](https://github.com/CPMpy/cpmpy/tree/master/examples/advanced) has the [MARCO algorithm](https://github.com/CPMpy/cpmpy/blob/master/examples/advanced/marco_musmss_enumeration.py) that can enumerate all MSS/MUSes.) ### Correcting an UNSAT program @@ -142,6 +143,8 @@ sat_cons = mss(model.constraints) # x[0] or x[1], x[2] -> x[1], ~x[0] cons_to_remove = (mcs(model.constraints)) # x[0] ``` +More information about these tools can be found in [their API documentation](api/tools/explain.html). + ## Debugging a satisfiable model, that does not contain an expected solution We will ignore the (possible) objective function here and focus on the feasibility part. @@ -158,5 +161,5 @@ This one is most annoying... Double check the printing of the model for oddities Try generating an explanation sequence for the solution... this requires a satisfaction problem, so remove the objective function or add a constraint that constraints the objective function to the value attained by the impossible solution. -As to generating the explanation sequence, check out our advanced example on [stepwise OCUS explanations](https://github.com/CPMpy/cpmpy/blob/master/examples/advanced/ocus_explanations.py) +As to generating the explanation sequence, check out our advanced example on [stepwise OCUS explanations](https://github.com/CPMpy/cpmpy/blob/master/examples/advanced/ocus_explanations.py). diff --git a/docs/index.rst b/docs/index.rst index 3d17e5df1..b26280641 100644 --- a/docs/index.rst +++ b/docs/index.rst @@ -71,4 +71,5 @@ Are you a solver developer? We are keen to integrate solvers that have a python api/expressions api/model api/solvers - api/transformations \ No newline at end of file + api/transformations + api/tools \ No newline at end of file diff --git a/docs/modeling.md b/docs/modeling.md index d19f79608..4014eb7d4 100644 --- a/docs/modeling.md +++ b/docs/modeling.md @@ -4,7 +4,7 @@ This page explains and demonstrates how to use CPMpy to model and solve combinat ## Installation -Installation is available through the `pip` python package manager. This will also install and use `ortools` as default solver: +Installation is available through the `pip` Python package manager. This will also install and use `ortools` as default solver (see how to use alternative solvers [here](./modeling.html#selecting-a-solver)): ```commandline pip install cpmpy @@ -14,16 +14,18 @@ See [installation instructions](installation_instructions.html) for more details ## Using the library -To conveniently use CPMpy in your python project, include it as follows: +To conveniently use CPMpy in your Python project, include it as follows: ```python from cpmpy import * ``` -This will overload the built-in `any()`, `all()`, `min()`, `max()`, `sum()` functions, such that they create CPMpy expressions when used on decision variables (see below). This convenience comes at the cost of some overhead for all uses of these functions in your code. +This will overload the built-in `any()`, `all()`, `min()`, `max()`, `sum()` functions, such that they create CPMpy expressions when used on decision variables (see below). This convenience comes at the cost of some overhead for all uses of these functions in your code (even when no CPMpy decision variables are involved). -You can also import it as a package, which does not overload the python built-ins: +You can also import it as a package, which does not overload the Python built-ins: ```python import cpmpy as cp ``` +The build-in operators can now be accessed through the package (for example `cp.any()`) without overloading Python's defaults. + We will use the latter in this document. ## Decision variables @@ -35,17 +37,19 @@ CPMpy supports discrete decision variables, namely Boolean and integer decision ```python import cpmpy as cp -b = cp.boolvar(name="b") -x = cp.intvar(1,10, name="x") +b = cp.boolvar(name="b") # a Boolean decision variable +x = cp.intvar(1,10, name="x") # an Integer decision variable ``` -Decision variables have a **domain**, a set of allowed values. For Boolean variables this is implicitly the values 'False' and 'True'. For Integer decision variables, you have to specify the lower-bound and upper-bound (`1` and `10` respectively above). +Decision variables have a **domain**, a set of allowed values. For Boolean variables this is implicitly the values 'False' and 'True'. For Integer decision variables, you have to specify the lower-bound and upper-bound (`1` and `10` respectively in the above example). + +If you want a **sparse domain**, containing only a few values, you can either define a suitable lower/upper bound and then forbid specific values, e.g. `x != 3, x != 5, x != 7`; or you can use the shorthand *InDomain* global constraint: `InDomain(x, [1,2,4,6,8,9])`. -Decision variables have a **unique name**. You can set it yourself, otherwise a unique name will automatically be assigned to it. If you print `print(b, x)` decision variables, it will print the name. Did we already say the name must be unique? Many solvers use the name as unique identifier, and it is near-impossible to debug with non-uniquely named variables. +Decision variables have a **unique name**. You can set it yourself, otherwise a unique name will automatically be assigned. If you print decision variables (`print(b)` or `print(x)`), it will display the name. Did we already say the name must be unique? Many solvers use the name as unique identifier, and it is near-impossible to debug with non-uniquely named decision variables. -A solver will set the **value** of the decision variables for which it solved, if it can find a solution. You can retrieve it with `v.value()`. Variables are not tied to a solver, so you can use the same variable in multiple models and solvers. When a solve call finishes, it will overwrite the value of all its decision variables. +A solver will set the **value** of the decision variables for which it solved, if it can find a solution. You can retrieve it with `v.value()`. Variables are not tied to a solver, so you can use the same variable across multiple models and solvers. When a solve call finishes, it will overwrite the value of all its decision variables. Before solving, this value will be `None`. After solving it has either taken a boolean or integer value, or it is still None. For example when the solver didn't find any solution or when the decision variable was never used in the model, i.e. never appeared in a constraint or the objective function (a stale decision variable). -Finally, by providing a **shape** you automatically create a **numpy n-dimensional array** of decision variables. They automatically get their index appended to their name to ensure it is unique: +Finally, by providing a **shape** you automatically create a **numpy n-dimensional array** of decision variables. These variables automatically get their index appended to their name (the name is provided on the array-level) as to ensure its uniqueness: ```python import cpmpy as cp @@ -57,15 +61,19 @@ x = cp.intvar(1,10, shape=(2,2), name="x") print(x) # [[x[0,0] x[0,1]] # [x[1,0] x[1,1]]] ``` -You can also call `v.value()` on these n-dimensional arrays, which will return an n-dimensional **numpy** array of values. And you can do vectorized operations and comparisons, like in regular numpy. As we will see below, this is very convenient and avoids having to write out many loops. It also makes it compatible with many existing scientific python tools, including machine learning and visualisation libraries, so a lot less glue code to write. +Similar to individual decision variables, you can call `v.value()` on these n-dimensional arrays. This will return an n-dimensional **numpy** array of values, one value for each of the included decision variables. + +Since the arrays of decision variables are based on numpy, you can do **vectorized operations** and **comparisons** on them. As we will see below, this is very convenient and avoids having to write out many loops. It also makes it compatible with many existing scientific Python tools, including machine learning and visualisation libraries. A lot less glue code will need to be written! See [the API documentation on variables](api/expressions/variables.html) for more detailed information. +Note that decision variables are not tied to a model. You can use the same variable in different models; its value() will be the one of the last succesful solve call. + ## Creating a model -A **model** is a collection of constraints over decision variables, optionally with an objective function. It represents a problem for which a solution must be found, e.g. by a solver. A solution is an assignment to the decision variables, such that each of the constraints is satisfied. +A **model** is a collection of constraints over decision variables, optionally with an objective function. It represents a problem for which a solution must be found, e.g. through the use of a solver. A solution is an assignment of values to the decision variables, such that the values lie within their respective variables' domain and such that each of the constraints is satisfied. -In CPMpy, the `Model()` object is a simple container that stores a list of CPMpy expressions representing constraints. It can also store a CPMpy expression representing an objective function that must be minimized or maximized. Constraints are added in the constructor, or using the built-in `+=` addition operator that corresponds to calling the `__add__()` function. +In CPMpy, the `Model()` object is a simple container that stores a list of CPMpy expressions representing constraints. In the case of an optimisation problem, it can also store a CPMpy expression representing an objective function that must be minimized or maximized. Constraints are added in the constructor, or using the built-in `+=` addition operator that corresponds to calling the `__add__()` function. Here is an example, where we explain how to express constraints in the next section: @@ -92,6 +100,8 @@ print(m) # pretty printing of the model, very useful for debugging ``` The `Model()` object has a number of other helpful functions, such as `to_file()` to store the model and `copy()` for creating a copy. +See [the API documentation on models](api/model.html) for more detailed information. + ## Expressing constraints A constraint is a relation between decision variables that restricts what values these variables can take. @@ -158,10 +168,10 @@ import cpmpy as cp a,b,c = cp.boolvar(shape=3) m = cp.Model( - a.implies(b), - b.implies(a), - a.implies(~c), - (~c).implies(a) + a.implies(b), # a -> b + b.implies(a), # b -> a + a.implies(~c), # a -> (~c) + (~c).implies(a) # (~c) -> a ) ``` For reverse implication, you switch the arguments yourself; it is difficult to read reverse implications out loud anyway. And as before, always use brackets around subexpressions to avoid surprises! @@ -170,7 +180,7 @@ For reverse implication, you switch the arguments yourself; it is difficult to r ### Simple comparison constraints -We overload Pythons comparison operators: `==, !=, <, <=, >, >=`. Comparisons are allowed between any CPMpy expressions as well as Boolean and integer constants. +We overload Python's comparison operators: `==, !=, <, <=, >, >=`. Comparisons are allowed between any CPMpy expressions as well as Boolean and integer constants. On a technical note, we treat Booleans as a subclass of integer expressions. This means that a Boolean (output) expression can be used anywhere a numeric expression can be used, where `True` is treated as `1` and `False` as `0`. But the inverse is not true: integers can NOT be used with Boolean operators, even when you intialise their domain to (0,1) they are still not Boolean: @@ -222,9 +232,9 @@ m = cp.Model( ### Arithmetic constraints -We overload Python's built-in arithmetic operators `+,-,*,//,%`. These can be used to built arbitrarily nested numeric expressions, which can then be turned into a constraint by adding a comparison to it. +We overload Python's built-in arithmetic operators `+,-,*,//,%`. These can be used to build arbitrarily nested numeric expressions, which can then be turned into a constraint by adding a comparison to it. -We also overwrite the built-in functions `abs(),sum(),min(),max()` which can be used to created numeric expressions. Some examples: +We also overwrite the built-in functions `abs(),sum(),min(),max()` which can be used for creating numeric expressions. Some examples: ```python import cpmpy as cp @@ -254,7 +264,7 @@ m = cp.Model( ) ``` -Note that because of our overloading of `+,-,*,//` some numpy functions like `np.sum(some_array)` will also create a CPMpy expression when used on CPMpy decision variables. However, this is not guaranteed, and other functions like `np.max(some_array)` will not. To **avoid surprises**, you should hence always take care to call the CPMpy functions `cp.sum(), `cp.max()` etc. We did overload `some_cpm_array.sum()` and `.min()/.max()` (including the axis= argument), so these are safe to use. +Note that because of our overloading of `+,-,*,//` some numpy functions like `np.sum(some_array)` will also create a CPMpy expression when used on CPMpy decision variables. However, this is not guaranteed, and other functions like `np.max(some_array)` will not. To **avoid surprises**, you should hence always take care to call the CPMpy functions `cp.sum()`, `cp.max()` etc. We did overload `some_cpm_array.sum()` and `.min()`/`.max()` (including the axis= argument), so these are safe to use. ### Global constraints @@ -265,15 +275,17 @@ In constraint solving, a global constraint is a function that expresses a relati A good example is the `AllDifferent()` global constraint that ensures all its arguments have distinct values. `AllDifferent(x,y,z)` can be decomposed into `[x!=y, x!=z, y!=z]`. For AllDifferent, the decomposition consists of _n*(n-1)_ pairwise inequalities, which are simpler constraints that most solvers support. -However, a solver that has specialised datastructures for this constraint specifically does not need to create the decomposition. Furthermore, for AllDifferent solvers can implement specialised algorithms that can propagate strictly stronger than the decomposed constraints can. +However, a solver that has specialised datastructures for this constraint specifically does not need to create the decomposition. Furthermore, solvers can implement specialised algorithms that can propagate strictly stronger than the decomposed constraints can. #### Global constraints -A non-exhaustive list of global constraints that are available in CPMpy is: `Xor(), AllDifferent(), AllDifferentExcept0(), Table(), Circuit(), Cumulative(), GlobalCardinalityCount()`. +Many global constraints are available in CPMpy. Some include `Xor(), AllDifferent(), AllDifferentExcept0(), Table(), Circuit(), Cumulative(), GlobalCardinalityCount()`. -For their meaning and more information on how to define your own global constraints, see [the API documentation on global constraints](api/expressions/globalconstraints.html). Global constraints can also be reified (e.g. used in an implication or equality constraint). +For a complete list of global constraints, their meaning and more information on how to define your own, see [the API documentation on global constraints](api/expressions/globalconstraints.html). + +Global constraints can also be reified (e.g. used in an implication or equality constraint). CPMpy will automatically decompose them if needed. If you want to see the decomposition yourself, you can call the `decompose()` function on them. @@ -283,24 +295,24 @@ x = cp.intvar(1,4, shape=4, name="x") b = cp.boolvar() cp.Model( cp.AllDifferent(x), - cp.AllDifferent(x).decompose(), # equivalent: [(x[0]) != (x[1]), (x[0]) != (x[2]), ... + cp.AllDifferent(x).decompose()[0], # equivalent: [(x[0]) != (x[1]), (x[0]) != (x[2]), ... b.implies(cp.AllDifferent(x)), cp.Xor(b, cp.AllDifferent(x)), # etc... ) ``` -`decompose()` returns two arguments, one that represents the constraints and an optional one that defines any new variables needed. This is technical, but important to make negation work, if you want to know more check the [the API documentation on global constraints](api/expressions/globalconstraints.html). +`decompose()` returns two arguments, one that represents the constraints and another that defines any newly created decision variables during the decomposition process. This is technical, but important to make negation work, if you want to know more check the [the API documentation on global constraints](api/expressions/globalconstraints.html). -#### Numeric global constraints +#### Global functions Coming back to the Python-builtin functions `min(),max(),abs()`, these are a bit special because they have a numeric return type. In fact, constraint solvers typically implement a global constraint `MinimumEq(args, var)` that represents `min(args) == var`, so it combines a numeric function with a comparison, where it will ensure that the bounds of the expressions on both sides satisfy the comparison relation. However, CPMpy also wishes to support the expressions `min(xs) > v` as well as `v + min(xs) != 4` and other nested expressions. -In CPMpy we do this by instantiating min/max/abs as **numeric global constraints**. E.g. `min([x,y,z])` becomes `Minimum([x,y,z])` which inherits from `GlobalFunction` because it has a numeric return type. Our library will transform the constraint model, including arbitrarly nested expressions, such that the numeric global constraint is used in a comparison with a variable. Then, the solver will either support it, or we will call `decompose_comparison()` on the numeric global constraint, which will decompose e.g. `min(xs) == v`. +In CPMpy we do this by instantiating min/max/abs as **global functions**. E.g. `min([x,y,z])` becomes `Minimum([x,y,z])` which inherits from `GlobalFunction` because it has a numeric return type. Our library will transform the constraint model, including arbitrarly nested expressions, such that the global function is used within a comparison with a variable. Then, the solver will either support it, or we will call `decompose_comparison()` ([link](api/expressions/globalfunctions.html#cpmpy.expressions.globalfunctions.Abs.decompose_comparison)) on the global function. A non-exhaustive list of **numeric global constraints** that are available in CPMpy is: `Minimum(), Maximum(), Count(), Element()`. -For their meaning and more information on how to define your own global constraints, see [the API documentation on global functions](api/expressions/globalfunctions.html). +For their meaning and more information on how to define your own global functions, see [the API documentation on global functions](api/expressions/globalfunctions.html). ```python import cpmpy as cp @@ -318,7 +330,7 @@ print(s.transform(cp.min(x) + cp.max(x) - 5 > 2*cp.Count(x, 2))) ``` -#### The Element numeric global constraint +#### The Element global function The `Element(Arr,Idx)` global function enforces that the result equals `Arr[Idx]` with `Arr` an array of constants or variables (the first argument) and `Idx` an integer decision variable, representing the index into the array. @@ -337,7 +349,7 @@ print(f"arr: {arr.value()}, idx: {idx.value()}, val: {arr[idx].value()}") # example output -- arr: [2 1 3 4], idx: 0, val: 2 ``` -The `arr[idx]` works because `arr` is a CPMpy `NDVarArray()` and we overloaded the `__getitem__()` python function. It even supports multi-dimensional access, e.g. `arr[idx1,idx2]`. +The `arr[idx]` works because `arr` is a CPMpy `NDVarArray()` and we overloaded the `__getitem__()` Python function. It even supports multi-dimensional access, e.g. `arr[idx1,idx2]`. This does not work on NumPy arrays though, as they don't know CPMpy. So you have to **wrap the array** in our `cpm_array()` or call `Element()` directly: @@ -367,7 +379,7 @@ print(f"arr: {arr.value()}, idx: {idx.value()}, val: {arr[idx].value()}") ## Objective functions -If a model has no objective function specified, then it is a satisfaction problem: the goal is to find out whether a solution, any solution, exists. When an objective function is added, this function needs to be minimized or maximized. +If a model has no objective function specified, then it represents a satisfaction problem: the goal is to find out whether a solution, any solution, exists. When an objective function is added, this function needs to be minimized or maximized. Any CPMpy expression can be added as objective function. Solvers are especially good in optimizing linear functions or the minimum/maximum of a set of expressions. Other (non-linear) expressions are supported too, just give it a try. @@ -429,21 +441,29 @@ n = m.solveAll() print("Nr of solutions:", n) # Nr of solutions: 6 ``` -When using `solveAll()`, a solver will use an optimized native implementation behind the scenes when that exists. +When using `solveAll()`, a solver will use an optimized native implementation behind the scenes when that exists. Otherwise it will be emulated with an iterative approach, resulting in a performance impact. -It has a `display=...` argument that can be used to display expressions or as a callback, as well as the `solution_limit=...` argument to set a solution limit. It also accepts any named argument, like `time_limit=...`, that the underlying solver accepts. +It has a `display=...` argument that can be used to display expressions (provide a list of expressions to be evaluated) or as a more generic callback (provide a function), both to be evaluated for each found solution. +It has a `solution_limit=...` argument to set a limit on the number of solutions to solve for. +It also accepts any named argument, like `time_limit=...`, that the underlying solver accepts. For more information about the available arguments, look at [the solver API documentation](api/solvers.html) for the solver in question. ```python +# Using list of expressions n = m.solveAll(display=[x,cp.sum(x)], solution_limit=3) # [array([1, 0]), 1] # [array([2, 0]), 2] # [array([3, 0]), 3] + +# Using callback function +n = m.solveAll(display=lambda: print([x.value(), cp.sum(x).value()]), solution_limit=3) +# ... ``` +(Note that the Exact solver, unlike other solvers, takes most of its arguments at construction time.) There is much more to say on enumerating solutions and the use of callbacks or blocking clauses. See the [the detailed documentation on finding multiple solutions](multiple_solutions.html). ## Debugging a model -If the solver is complaining about your model, then a good place to start debugging is to **print** the model you have created, or the individual constraints. If they look fine (e.g. no integers, or shorter or longer expressions then what you intended) and you don't know which constraint specifically is causing the error, then you can feed the constraints incrementally to the solver you are using: +If the solver is complaining about your model, then a good place to start debugging is to **print** the model (or individual constraints) that you have created. If they look fine (e.g. no integers, or shorter or longer expressions then what you intended) and you don't know which constraint specifically is causing the error, then you can feed the constraints incrementally to the solver you are using: ```python import cpmpy as cp @@ -456,7 +476,7 @@ m = cp.Model(cons) # any model created # e.g. if you wrote `all(x)` instead of `cp.all(x)` it will contain 'True' instead of the conjunction print(m) -s = cp.SolverLookup.get("ortools") +s = cp.SolverLookup.get("ortools") # will be explained later # feed the constraints one-by-one for c in m.constraints: s += c # add the constraints incrementally until you hit the error @@ -466,8 +486,9 @@ If that is not sufficient or you want to debug an unexpected (non)solution, have ## Selecting a solver -The default solver is OR-Tools CP-SAT, an award winning constraint solver. But CPMpy supports multiple other solvers: a MIP solver (gurobi), SAT solvers (those in PySAT), the Z3 SMT solver, even a knowledge compiler (PySDD) and any CP solver supported by the text-based MiniZinc language. +The default solver is [OR-Tools CP-SAT](https://developers.google.com/optimization), an award winning constraint solver. But CPMpy supports multiple other solvers: a MIP solver (gurobi), SAT solvers (those in PySAT), the Z3 SMT solver, even a knowledge compiler (PySDD) and any CP solver supported by the text-based MiniZinc language. +The list of supported solver interfaces can be found in [the API documentation](api/solvers.html). See the full list of solvers known by CPMpy with: ```python @@ -475,7 +496,7 @@ import cpmpy as cp cp.SolverLookup.solvernames() ``` -On my system, with pysat and minizinc installed, this gives `['ortools', 'minizinc', 'minizinc:chuffed', 'minizinc:coin-bc', ..., 'pysat:minicard', 'pysat:minisat22', 'pysat:minisat-gh'] +On a system with pysat and minizinc installed, this for example gives `['ortools', 'minizinc', 'minizinc:chuffed', 'minizinc:coin-bc', ..., 'pysat:minicard', 'pysat:minisat22', 'pysat:minisat-gh']` You can specify a solvername when calling `solve()` on a model: @@ -487,11 +508,11 @@ m = cp.Model(cp.sum(x) <= 5) m.solve(solver="minizinc:chuffed") ``` -Note that for solvers other than "ortools", you will need to **install additional package(s)**. You can check if a solver, e.g. "minizinc", is supported by calling `cp.SolverLookup.get("gurobi")` and it will raise a helpful error if it is not yet installed on your system. See [the API documentation](api/solvers.html) of the solver for detailed installation instructions. +Note that for solvers other than "ortools", you will need to **install additional package(s)**. You can check if a solver, e.g. "gurobi", is supported by calling `cp.SolverLookup.get("gurobi")` and it will raise a helpful error if it is not yet installed on your system. See [the API documentation](api/solvers.html) of the solver for detailed installation instructions. ## Model versus solver interface -A `Model()` is a **lazy container**. It simply stores the constraints. Only when `solve()` is called will it instantiate a solver, and send the entire model to it at once. So `m.solve("ortools")` is equivalent to: +A `Model()` is a **lazy container**. It simply stores the constraints. Only when `solve()` is called will it instantiate a solver instance, and send the entire model to it at once. So `m.solve("ortools")` is equivalent to: ```python s = SolverLookup.get("ortools", m) s.solve() @@ -517,7 +538,7 @@ On a technical note, remark that a solver object does not modify the Model objec ## Setting solver parameters Now lets use our solver-specific powers. -For example, with `m` a CPMpy Model(), you can do the following to make or-tools use 8 parallel cores and print search progress: +For example, with `m` a CPMpy Model(), you can do the following to make OR-Tools use 8 parallel cores and print search progress: ```python import cpmpy as cp @@ -526,7 +547,7 @@ s = cp.SolverLookup.get("ortools", m) s.solve(num_search_workers=8, log_search_progress=True) ``` -Modern CP-solvers support a variety of hyperparameters. (See the full list of [OR-tools parameters](https://github.com/google/or-tools/blob/stable/ortools/sat/sat_parameters.proto) for example). +Modern CP-solvers support a variety of hyperparameters. (See the full list of [OR-Tools parameters](https://github.com/google/or-tools/blob/stable/ortools/sat/sat_parameters.proto) for example). Using the solver interface, any parameter that the solver supports can be passed using the `.solve()` call. These parameters will be posted to the native solver before solving the model. @@ -541,7 +562,7 @@ See [the API documentation of the solvers](api/solvers.html) for information and ## Accessing the underlying solver object After solving, we can access the underlying solver object to retrieve some information about the solve. -For example in ortools we can find the number of search branches like this, expanding our previous example. +For example in OR-Tools we can find the number of search branches like this (expanding on our previous example): ```python import cpmpy as cp x = cp.intvar(0,10, shape=3) @@ -552,7 +573,7 @@ s.solve() print(s.ort_solver.NumBranches()) ``` -Other solvers, like Minizinc might have other native objects stored. +Other solvers, like Minizinc, might have other native objects stored. You can see which solver native objects are initialized for each solver in [the API documentation](api/solvers.html) of the solver. We can access the solver statistics from the mzn_result object like this: @@ -564,7 +585,7 @@ s = cp.SolverLookup.get("minizinc") s += (cp.sum(x) <= 5) s.solve() print(s.mzn_result.statistics) -print(s.mzn_result.statistics['nodes']) #if we are only interested in the nb of search nodes +print(s.mzn_result.statistics['nodes']) # if we are only interested in the nb of search nodes ``` ## Incremental solving @@ -573,18 +594,18 @@ It is important to realize that a CPMpy solver interface is _eager_. That means This has two potential benefits for incremental solving, whereby you add more constraints and variables inbetween solve calls: 1) CPMpy only translates and posts each constraint once, even if the model is solved multiple times; and - 2) if the solver itself is incremental then it can reuse any information from call to call, as the state of the native solver object is kept between solver calls and can therefore rely on information derived during a previous `solve` call. + 2) if the solver itself is incremental then it can reuse any information from call to call, as the state of the native solver object is kept between solver calls and can therefore rely on information derived during a previous `.solve()` call. ```python -gs = SolverLookup.get("gurobi") +s = SolverLookup.get("gurobi") -gs += sum(ivar) <= 5 -gs.solve() +s += sum(ivar) <= 5 +s.solve() -gs += sum(ivar) == 3 +s += sum(ivar) == 3 # the underlying gurobi instance is reused, only the new constraint is added to it. # gurobi is an incremental solver and will look for solutions starting from the previous one. -gs.solve() +s.solve() ``` _Technical note_: OR-Tools its model representation is incremental but its solving itself is not (yet?). Gurobi and the PySAT solvers are fully incremental, as is Z3. The text-based MiniZinc language is not incremental. @@ -593,7 +614,7 @@ _Technical note_: OR-Tools its model representation is incremental but its solvi SAT and CP-SAT solvers oftentimes support solving under assumptions, which is also supported by their CPMpy interface. Assumption variables are usefull for incremental solving when you want to activate/deactivate different subsets of constraints without copying (parts of) the model or removing constraints and re-solving. By relying on the solver interface directly as in the previous section, the state of the solver is kept in between solve-calls. -Many explanation-generation algorithms (see `cpmpy.tools.explain`) make use of this feature to speed up the solving. +Many explanation-generation algorithms ([see](api/tools/explain.html) `cpmpy.tools.explain`) make use of this feature to speed up the solving. ```python import cpmpy as cp @@ -608,7 +629,7 @@ cp.Model([c1,c2,c3]).solve() # Will be UNSAT s = cp.SolverLookup.get("exact") # OR-tools, PySAT and Exact support solving under assumptions assump = cp.boolvar(shape=3, name="assump") -s += assump.implies([c1,c2,c3]) +s += assump.implies([c1,c2,c3]) # assump[0] -> c1, assump[1] -> c2, assump[2] -> c3 # Underlying solver state will be kept inbetween solve-calls s.solve(assumptions=assump[0,1]) # Will be SAT @@ -620,15 +641,15 @@ s.solve(assumptions=assump[1,2]) # Will be SAT ## Using solver-specific CPMpy features -We sometimes add solver-specific functions to the CPMpy interface, for convenient access. Two examples of this are `solution_hint()` and `get_core()` which is supported by the OR-Tools and PySAT solvers and interfaces. Other solvers may work differently and not have these concepts. +We sometimes add solver-specific functions to the CPMpy interface, for convenient access. Two examples of this are `solution_hint()` and `get_core()` which is supported by the OR-Tools, PySAT and Exact solvers and interfaces. Other solvers may work differently and not have these concepts. -`solution_hint()` tells the solver that it could use these variable-values first during search, e.g. typically from a previous solution: +`solution_hint()` tells the solver that it could start from these value-to-variable assignments during search, e.g. start from a previous solution: ```python import cpmpy as cp x = cp.intvar(0,10, shape=3) s = cp.SolverLookup.get("ortools") s += cp.sum(x) <= 5 -# we are operating on a ortools' interface here +# we are operating on an ortools' interface here s.solution_hint(x, [1,2,3]) s.solve() print(x.value()) @@ -679,7 +700,10 @@ s += cp.DirectConstraint("AddAutomaton", (trans_vars, 0, [2], trans_tabl), A minimal example of the DirectConstraint for every supported solver is [in the test suite](https://github.com/CPMpy/cpmpy/tree/master/tests/test_direct.py). -The `DirectConstraint` is a very powerful primitive to get the most out of specific solvers. See the following examples: [nonogram_ortools.ipynb](https://github.com/CPMpy/cpmpy/tree/master/examples/nonogram_ortools.ipynb) which uses a helper function that generates automatons with DirectConstraints; [vrp_ortools.py](https://github.com/CPMpy/cpmpy/tree/master/examples/vrp_ortools.ipynb) demonstrating ortools' newly introduced multi-circuit global constraint through DirectConstraint; and [pctsp_ortools.py](https://github.com/CPMpy/cpmpy/tree/master/examples/pctsp_ortools.ipynb) that uses a DirectConstraint to use OR-Tools circuit to post a sub-circuit constraint as needed for this price-collecting TSP variant. +The `DirectConstraint` is a very powerful primitive to get the most out of specific solvers. See the following examples: +- [nonogram_ortools.ipynb](https://github.com/CPMpy/cpmpy/tree/master/examples/nonogram_ortools.ipynb) which uses a helper function that generates automatons with DirectConstraints; +- [vrp_ortools.py](https://github.com/CPMpy/cpmpy/blob/master/examples/vrp_ortools.py) demonstrating OR-Tools' newly introduced multi-circuit global constraint through DirectConstraint; and +- [pctsp_ortools.py](https://github.com/CPMpy/cpmpy/blob/master/examples/pctsp_ortools.py) that uses a DirectConstraint to use OR-Tools circuit to post a sub-circuit constraint as needed for this price-collecting TSP variant. ### Directly accessing the underlying solver @@ -705,7 +729,7 @@ Because CPMpy offers programmatic access to the solver API, hyperparameter searc ### Built-in tuners -The tools directory contains a utility to efficiently search through the hyperparameter space defined by the solvers `tunable_params`. +The tools directory contains a [utility](https://github.com/CPMpy/cpmpy/blob/master/cpmpy/tools/tune_solver.py) to efficiently search through the hyperparameter space defined by the solvers' `tunable_params`. Solver interfaces not providing the set of tunable parameters can still be tuned by using this utility and providing the parameter (values) yourself. diff --git a/docs/multiple_solutions.md b/docs/multiple_solutions.md index 06d8b2b09..b0234df98 100644 --- a/docs/multiple_solutions.md +++ b/docs/multiple_solutions.md @@ -6,8 +6,8 @@ When using `solveAll()`, a solver will use an optimized native implementation be It has two special named optional arguments: - * `display=...`: accepts a CPMpy expression, a list of CPMpy expressions or a callback function that will be called on every solution found (default: None) - * `solution_limit=...`: stop after this many solutions (default: None) + - `display=...` : accepts a CPMpy expression, a list of CPMpy expressions or a callback function that will be called on every solution found (default: None) + - `solution_limit=...` : stop after this many solutions (default: None) It also accepts named argument `time_limit=...` and any other keyword argument is passed on to the solver just like `solve()` does. @@ -93,7 +93,7 @@ while s.solve(): s += ~all(x == x.value()) ``` -In case of multiple variables you should put them in one long python-native list, as such: +In case of multiple variables you should put them in one long Python-native list, as such: ```python x = intvar(0,3, shape=2) b = boolvar() @@ -111,7 +111,7 @@ while s.solve(): ## Diverse solution search A better, more complex example of repeated solving is when searching for diverse solutions. -The goal is to iteratively find solutions that are as diverse as possible with the previous solutions. Many definitions of diversity between solutions exist. We can for example measure the difference between two solutions with the Hamming distance (comparing the number of different values) or the Euclidian distance (compare the absolute difference in value for the variables). +The goal is to iteratively find solutions that are as diverse as possible with the previous solutions. Many definitions of diversity between solutions exist. We can for example measure the difference between two solutions with the [Hamming distance](https://en.wikipedia.org/wiki/Hamming_distance#:~:text=In%20information%20theory%2C%20the%20Hamming,the%20corresponding%20symbols%20are%20different.) (comparing the number of different values) or the [Euclidian distance](https://en.wikipedia.org/wiki/Euclidean_distance) (compare the absolute difference in value for the variables). Here is the example code for enumerating K diverse solutions with Hamming distance, which overwrites the objective function in each iteration: @@ -130,7 +130,7 @@ while len(store) < 3 and s.solve(): s.maximize(sum([sum(x != sol) for sol in store])) ``` -As a fun fact, observe how `x != sol` works, even though one is a vector of Boolean variables and sol is Numpy array. However, both have the same length, so this is automatically translated into a pairwise comparison constraint by CPMpy. These numpy-style vectorized operations mean we have to write much less loops while modelling. +As a fun fact, observe how `x != sol` works, even though one is a vector of Boolean variables and sol is Numpy array. However, both have the same length, so this is automatically translated into a pairwise comparison constraint by CPMpy. These numpy-style vectorized operations mean we have to write much less loops while modeling. To use the Euclidian distance, only the last line needs to be changed. We again use vectorized operations and obtain succinct models. The creation of intermediate variables (with appropriate domains) is done by CPMpy behind the scenes. @@ -157,3 +157,4 @@ cb = OrtSolutionPrinter() s.solve(enumerate_all_solutions=True, solution_callback=cb) print("Nr of solutions:",cb.solution_count()) ``` +Have a look at `OrtSolutionPrinter`'s [implementation](https://github.com/CPMpy/cpmpy/blob/master/cpmpy/solvers/ortools.py#L640). diff --git a/docs/requirements.txt b/docs/requirements.txt index bfbbd40c7..2b3a4c6b3 100644 --- a/docs/requirements.txt +++ b/docs/requirements.txt @@ -1,7 +1,6 @@ numpy>=1.17.0 m2r2>=0.2.7 sphinx-automodapi -sphinx_rtd_theme -sphinx==5.3.0 -sphinx_rtd_theme==1.1.1 -readthedocs-sphinx-search==0.1.1 \ No newline at end of file +sphinx>=5.3.0 # cannot be updated due to m2r2 no longer being maintained -> switch to MyST? +sphinx_rtd_theme>=2.0.0 +readthedocs-sphinx-search>=0.3.2 \ No newline at end of file diff --git a/docs/solvers.md b/docs/solvers.md index 11ecd309a..76035b576 100644 --- a/docs/solvers.md +++ b/docs/solvers.md @@ -2,7 +2,7 @@ CPMpy can be used as a declarative modeling language: you create a `Model()`, add constraints and call `solve()` on it. -The default solver is ortools CP-SAT, an award winning constraint solver. But CPMpy supports multiple other solvers: a MIP solver (gurobi), SAT solvers (those in PySAT), the Z3 SMT solver, a conflict-driven cutting-planes solver (Exact), even a knowledge compiler (PySDD) and any CP solver supported by the text-based MiniZinc language. +The default solver is OR-Tools CP-SAT, an award winning constraint solver. But CPMpy supports multiple other solvers: a MIP solver (gurobi), SAT solvers (those in PySAT), the Z3 SMT solver, a conflict-driven cutting-planes solver (Exact), even a knowledge compiler (PySDD) and any CP solver supported by the text-based MiniZinc language. See the list of solvers known by CPMpy with: diff --git a/docs/unsat_core_extraction.md b/docs/unsat_core_extraction.md index fd8228097..aadd2b034 100644 --- a/docs/unsat_core_extraction.md +++ b/docs/unsat_core_extraction.md @@ -4,9 +4,9 @@ When a model is unsatisfiable, it can be desirable to get a better idea of which In the SATisfiability literature, the Boolean variables of interests are called _assumption_ variables and the solver will assume they are true. The subset of these variables that, when true, make the model unsatisfiable is called an unsatisfiable _core_. -Lazy Clause Generation solvers, like or-tools, are built on SAT solvers and hence can inherit the ability to define assumption variables and extract an unsatisfiable core. +Lazy Clause Generation solvers, like OR-Tools, are built on SAT solvers and hence can inherit the ability to define assumption variables and extract an unsatisfiable core. -Since version 8.2, or-tools supports declaring assumption variables, and extracting an unsat core. We also implement this functionality in CPMpy, using PySAT-like `s.solve(assumptions=[...])` and `s.get_core()`: +Since version 8.2, OR-Tools supports declaring assumption variables, and extracting an unsat core. We also implement this functionality in CPMpy, using PySAT-like `s.solve(assumptions=[...])` and `s.get_core()`: ```python from cpmpy import * @@ -39,8 +39,10 @@ from cpmpy.tools import mus print(mus(m.constraints)) ``` -We welcome any additional examples that use CPMpy in this way!! Here is one example: the [MARCO algorithm for enumerating all MUS/MSSes](http://github.com/tias/cppy/tree/master/examples/advanced/marco_musmss_enumeration.py). Here is another: a [stepwise explanation algorithm](https://github.com/CPMpy/cpmpy/blob/master/examples/advanced/ocus_explanations.py) for SAT problems (implicit hitting-set based) +We welcome any additional examples that use CPMpy in this way!! Here is one example: the [MARCO algorithm for enumerating all MUS/MSSes](http://github.com/tias/cppy/tree/master/examples/advanced/marco_musmss_enumeration.py). Here is another: a [stepwise explanation algorithm](https://github.com/CPMpy/cpmpy/blob/master/examples/advanced/ocus_explanations.py) for SAT problems (implicit hitting-set based). -One OR-TOOLS specific caveat is that this particular (default) solver its Python interface is by design _stateless_. That means that, unlike in PySAT, calling `s.solve(assumptions=bv)` twice for a different `bv` array does NOT REUSE anything from the previous run: no warm-starting, no learnt clauses that are kept, no incrementality, so there will be some pre-processing overhead. If you know of another CP solver with a (Python) assumption interface that is incremental, let us know!! +More information on how to use these tools can be found in [the tools API documentation](api/tools.html) -A final-final note is that you can manually warm-start or-tools with a previously found solution with s.solution\_hint(); see also the MARCO code linked above. +One OR-Tools specific caveat is that this particular (default) solver its Python interface is by design _stateless_. That means that, unlike in PySAT, calling `s.solve(assumptions=bv)` twice for a different `bv` array does NOT REUSE anything from the previous run: no warm-starting, no learnt clauses that are kept, no incrementality, so there will be some pre-processing overhead. If you know of another CP solver with a (Python) assumption interface that is incremental, let us know!! + +A final-final note is that you can manually warm-start OR-Tools with a previously found solution through `s.solution_hint()`; see also the MARCO code linked above. diff --git a/examples/advanced/exact_maximal_propagate.py b/examples/advanced/exact_maximal_propagate.py new file mode 100644 index 000000000..66ff60c6e --- /dev/null +++ b/examples/advanced/exact_maximal_propagate.py @@ -0,0 +1,81 @@ +""" +Example showing direct solver access to Exact solver object. + +Solver-native implementation of the Maximal-Propagation algorithm as found in cpmpy.tools.maximal_propagate +Can be used for finding the maximal consequence of a set of constraints under assumptions. +""" + +from cpmpy import * +from cpmpy.expressions.variables import NegBoolView +from cpmpy.solvers import CPM_exact +from cpmpy.expressions.utils import is_any_list, is_num, flatlist + + +class PropagationSolver(CPM_exact): + + def __init__(self, cpm_model=None, subsolver=None): + # first need to set the encoding before adding constraints and variables + super(PropagationSolver, self).__init__(cpm_model=None, subsolver=subsolver) + self.encoding = "onehot" # required for use of pruneDomains function + + if cpm_model is not None: + # post all constraints at once, implemented in __add__() + self += cpm_model.constraints + + # post objective + if cpm_model.objective_ is not None: + if cpm_model.objective_is_min: + self.minimize(cpm_model.objective_) + else: + self.maximize(cpm_model.objective_) + + def maximal_propagate(self, assumptions=[], vars=None): + """ + Wrapper for the `pruneDomains` method of Exact + Automatically converts CPMpy variables to Exact variables and supports assumptions + :param: assumptions: list of Boolean variables (or their negation) assumed to be True + :param: vars: list of variables to propagate + """ + if self.solver_is_initialized is False: + self.solve() # this initializes solver and necessary datastructures + + self.xct_solver.clearAssumptions() # clear any old assumptions + + assump_vals = [int(not isinstance(v, NegBoolView)) for v in assumptions] + assump_vars = [self.solver_var(v._bv if isinstance(v, NegBoolView) else v) for v in assumptions] + for var, val in zip(assump_vars, assump_vals): + self.xct_solver.setAssumption(var, [val]) + + if vars is None: + vars = flatlist(self.user_vars) + + domains = self.xct_solver.pruneDomains(flatlist(self.solver_vars(vars))) + + return {var : {int(v) for v in dom} for var, dom in zip(vars, domains)} + + + +if __name__ == "__main__": + import cpmpy as cp + from cpmpy.tools.explain.utils import make_assump_model + + x = cp.intvar(1, 5, shape=5, name="x") + + c1 = cp.AllDifferent(x) + c2 = x[0] == cp.min(x) + c3 = x[-1] == 1 + + model, cons, assump = make_assump_model([c1,c2,c3]) + + propsolver = PropagationSolver(model) + + print(f"Propagating constraints {c1} and {c2}:") + possible_vals = propsolver.maximal_propagate(assumptions=assump[:2]) + for var, values in possible_vals.items(): + print(f"{var}: {sorted(values)}") + + print("\n\n") + print(f"Propagating constraints {c1}, {c2} and {c3}:") + possible_vals = propsolver.maximal_propagate(assumptions=assump) + for var, values in possible_vals.items(): + print(f"{var}: {sorted(values)}") diff --git a/examples/advanced/maximal_propagate.py b/examples/advanced/maximal_propagate.py deleted file mode 100644 index 1da9d2fe8..000000000 --- a/examples/advanced/maximal_propagate.py +++ /dev/null @@ -1,20 +0,0 @@ -""" -Example showing repeated solving for deriving the maximal consequence of set of constraints. - -Iteratively finds a new solutions by forcing at least one variable to have a value not seen in any other solution. -""" - -from cpmpy import * -from cpmpy.tools.maximal_propagate import maximal_propagate - - -if __name__ == "__main__": - x,y,z = [intvar(lb=0,ub=5, name=n) for n in "xyz"] - - model = Model([x < y, y != 3]) - print("Propagating constraints in model:", model, sep="\n", end="\n\n") - - possible_vals = maximal_propagate(model.constraints) - - for var, values in possible_vals.items(): - print(f"{var}: {sorted(values)}") diff --git a/requirements.txt b/requirements.txt index 5b2578c91..b8e491892 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,2 +1,2 @@ -numpy -ortools >= 9.6 +numpy < 2.0 +ortools >= 9.9 diff --git a/tests/test_constraints.py b/tests/test_constraints.py index 895d675bf..bdcc3fa0c 100644 --- a/tests/test_constraints.py +++ b/tests/test_constraints.py @@ -1,3 +1,6 @@ +import inspect + +import cpmpy from cpmpy import Model, SolverLookup, BoolVal from cpmpy.expressions.globalconstraints import * from cpmpy.expressions.globalfunctions import * @@ -8,15 +11,31 @@ # make sure that `SolverLookup.get(solver)` works # also add exclusions to the 3 EXCLUDE_* below as needed SOLVERNAMES = [name for name, solver in SolverLookup.base_solvers() if solver.supported()] +ALL_SOLS = False # test wheter all solutions returned by the solver satisfy the constraint # Exclude some global constraints for solvers -# Can be used when .value() method is not implemented/contains bugs -EXCLUDE_GLOBAL = {"ortools": {}, - "gurobi": {}, - "minizinc": {"circuit"}, - "pysat": {"circuit", "element","min","max","allequal","alldifferent","cumulative"}, - "pysdd": {"circuit", "element","min","max","allequal","alldifferent","cumulative",'xor'}, - "exact": {}, +NUM_GLOBAL = { + "AllEqual", "AllDifferent", "AllDifferentExcept0", + "AllDifferentExceptN", "AllEqualExceptN", + "GlobalCardinalityCount", "InDomain", "Inverse", "Table", 'NegativeTable', "Circuit", + "Increasing", "IncreasingStrict", "Decreasing", "DecreasingStrict", + "Precedence", "Cumulative", "NoOverlap", + "LexLess", "LexLessEq", "LexChainLess", "LexChainLessEq", + # also global functions + "Abs", "Element", "Minimum", "Maximum", "Count", "Among", "NValue", "NValueExcept" +} + +# Solvers not supporting arithmetic constraints +SAT_SOLVERS = {"pysat", "pysdd"} + +EXCLUDE_GLOBAL = {"pysat": NUM_GLOBAL, + "pysdd": NUM_GLOBAL | {"Xor"}, + "z3": {"Inverse"}, + "choco": {"Inverse"}, + "ortools":{"Inverse"}, + "exact": {"Inverse"}, + "minizinc": {"IncreasingStrict"}, # bug #813 reported on libminizinc + "gcs": {} } # Exclude certain operators for solvers. @@ -27,16 +46,6 @@ "exact": {"mod","pow","div","mul"}, } -# Some solvers only support a subset of operators in imply-constraints -# This subset can differ between left and right hand side of the implication -EXCLUDE_IMPL = {"ortools": {}, - "minizinc": {}, - "z3": {}, - "pysat": {}, - "pysdd": {}, - "exact": {"mod","pow","div","mul"}, - } - # Variables to use in the rest of the test script NUM_ARGS = [intvar(-3, 5, name=n) for n in "xyz"] # Numerical variables NN_VAR = intvar(0, 10, name="n_neg") # Non-negative variable, needed in power functions @@ -59,8 +68,7 @@ def numexprs(solver): Numexpr: - Operator (non-Boolean) with all args Var/constant (examples: +,*,/,mod,wsum) (CPMpy class 'Operator', not is_bool()) - - Global constraint (non-Boolean) (examples: Max,Min,Element) - (CPMpy class 'GlobalConstraint', not is_bool())) + - Global functions (examples: Max,Min,Element) (CPMpy class 'GlobalFunction') """ names = [(name, arity) for name, (arity, is_bool) in Operator.allowed.items() if not is_bool] if solver in EXCLUDE_OPERATORS: @@ -79,6 +87,35 @@ def numexprs(solver): yield Operator(name, operator_args) + # boolexprs are also numeric + for expr in bool_exprs(solver): + yield expr + + # also global functions + classes = inspect.getmembers(cpmpy.expressions.globalfunctions, inspect.isclass) + classes = [(name, cls) for name, cls in classes if issubclass(cls, GlobalFunction) and name != "GlobalFunction"] + classes = [(name, cls) for name, cls in classes if name not in EXCLUDE_GLOBAL.get(solver, {})] + + for name, cls in classes: + if name == "Abs": + expr = cls(NUM_ARGS[0]) + elif name == "Count": + expr = cls(NUM_ARGS, NUM_VAR) + elif name == "Element": + expr = cls(NUM_ARGS, POS_VAR) + elif name == "NValueExcept": + expr = cls(NUM_ARGS, 3) + elif name == "Among": + expr = cls(NUM_ARGS, [1,2]) + else: + expr = cls(NUM_ARGS) + + if solver in EXCLUDE_GLOBAL and expr.name in EXCLUDE_GLOBAL[solver]: + continue + else: + yield expr + + # Generate all possible comparison constraints def comp_constraints(solver): @@ -91,28 +128,19 @@ def comp_constraints(solver): - Numeric inequality (>=,>,<,<=): Numexpr >=< Var (CPMpy class 'Comparison') """ for comp_name in Comparison.allowed: + for numexpr in numexprs(solver): - for rhs in [NUM_VAR, BOOL_VAR, 1, BoolVal(True)]: + # numeric vs bool/num var/val (incl global func) + lb, ub = get_bounds(numexpr) + for rhs in [NUM_VAR, BOOL_VAR, BoolVal(True), 1]: + if solver in SAT_SOLVERS and not is_num(rhs): + continue + if comp_name == ">" and ub <= get_bounds(rhs)[1]: + continue + if comp_name == "<" and lb >= get_bounds(rhs)[0]: + continue yield Comparison(comp_name, numexpr, rhs) - for comp_name in Comparison.allowed: - for glob_expr in global_constraints(solver): - if not glob_expr.is_bool(): - for rhs in [NUM_VAR, BOOL_VAR, 1, BoolVal(True)]: - yield Comparison(comp_name, glob_expr, rhs) - - if solver == "z3": - for comp_name in Comparison.allowed: - for boolexpr in bool_exprs(solver): - for rhs in [NUM_VAR, BOOL_VAR, 1, BoolVal(True)]: - if comp_name == '>': - # >1 is unsat for boolean expressions, so change it to 0 - if isinstance(rhs, int) and rhs == 1: - rhs = 0 - if isinstance(rhs, BoolVal) and rhs.args[0] == True: - rhs = BoolVal(False) - yield Comparison(comp_name, boolexpr, rhs) - # Generate all possible boolean expressions def bool_exprs(solver): @@ -120,6 +148,7 @@ def bool_exprs(solver): Generate all boolean expressions: - Boolean operators: and([Var]), or([Var]) (CPMpy class 'Operator', is_bool()) - Boolean equality: Var == Var (CPMpy class 'Comparison') + - Global constraints """ names = [(name, arity) for name, (arity, is_bool) in Operator.allowed.items() if is_bool] @@ -140,36 +169,83 @@ def bool_exprs(solver): yield Comparison(eq_name, *BOOL_ARGS[:2]) for cpm_cons in global_constraints(solver): - if cpm_cons.is_bool(): - yield cpm_cons + yield cpm_cons def global_constraints(solver): """ Generate all global constraints - AllDifferent, AllEqual, Circuit, Minimum, Maximum, Element, - Xor, Cumulative + Xor, Cumulative, NValue, Count, Increasing, Decreasing, IncreasingStrict, DecreasingStrict, LexLessEq, LexLess """ - global_cons = [AllDifferent, AllEqual, Minimum, Maximum] - for global_type in global_cons: - cons = global_type(NUM_ARGS) - if solver not in EXCLUDE_GLOBAL or cons.name not in EXCLUDE_GLOBAL[solver]: - yield cons - - # "special" constructors - if solver not in EXCLUDE_GLOBAL or "element" not in EXCLUDE_GLOBAL[solver]: - yield cpm_array(NUM_ARGS)[NUM_VAR] - - if solver not in EXCLUDE_GLOBAL or "xor" not in EXCLUDE_GLOBAL[solver]: - yield Xor(BOOL_ARGS) - - if solver not in EXCLUDE_GLOBAL or "cumulative" not in EXCLUDE_GLOBAL[solver]: - s = intvar(0,10,shape=3,name="start") - e = intvar(0,10,shape=3,name="end") - dur = [1,4,3] - demand = [4,5,7] - cap = 10 - yield Cumulative(s, dur, e, demand, cap) - + classes = inspect.getmembers(cpmpy.expressions.globalconstraints, inspect.isclass) + classes = [(name, cls) for name, cls in classes if issubclass(cls, GlobalConstraint) and name != "GlobalConstraint"] + classes = [(name, cls) for name, cls in classes if name not in EXCLUDE_GLOBAL.get(solver, {})] + + for name, cls in classes: + if solver in EXCLUDE_GLOBAL and name in EXCLUDE_GLOBAL[solver]: + continue + + if name == "Xor": + expr = cls(BOOL_ARGS) + elif name == "Inverse": + expr = cls(NUM_ARGS, [1,0,2]) + elif name == "Table": + expr = cls(NUM_ARGS, [[0,1,2],[1,2,0],[1,0,2]]) + elif name == "NegativeTable": + expr = cls(NUM_ARGS, [[0, 1, 2], [1, 2, 0], [1, 0, 2]]) + elif name == "IfThenElse": + expr = cls(*BOOL_ARGS) + elif name == "InDomain": + expr = cls(NUM_VAR, [0,1,6]) + elif name == "Cumulative": + s = intvar(0, 10, shape=3, name="start") + e = intvar(0, 10, shape=3, name="end") + dur = [1, 4, 3] + demand = [4, 5, 7] + cap = 10 + expr = Cumulative(s, dur, e, demand, cap) + elif name == "GlobalCardinalityCount": + vals = [1, 2, 3] + cnts = intvar(0,10,shape=3) + expr = cls(NUM_ARGS, vals, cnts) + elif name == "AllDifferentExceptN": + expr = cls(NUM_ARGS, NUM_VAR) + elif name == "AllEqualExceptN": + expr = cls(NUM_ARGS, NUM_VAR) + elif name == "Precedence": + x = intvar(0,5, shape=3, name="x") + expr = cls(x, [3,1,0]) + elif name == "NoOverlap": + s = intvar(0, 10, shape=3, name="start") + e = intvar(0, 10, shape=3, name="end") + dur = [1,4,3] + expr = cls(s, dur, e) + elif name == "GlobalCardinalityCount": + vals = [1, 2, 3] + cnts = intvar(0,10,shape=3) + expr = cls(NUM_ARGS, vals, cnts) + elif name == "LexLessEq": + X = intvar(0, 3, shape=3) + Y = intvar(0, 3, shape=3) + expr = LexLessEq(X, Y) + elif name == "LexLess": + X = intvar(0, 3, shape=3) + Y = intvar(0, 3, shape=3) + expr = LexLess(X, Y) + elif name == "LexChainLess": + X = intvar(0, 3, shape=(3,3)) + expr = LexChainLess(X) + elif name == "LexChainLessEq": + X = intvar(0, 3, shape=(3,3)) + expr = LexChainLess(X) + else: # default constructor, list of numvars + expr= cls(NUM_ARGS) + + if solver in EXCLUDE_GLOBAL and name in EXCLUDE_GLOBAL[solver]: + continue + else: + yield expr + def reify_imply_exprs(solver): """ @@ -178,44 +254,59 @@ def reify_imply_exprs(solver): Var -> Boolexpr (CPMpy class 'Operator', is_bool()) """ for bool_expr in bool_exprs(solver): - if solver not in EXCLUDE_IMPL or \ - bool_expr.name not in EXCLUDE_IMPL[solver]: - yield bool_expr.implies(BOOL_VAR) - yield BOOL_VAR.implies(bool_expr) - yield bool_expr == BOOL_VAR + yield bool_expr.implies(BOOL_VAR) + yield BOOL_VAR.implies(bool_expr) + yield bool_expr == BOOL_VAR for comp_expr in comp_constraints(solver): lhs, rhs = comp_expr.args - if solver not in EXCLUDE_IMPL or \ - (not isinstance(lhs, Expression) or lhs.name not in EXCLUDE_IMPL[solver]) and \ - (not isinstance(rhs, Expression) or rhs.name not in EXCLUDE_IMPL[solver]): - yield comp_expr.implies(BOOL_VAR) - yield BOOL_VAR.implies(comp_expr) - yield comp_expr == BOOL_VAR + yield comp_expr.implies(BOOL_VAR) + yield BOOL_VAR.implies(comp_expr) + yield comp_expr == BOOL_VAR + + +def verify(cons): + assert argval(cons) + assert cons.value() -@pytest.mark.parametrize(("solver","constraint"),_generate_inputs(bool_exprs), ids=str) +@pytest.mark.parametrize(("solver","constraint"),list(_generate_inputs(bool_exprs)), ids=str) def test_bool_constaints(solver, constraint): """ Tests boolean constraint by posting it to the solver and checking the value after solve. """ - assert SolverLookup.get(solver, Model(constraint)).solve() - assert constraint.value() + if ALL_SOLS: + n_sols = SolverLookup.get(solver, Model(constraint)).solveAll(display=lambda: verify(constraint)) + assert n_sols >= 1 + else: + assert SolverLookup.get(solver, Model(constraint)).solve() + assert argval(constraint) + assert constraint.value() -@pytest.mark.parametrize(("solver","constraint"), _generate_inputs(comp_constraints), ids=str) +@pytest.mark.parametrize(("solver","constraint"), list(_generate_inputs(comp_constraints)), ids=str) def test_comparison_constraints(solver, constraint): """ Tests comparison constraint by posting it to the solver and checking the value after solve. """ - assert SolverLookup.get(solver,Model(constraint)).solve() - assert constraint.value() + if ALL_SOLS: + n_sols = SolverLookup.get(solver, Model(constraint)).solveAll(display= lambda: verify(constraint)) + assert n_sols >= 1 + else: + assert SolverLookup.get(solver,Model(constraint)).solve() + assert argval(constraint) + assert constraint.value() -@pytest.mark.parametrize(("solver","constraint"), _generate_inputs(reify_imply_exprs), ids=str) +@pytest.mark.parametrize(("solver","constraint"), list(_generate_inputs(reify_imply_exprs)), ids=str) def test_reify_imply_constraints(solver, constraint): """ Tests boolean expression by posting it to solver and checking the value after solve. """ - assert SolverLookup.get(solver, Model(constraint)).solve() - assert constraint.value() + if ALL_SOLS: + n_sols = SolverLookup.get(solver, Model(constraint)).solveAll(display=lambda: verify(constraint)) + assert n_sols >= 1 + else: + assert SolverLookup.get(solver, Model(constraint)).solve() + assert argval(constraint) + assert constraint.value() diff --git a/tests/test_direct.py b/tests/test_direct.py index 63da961ed..71b94d78b 100644 --- a/tests/test_direct.py +++ b/tests/test_direct.py @@ -4,7 +4,7 @@ import pytest from cpmpy import * -from cpmpy.solvers import CPM_gurobi, CPM_pysat, CPM_minizinc, CPM_pysdd, CPM_z3, CPM_exact +from cpmpy.solvers import CPM_gurobi, CPM_pysat, CPM_minizinc, CPM_pysdd, CPM_z3, CPM_exact, CPM_choco class TestDirectORTools(unittest.TestCase): @@ -127,4 +127,17 @@ def test_direct_poly(self): self.assertEqual(y.value(), poly_val) +@pytest.mark.skipif(not CPM_choco.supported(), + reason="pychoco not installed") +class TestDirectChoco(unittest.TestCase): + + def test_direct_global(self): + iv = intvar(1,9, shape=3) + + model = SolverLookup.get("choco") + + model += DirectConstraint("increasing", iv) + model += iv[1] < iv[0] + + self.assertFalse(model.solve()) diff --git a/tests/test_expressions.py b/tests/test_expressions.py index 8921aa5ea..7ff22e88e 100644 --- a/tests/test_expressions.py +++ b/tests/test_expressions.py @@ -6,7 +6,7 @@ from cpmpy.expressions import * from cpmpy.expressions.variables import NDVarArray from cpmpy.expressions.core import Operator, Expression -from cpmpy.expressions.utils import get_bounds +from cpmpy.expressions.utils import get_bounds, argval class TestComparison(unittest.TestCase): def test_comps(self): @@ -437,6 +437,7 @@ def test_incomplete_func(self): if cp.SolverLookup.lookup("z3").supported(): self.assertTrue(m.solve(solver="z3")) # ortools does not support divisor spanning 0 work here self.assertRaises(IncompleteFunctionError, cons.value) + self.assertFalse(argval(cons)) # mayhem cons = (arr[10 // (a - b)] == 1).implies(p) diff --git a/tests/test_globalconstraints.py b/tests/test_globalconstraints.py index 07cd8184e..6dc630125 100644 --- a/tests/test_globalconstraints.py +++ b/tests/test_globalconstraints.py @@ -5,7 +5,7 @@ import cpmpy as cp from cpmpy.expressions.globalfunctions import GlobalFunction -from cpmpy.exceptions import TypeError +from cpmpy.exceptions import TypeError, NotSupportedError from cpmpy.solvers import CPM_minizinc @@ -94,12 +94,57 @@ def test_alldifferent_except0(self): bv = cp.boolvar() self.assertTrue(cp.Model([cp.AllDifferentExcept0(iv[0], bv)]).solve()) + def test_alldifferent_except_n(self): + # test known input/outputs + tuples = [ + ((1, 2, 3), True), + ((1, 2, 1), False), + ((0, 1, 2), True), + ((2, 0, 3), True), + ((2, 0, 2), True), + ((0, 0, 2), False), + ] + iv = cp.intvar(0, 4, shape=3) + c = cp.AllDifferentExceptN(iv, 2) + for (vals, oracle) in tuples: + ret = cp.Model(c, iv == vals).solve() + assert (ret == oracle), f"Mismatch solve for {vals, oracle}" + # don't try this at home, forcibly overwrite variable values (so even when ret=false) + for (var, val) in zip(iv, vals): + var._value = val + assert (c.value() == oracle), f"Wrong value function for {vals, oracle}" + + # and some more + iv = cp.intvar(-8, 8, shape=3) + self.assertTrue(cp.Model([cp.AllDifferentExceptN(iv,2)]).solve()) + self.assertTrue(cp.AllDifferentExceptN(iv,4).value()) + self.assertTrue(cp.Model([cp.AllDifferentExceptN(iv,7), iv == [7, 7, 1]]).solve()) + self.assertTrue(cp.AllDifferentExceptN(iv,7).value()) + + # test with mixed types + bv = cp.boolvar() + self.assertTrue(cp.Model([cp.AllDifferentExceptN([iv[0], bv],4)]).solve()) + + # test with list of n + iv = cp.intvar(0, 4, shape=7) + self.assertFalse(cp.Model([cp.AllDifferentExceptN([iv], [7,8])]).solve()) + self.assertTrue(cp.Model([cp.AllDifferentExceptN([iv], [4, 1])]).solve()) + def test_not_alldifferentexcept0(self): iv = cp.intvar(-8, 8, shape=3) self.assertTrue(cp.Model([~cp.AllDifferentExcept0(iv)]).solve()) self.assertFalse(cp.AllDifferentExcept0(iv).value()) self.assertFalse(cp.Model([~cp.AllDifferentExcept0(iv), iv == [0, 0, 1]]).solve()) + def test_alldifferent_onearg(self): + iv = cp.intvar(0,10) + for s, cls in cp.SolverLookup.base_solvers(): + print(s) + if cls.supported(): + try: + self.assertTrue(cp.Model(cp.AllDifferent([iv])).solve(solver=s)) + except (NotImplementedError, NotSupportedError): + pass def test_circuit(self): """ @@ -126,6 +171,7 @@ def test_circuit(self): self.assertTrue(model.solve()) self.assertTrue(cp.Circuit(x).value()) + def test_not_circuit(self): x = cp.intvar(lb=0, ub=2, shape=3) circuit = cp.Circuit(x) @@ -139,25 +185,25 @@ def test_not_circuit(self): self.assertFalse(cp.Model([circuit, ~circuit]).solve()) - circuit_sols = set() - not_circuit_sols = set() + all_sols = set() + not_all_sols = set() - circuit_models = cp.Model(circuit).solveAll(display=lambda : circuit_sols.add(tuple(x.value()))) - not_circuit_models = cp.Model(~circuit).solveAll(display=lambda : not_circuit_sols.add(tuple(x.value()))) + circuit_models = cp.Model(circuit).solveAll(display=lambda : all_sols.add(tuple(x.value()))) + not_circuit_models = cp.Model(~circuit).solveAll(display=lambda : not_all_sols.add(tuple(x.value()))) total = cp.Model(x == x).solveAll() - for sol in circuit_sols: + for sol in all_sols: for var,val in zip(x, sol): var._value = val self.assertTrue(circuit.value()) - for sol in not_circuit_sols: + for sol in not_all_sols: for var,val in zip(x, sol): var._value = val self.assertFalse(circuit.value()) - self.assertEqual(total, len(circuit_sols) + len(not_circuit_sols)) + self.assertEqual(total, len(all_sols) + len(not_all_sols)) def test_inverse(self): @@ -185,6 +231,17 @@ def test_inverse(self): # constraint can be used as value self.assertTrue(inv.value()) + def test_inverse_onearg(self): + iv = cp.intvar(0,10) + for s, cls in cp.SolverLookup.base_solvers(): + print(s) + if cls.supported(): + try: + self.assertTrue(cp.Model(cp.Inverse([iv], [0])).solve(solver=s)) + except (NotImplementedError, NotSupportedError): + pass + + def test_InDomain(self): iv = cp.intvar(-8, 8) iv_arr = cp.intvar(-8, 8, shape=5) @@ -213,6 +270,107 @@ def test_InDomain(self): model = cp.Model(cons) self.assertTrue(model.solve()) self.assertIn(iv.value(), vals) + vals = [1, 5, 8, -4] + bv = cp.boolvar() + cons = [cp.InDomain(bv, vals)] + model = cp.Model(cons) + self.assertTrue(model.solve()) + self.assertIn(bv.value(), vals) + vals = [iv2, 5, 8, -4] + bv = cp.boolvar() + cons = [cp.InDomain(bv, vals)] + model = cp.Model(cons) + self.assertTrue(model.solve()) + self.assertIn(bv.value(), vals) + vals = [bv & bv, 5, 8, -4] + bv = cp.boolvar() + cons = [cp.InDomain(bv, vals)] + model = cp.Model(cons) + self.assertTrue(model.solve()) + self.assertIn(bv.value(), vals) + + def test_lex_lesseq(self): + from cpmpy import BoolVal + X = cp.intvar(0, 3, shape=10) + c1 = X[:-1] == 1 + Y = [1, 1, 1, 1, 1, 1, 1, 1, 1, 1] + c = cp.LexLessEq(X, Y) + c2 = c != (BoolVal(True)) + m = cp.Model([c1, c2]) + self.assertTrue(m.solve()) + self.assertTrue(c2.value()) + self.assertFalse(c.value()) + + Y = cp.intvar(0, 0, shape=10) + c = cp.LexLessEq(X, Y) + m = cp.Model(c) + self.assertTrue(m.solve("ortools")) + from cpmpy.expressions.utils import argval + self.assertTrue(sum(argval(X)) == 0) + + def test_lex_less(self): + from cpmpy import BoolVal + X = cp.intvar(0, 3, shape=10) + c1 = X[:-1] == 1 + Y = [1, 1, 1, 1, 1, 1, 1, 1, 1, 1] + c = cp.LexLess(X, Y) + c2 = c != (BoolVal(True)) + m = cp.Model([c1, c2]) + self.assertTrue(m.solve()) + self.assertTrue(c2.value()) + self.assertFalse(c.value()) + + Y = cp.intvar(0, 0, shape=10) + c = cp.LexLess(X, Y) + m = cp.Model(c) + self.assertFalse(m.solve("ortools")) + + Z = [0, 0, 0, 0, 0, 0, 0, 0, 0, 1] + c = cp.LexLess(X, Z) + m = cp.Model(c) + self.assertTrue(m.solve("ortools")) + from cpmpy.expressions.utils import argval + self.assertTrue(sum(argval(X)) == 0) + + + def test_lex_chain(self): + from cpmpy import BoolVal + X = cp.intvar(0, 3, shape=10) + c1 = X[:-1] == 1 + Y = [1, 1, 1, 1, 1, 1, 1, 1, 1, 1] + c = cp.LexChainLess([X, Y]) + c2 = c != (BoolVal(True)) + m = cp.Model([c1, c2]) + self.assertTrue(m.solve()) + self.assertTrue(c2.value()) + self.assertFalse(c.value()) + + Y = cp.intvar(0, 0, shape=10) + c = cp.LexChainLessEq([X, Y]) + m = cp.Model(c) + self.assertTrue(m.solve("ortools")) + from cpmpy.expressions.utils import argval + self.assertTrue(sum(argval(X)) == 0) + + Z = cp.intvar(0, 1, shape=(3,2)) + c = cp.LexChainLess(Z) + m = cp.Model(c) + self.assertTrue(m.solve()) + self.assertTrue(sum(argval(Z[0])) == 0) + self.assertTrue(sum(argval(Z[1])) == 1) + self.assertTrue(sum(argval(Z[2])) >= 1) + + + def test_indomain_onearg(self): + + iv = cp.intvar(0, 10) + for s, cls in cp.SolverLookup.base_solvers(): + print(s) + if cls.supported(): + try: + self.assertTrue(cp.Model(cp.InDomain(iv, [2])).solve(solver=s)) + except (NotImplementedError, NotSupportedError): + pass def test_table(self): iv = cp.intvar(-8,8,3) @@ -242,6 +400,59 @@ def test_table(self): model = cp.Model(constraints[0].decompose()) self.assertFalse(model.solve()) + def test_negative_table(self): + iv = cp.intvar(-8,8,3) + + constraints = [cp.NegativeTable([iv[0], iv[1], iv[2]], [ (5, 2, 2)])] + model = cp.Model(constraints) + self.assertTrue(model.solve()) + + model = cp.Model(constraints[0].decompose()) + self.assertTrue(model.solve()) + + constraints = [cp.NegativeTable(iv, [[10, 8, 2], [5, 2, 2]])] + model = cp.Model(constraints) + self.assertTrue(model.solve()) + + model = cp.Model(constraints[0].decompose()) + self.assertTrue(model.solve()) + + self.assertTrue(cp.NegativeTable(iv, [[10, 8, 2], [5, 2, 2]]).value()) + + constraints = [~cp.NegativeTable(iv, [[10, 8, 2], [5, 2, 2]])] + model = cp.Model(constraints) + self.assertTrue(model.solve()) + self.assertFalse(cp.NegativeTable(iv, [[10, 8, 2], [5, 2, 2]]).value()) + self.assertTrue(cp.Table(iv, [[10, 8, 2], [5, 2, 2]]).value()) + + constraints = [cp.NegativeTable(iv, [[10, 8, 2], [5, 9, 2]])] + model = cp.Model(constraints) + self.assertTrue(model.solve()) + + constraints = [cp.NegativeTable(iv, [[10, 8, 2], [5, 9, 2]])] + model = cp.Model(constraints[0].decompose()) + self.assertTrue(model.solve()) + + constraints = [cp.NegativeTable(iv, [[10, 8, 2], [5, 9, 2]]), cp.Table(iv, [[10, 8, 2], [5, 9, 2]])] + model = cp.Model(constraints) + self.assertFalse(model.solve()) + + constraints = [cp.NegativeTable(iv, [[10, 8, 2], [5, 9, 2]]), cp.Table(iv, [[10, 8, 2], [5, 9, 2]])] + model = cp.Model(constraints[0].decompose()) + model += constraints[1].decompose() + self.assertFalse(model.solve()) + + def test_table_onearg(self): + + iv = cp.intvar(0, 10) + for s, cls in cp.SolverLookup.base_solvers(): + print(s) + if cls.supported(): + try: + self.assertTrue(cp.Model(cp.Table([iv], [[0]])).solve(solver=s)) + except (NotImplementedError, NotSupportedError): + pass + def test_minimum(self): iv = cp.intvar(-8, 8, 3) constraints = [cp.Minimum(iv) + 9 == 8] @@ -253,6 +464,17 @@ def test_minimum(self): self.assertTrue(model.solve()) self.assertEqual(str(min(iv.value())), '4') + def test_minimum_onearg(self): + + iv = cp.intvar(0, 10) + for s, cls in cp.SolverLookup.base_solvers(): + print(s) + if cls.supported(): + try: + self.assertTrue(cp.Model(cp.min([iv]) == 0).solve(solver=s)) + except (NotImplementedError, NotSupportedError): + pass + def test_maximum(self): iv = cp.intvar(-8, 8, 3) constraints = [cp.Maximum(iv) + 9 <= 8] @@ -264,6 +486,17 @@ def test_maximum(self): self.assertTrue(model.solve()) self.assertNotEqual(str(max(iv.value())), '4') + def test_maximum_onearg(self): + + iv = cp.intvar(0, 10) + for s, cls in cp.SolverLookup.base_solvers(): + print(s) + if cls.supported(): + try: + self.assertTrue(cp.Model(cp.max([iv]) == 0).solve(solver=s)) + except (NotImplementedError, NotSupportedError): + pass + def test_abs(self): from cpmpy.transformations.decompose_global import decompose_in_tree iv = cp.intvar(-8, 8) @@ -298,7 +531,7 @@ def test_element(self): self.assertEqual(iv.value()[idx.value()], 8) # test 2-D iv = cp.intvar(-8, 8, shape=(3, 3)) - a,b = cp.intvar(0, 3, shape=2) + a,b = cp.intvar(0, 2, shape=2) cons = iv[a,b] == 8 model = cp.Model(cons) self.assertTrue(model.solve()) @@ -311,6 +544,18 @@ def test_element(self): self.assertTrue(cons.value()) self.assertEqual(arr[a.value(), b.value()], 1) + def test_element_onearg(self): + + iv = cp.intvar(0, 10) + idx = cp.intvar(0,0) + for s, cls in cp.SolverLookup.base_solvers(): + print(s) + if cls.supported(): + try: + self.assertTrue(cp.Model(cp.Element([iv],idx) == 0).solve(solver=s)) + except (NotImplementedError, NotSupportedError): + pass + def test_xor(self): bv = cp.boolvar(5) self.assertTrue(cp.Model(cp.Xor(bv)).solve()) @@ -401,6 +646,24 @@ def test_cumulative_single_demand(self): self.assertTrue(m.solve(solver="ortools")) self.assertTrue(m.solve(solver="minizinc")) + @pytest.mark.skipif(not CPM_minizinc.supported(), + reason="Minizinc not installed") + def test_cumulative_nested(self): + start = cp.intvar(0, 10, name="start", shape=3) + dur = [5,5,5] + end = cp.intvar(0, 10, name="end", shape=3) + demand = [5,5,9] + capacity = 10 + bv = cp.boolvar() + + cons = cp.Cumulative([start], [dur], [end], [demand], capacity) + + m = cp.Model(bv.implies(cons), start + dur != end) + + self.assertTrue(m.solve(solver="ortools")) + self.assertTrue(m.solve(solver="minizinc")) + + def test_cumulative_no_np(self): start = cp.intvar(0, 10, 4, "start") @@ -480,6 +743,17 @@ def test_not_global_cardinality_count(self): self.assertFalse(all(cp.Count(iv, val[i]).value() == occ[i] for i in range(len(val)))) self.assertTrue(~cp.GlobalCardinalityCount([iv[0],iv[2],iv[1],iv[4],iv[3]], val, occ).value()) + def test_gcc_onearg(self): + + iv = cp.intvar(0, 10) + for s, cls in cp.SolverLookup.base_solvers(): + print(s) + if cls.supported(): + try: + self.assertTrue(cp.Model(cp.GlobalCardinalityCount([iv], [3],[1])).solve(solver=s)) + except (NotImplementedError, NotSupportedError): + pass + def test_count(self): iv = cp.intvar(-8, 8, shape=3) self.assertTrue(cp.Model([iv[0] == 0, iv[1] != 1, iv[2] != 2, cp.Count(iv, 0) == 3]).solve()) @@ -497,6 +771,118 @@ def test_count(self): self.assertTrue(cp.Model(cp.Count([iv[0],iv[2],iv[1]], x) > y).solve()) + def test_count_onearg(self): + + iv = cp.intvar(0, 10) + for s, cls in cp.SolverLookup.base_solvers(): + print(s) + if cls.supported(): + try: + self.assertTrue(cp.Model(cp.Count([iv], 1) == 0).solve(solver=s)) + except (NotImplementedError, NotSupportedError): + pass + + def test_nvalue(self): + + iv = cp.intvar(-8, 8, shape=3) + cnt = cp.intvar(0,10) + + self.assertFalse(cp.Model(cp.all(iv == 1), cp.NValue(iv) > 1).solve()) + self.assertTrue(cp.Model(cp.all(iv == 1), cp.NValue(iv) > cnt).solve()) + self.assertGreater(len(set(iv.value())), cnt.value()) + + self.assertTrue(cp.Model(cp.NValue(iv) != cnt).solve()) + self.assertTrue(cp.Model(cp.NValue(iv) >= cnt).solve()) + self.assertTrue(cp.Model(cp.NValue(iv) <= cnt).solve()) + self.assertTrue(cp.Model(cp.NValue(iv) < cnt).solve()) + self.assertTrue(cp.Model(cp.NValue(iv) > cnt).solve()) + + # test nested + bv = cp.boolvar() + cons = bv == (cp.NValue(iv) <= 2) + + def check_true(): + self.assertTrue(cons.value()) + cp.Model(cons).solveAll(display=check_true) + + def test_nvalue_except(self): + + iv = cp.intvar(-8, 8, shape=3) + cnt = cp.intvar(0, 10) + + + self.assertFalse(cp.Model(cp.all(iv == 1), cp.NValueExcept(iv, 6) > 1).solve()) + self.assertTrue(cp.Model(cp.NValueExcept(iv, 10) > 1).solve()) + self.assertTrue(cp.Model(cp.all(iv == 1), cp.NValueExcept(iv, 1) == 0).solve()) + self.assertTrue(cp.Model(cp.all(iv == 1), cp.NValueExcept(iv, 6) > cnt).solve()) + self.assertGreater(len(set(iv.value())), cnt.value()) + + val = 6 + self.assertTrue(cp.Model(cp.NValueExcept(iv, val) != cnt).solve()) + self.assertTrue(cp.Model(cp.NValueExcept(iv, val) >= cnt).solve()) + self.assertTrue(cp.Model(cp.NValueExcept(iv, val) <= cnt).solve()) + self.assertTrue(cp.Model(cp.NValueExcept(iv, val) < cnt).solve()) + self.assertTrue(cp.Model(cp.NValueExcept(iv, val) > cnt).solve()) + + # test nested + bv = cp.boolvar() + cons = bv == (cp.NValueExcept(iv, val) <= 2) + + def check_true(): + self.assertTrue(cons.value()) + + cp.Model(cons).solveAll(display=check_true) + + @pytest.mark.skipif(not CPM_minizinc.supported(), + reason="Minizinc not installed") + def test_nvalue_minizinc(self): + iv = cp.intvar(-8, 8, shape=3) + cnt = cp.intvar(0, 10) + + self.assertFalse(cp.Model(cp.all(iv == 1), cp.NValue(iv) > 1).solve('minizinc')) + self.assertTrue(cp.Model(cp.all(iv == 1), cp.NValue(iv) > cnt).solve('minizinc')) + self.assertGreater(len(set(iv.value())), cnt.value()) + + self.assertTrue(cp.Model(cp.NValue(iv) != cnt).solve('minizinc')) + self.assertTrue(cp.Model(cp.NValue(iv) >= cnt).solve('minizinc')) + self.assertTrue(cp.Model(cp.NValue(iv) <= cnt).solve('minizinc')) + self.assertTrue(cp.Model(cp.NValue(iv) < cnt).solve('minizinc')) + self.assertTrue(cp.Model(cp.NValue(iv) > cnt).solve('minizinc')) + + # test nested + bv = cp.boolvar() + cons = bv == (cp.NValue(iv) <= 2) + + def check_true(): + self.assertTrue(cons.value()) + + cp.Model(cons).solveAll(solver='minizinc') + + + def test_precedence(self): + iv = cp.intvar(0,5, shape=6, name="x") + + cons = cp.Precedence(iv, [0,2,1]) + self.assertTrue(cp.Model([cons, iv == [5,0,2,0,0,1]]).solve()) + self.assertTrue(cons.value()) + self.assertTrue(cp.Model([cons, iv == [0,0,0,0,0,0]]).solve()) + self.assertTrue(cons.value()) + self.assertFalse(cp.Model([cons, iv == [0,1,2,0,0,0]]).solve()) + + + def test_no_overlap(self): + start = cp.intvar(0,5, shape=3) + end = cp.intvar(0,5, shape=3) + cons = cp.NoOverlap(start, [2,1,1], end) + self.assertTrue(cp.Model(cons).solve()) + self.assertTrue(cons.value()) + self.assertTrue(cp.Model(cons.decompose()).solve()) + self.assertTrue(cons.value()) + + def check_val(): + assert cons.value() is False + + cp.Model(~cons).solveAll(display=check_val) class TestBounds(unittest.TestCase): def test_bounds_minimum(self): @@ -586,7 +972,103 @@ def test_allEqual(self): a = cp.boolvar() self.assertTrue(cp.Model([cp.AllEqual(x,y,-1)]).solve()) self.assertTrue(cp.Model([cp.AllEqual(a,b,False, a | b)]).solve()) - #self.assertTrue(cp.Model([cp.AllEqual(x,y,b)]).solve()) + self.assertFalse(cp.Model([cp.AllEqual(x,y,b)]).solve()) + + def test_allEqualExceptn(self): + x = cp.intvar(-8, 8) + y = cp.intvar(-7, -1) + b = cp.boolvar() + a = cp.boolvar() + self.assertTrue(cp.Model([cp.AllEqualExceptN([x,y,-1],211)]).solve()) + self.assertTrue(cp.Model([cp.AllEqualExceptN([x,y,-1,4],4)]).solve()) + self.assertTrue(cp.Model([cp.AllEqualExceptN([x,y,-1,4],-1)]).solve()) + self.assertTrue(cp.Model([cp.AllEqualExceptN([a,b,False, a | b], 4)]).solve()) + self.assertTrue(cp.Model([cp.AllEqualExceptN([a,b,False, a | b], 0)]).solve()) + self.assertTrue(cp.Model([cp.AllEqualExceptN([a,b,False, a | b, y], -1)]).solve()) + + # test with list of n + iv = cp.intvar(0, 4, shape=7) + self.assertFalse(cp.Model([cp.AllEqualExceptN([iv], [7,8]), iv[0] != iv[1]]).solve()) + self.assertTrue(cp.Model([cp.AllEqualExceptN([iv], [4, 1]), iv[0] != iv[1]]).solve()) + + def test_not_allEqualExceptn(self): + x = cp.intvar(lb=0, ub=3, shape=3) + n = 2 + constr = cp.AllEqualExceptN(x,n) + + model = cp.Model([~constr, x == [1, 2, 1]]) + self.assertFalse(model.solve()) + + model = cp.Model([~constr]) + self.assertTrue(model.solve()) + self.assertFalse(constr.value()) + + self.assertFalse(cp.Model([constr, ~constr]).solve()) + + all_sols = set() + not_all_sols = set() + + circuit_models = cp.Model(constr).solveAll(display=lambda: all_sols.add(tuple(x.value()))) + not_circuit_models = cp.Model(~constr).solveAll(display=lambda: not_all_sols.add(tuple(x.value()))) + + total = cp.Model(x == x).solveAll() + + for sol in all_sols: + for var, val in zip(x, sol): + var._value = val + self.assertTrue(constr.value()) + + for sol in not_all_sols: + for var, val in zip(x, sol): + var._value = val + self.assertFalse(constr.value()) + + self.assertEqual(total, len(all_sols) + len(not_all_sols)) + + + def test_increasing(self): + x = cp.intvar(-8, 8) + y = cp.intvar(-7, -1) + b = cp.boolvar() + a = cp.boolvar() + self.assertTrue(cp.Model([cp.Increasing(x,y)]).solve()) + self.assertTrue(cp.Model([cp.Increasing(a,b)]).solve()) + self.assertTrue(cp.Model([cp.Increasing(x,y,b)]).solve()) + z = cp.intvar(2,5) + self.assertFalse(cp.Model([cp.Increasing(z,b)]).solve()) + + def test_decreasing(self): + x = cp.intvar(-8, 8) + y = cp.intvar(-7, -1) + b = cp.boolvar() + a = cp.boolvar() + self.assertTrue(cp.Model([cp.Decreasing(x,y)]).solve()) + self.assertTrue(cp.Model([cp.Decreasing(a,b)]).solve()) + self.assertFalse(cp.Model([cp.Decreasing(x,y,b)]).solve()) + z = cp.intvar(2,5) + self.assertTrue(cp.Model([cp.Decreasing(z,b)]).solve()) + + def test_increasing_strict(self): + x = cp.intvar(-8, 8) + y = cp.intvar(-7, -1) + b = cp.boolvar() + a = cp.boolvar() + self.assertTrue(cp.Model([cp.IncreasingStrict(x,y)]).solve()) + self.assertTrue(cp.Model([cp.IncreasingStrict(a,b)]).solve()) + self.assertTrue(cp.Model([cp.IncreasingStrict(x,y,b)]).solve()) + z = cp.intvar(1,5) + self.assertFalse(cp.Model([cp.IncreasingStrict(z,b)]).solve()) + + def test_decreasing_strict(self): + x = cp.intvar(-8, 8) + y = cp.intvar(-7, 0) + b = cp.boolvar() + a = cp.boolvar() + self.assertTrue(cp.Model([cp.DecreasingStrict(x,y)]).solve()) + self.assertTrue(cp.Model([cp.DecreasingStrict(a,b)]).solve()) + self.assertFalse(cp.Model([cp.DecreasingStrict(x,y,b)]).solve()) + z = cp.intvar(1,5) + self.assertTrue(cp.Model([cp.DecreasingStrict(z,b)]).solve()) def test_circuit(self): x = cp.intvar(-8, 8) @@ -681,7 +1163,7 @@ def test_gcc(self): b = cp.boolvar() a = cp.boolvar() - self.assertTrue(cp.Model([cp.GlobalCardinalityCount([x,y], [z,q], [h,v])]).solve()) + # type checks self.assertRaises(TypeError, cp.GlobalCardinalityCount, [x,y], [x,False], [h,v]) self.assertRaises(TypeError, cp.GlobalCardinalityCount, [x,y], [z,b], [h,v]) self.assertRaises(TypeError, cp.GlobalCardinalityCount, [b,a], [a,b], [h,v]) @@ -689,6 +1171,13 @@ def test_gcc(self): self.assertRaises(TypeError, cp.GlobalCardinalityCount, [x, y], [x, h], [True, v]) self.assertRaises(TypeError, cp.GlobalCardinalityCount, [x, y], [x, h], [v, a]) + iv = cp.intvar(0,10, shape=3) + SOLVERNAMES = [name for name, solver in cp.SolverLookup.base_solvers() if solver.supported()] + for name in SOLVERNAMES: + if name in ("pysat", "pysdd"): continue + self.assertTrue(cp.Model([cp.GlobalCardinalityCount(iv, [1,4], [1,1])]).solve(solver=name)) + # test closed version + self.assertFalse(cp.Model(cp.GlobalCardinalityCount(iv, [1,4], [0,0], closed=True)).solve(solver=name)) def test_count(self): x = cp.intvar(0, 1) @@ -699,4 +1188,32 @@ def test_count(self): a = cp.boolvar() self.assertTrue(cp.Model([cp.Count([x,y],z) == 1]).solve()) - self.assertRaises(TypeError, cp.Count, [x,y],[x,False]) \ No newline at end of file + self.assertRaises(TypeError, cp.Count, [x,y],[x,False]) + + def test_among(self): + + iv = cp.intvar(0,10, shape=3, name="x") + + for name, cls in cp.SolverLookup.base_solvers(): + + if cls.supported() is False: + continue + try: + self.assertTrue(cp.Model([cp.Among(iv, [1,2]) == 3]).solve(solver=name)) + self.assertTrue(all(x.value() in [1,2] for x in iv)) + self.assertTrue(cp.Model([cp.Among(iv, [1,100]) > 2]).solve(solver=name)) + self.assertTrue(all(x.value() == 1 for x in iv)) + except NotSupportedError: + continue + + + def test_table(self): + iv = cp.intvar(-8,8,3) + + constraints = [cp.Table([iv[0], [iv[1], iv[2]]], [ (5, 2, 2)])] # not flatlist, should work + model = cp.Model(constraints) + self.assertTrue(model.solve()) + + self.assertRaises(TypeError, cp.Table, [iv[0], iv[1], iv[2], 5], [(5, 2, 2)]) + self.assertRaises(TypeError, cp.Table, [iv[0], iv[1], iv[2], [5]], [(5, 2, 2)]) + self.assertRaises(TypeError, cp.Table, [iv[0], iv[1], iv[2], ['a']], [(5, 2, 2)]) diff --git a/tests/test_model.py b/tests/test_model.py index 9667e8f8e..9b917d8fa 100644 --- a/tests/test_model.py +++ b/tests/test_model.py @@ -108,4 +108,11 @@ def test_deepcopy(self): m2 = copy.deepcopy(m) for cons in flatlist(m2.constraints): - self.assertTrue(cons.value()) \ No newline at end of file + self.assertTrue(cons.value()) + + + def test_unknown_solver(self): + + model = cp.Model(cp.any(cp.boolvar(shape=3))) + + self.assertRaises(ValueError, lambda : model.solve(solver="notasolver")) \ No newline at end of file diff --git a/tests/test_solveAll.py b/tests/test_solveAll.py index 88b31a3e5..a8d0dd82b 100644 --- a/tests/test_solveAll.py +++ b/tests/test_solveAll.py @@ -29,8 +29,7 @@ def test_solveall_no_obj(self): count = solver.solveAll(solution_limit=1000, display=add_sol) self.assertEqual(3, count) self.assertSetEqual(sols, - {"[True, True]","[True, False]","[False, True]"}) - + {"[True, True]", "[True, False]", "[False, True]"}) def test_solveall_with_obj(self): diff --git a/tests/test_solvers.py b/tests/test_solvers.py index e8f06849b..55108c91d 100644 --- a/tests/test_solvers.py +++ b/tests/test_solvers.py @@ -8,6 +8,8 @@ from cpmpy.solvers.minizinc import CPM_minizinc from cpmpy.solvers.gurobi import CPM_gurobi from cpmpy.solvers.exact import CPM_exact +from cpmpy.solvers.choco import CPM_choco + from cpmpy.exceptions import MinizincNameException @@ -526,6 +528,118 @@ def test_false(self): if cls.supported(): self.assertFalse(m.solve(solver=name)) + @pytest.mark.skipif(not CPM_choco.supported(), + reason="pychoco not installed") + def test_choco(self): + bv = cp.boolvar(shape=3) + iv = cp.intvar(0, 9, shape=3) + # circular 'bigger then', UNSAT + m = cp.Model([ + bv[0].implies(iv[0] > iv[1]), + bv[1].implies(iv[1] > iv[2]), + bv[2].implies(iv[2] > iv[0]) + ]) + m += sum(bv) == len(bv) + s = cp.SolverLookup.get("choco", m) + + self.assertFalse(s.solve()) + + m = cp.Model(~(iv[0] != iv[1])) + s = cp.SolverLookup.get("choco", m) + self.assertTrue(s.solve()) + + m = cp.Model((iv[0] == 0) & ((iv[0] != iv[1]) == 0)) + s = cp.SolverLookup.get("choco", m) + self.assertTrue(s.solve()) + + m = cp.Model([~bv, ~((iv[0] + abs(iv[1])) == sum(iv))]) + s = cp.SolverLookup.get("choco", m) + self.assertTrue(s.solve()) + + @pytest.mark.skipif(not CPM_choco.supported(), + reason="pychoco not installed") + def test_choco_element(self): + + # test 1-D + iv = cp.intvar(-8, 8, 3) + idx = cp.intvar(-8, 8) + # test directly the constraint + constraints = [cp.Element(iv, idx) == 8] + model = cp.Model(constraints) + s = cp.SolverLookup.get("choco", model) + self.assertTrue(s.solve()) + self.assertTrue(iv.value()[idx.value()] == 8) + self.assertTrue(cp.Element(iv, idx).value() == 8) + # test through __get_item__ + constraints = [iv[idx] == 8] + model = cp.Model(constraints) + s = cp.SolverLookup.get("choco", model) + self.assertTrue(s.solve()) + self.assertTrue(iv.value()[idx.value()] == 8) + self.assertTrue(cp.Element(iv, idx).value() == 8) + # test 2-D + iv = cp.intvar(-8, 8, shape=(3, 3)) + idx = cp.intvar(0, 3) + idx2 = cp.intvar(0, 3) + constraints = [iv[idx, idx2] == 8] + model = cp.Model(constraints) + s = cp.SolverLookup.get("choco", model) + self.assertTrue(s.solve()) + self.assertTrue(iv.value()[idx.value(), idx2.value()] == 8) + + @pytest.mark.skipif(not CPM_choco.supported(), + reason="pychoco not installed") + def test_choco_gcc_alldiff(self): + + iv = cp.intvar(-8, 8, shape=5) + occ = cp.intvar(0, len(iv), shape=3) + val = [1, 4, 5] + model = cp.Model([cp.GlobalCardinalityCount(iv, val, occ)]) + solver = cp.SolverLookup.get("choco", model) + self.assertTrue(solver.solve()) + self.assertTrue(cp.GlobalCardinalityCount(iv, val, occ).value()) + self.assertTrue(all(cp.Count(iv, val[i]).value() == occ[i].value() for i in range(len(val)))) + occ = [2, 3, 0] + model = cp.Model([cp.GlobalCardinalityCount(iv, val, occ), cp.AllDifferent(val)]) + solver = cp.SolverLookup.get("choco", model) + self.assertTrue(solver.solve()) + self.assertTrue(cp.GlobalCardinalityCount(iv, val, occ).value()) + self.assertTrue(all(cp.Count(iv, val[i]).value() == occ[i] for i in range(len(val)))) + self.assertTrue(cp.GlobalCardinalityCount([iv[0],iv[2],iv[1],iv[4],iv[3]], val, occ).value()) + + @pytest.mark.skipif(not CPM_choco.supported(), + reason="pychoco not installed") + def test_choco_inverse(self): + from cpmpy.solvers.ortools import CPM_ortools + + fwd = cp.intvar(0, 9, shape=10) + rev = cp.intvar(0, 9, shape=10) + + # Fixed value for `fwd` + fixed_fwd = [9, 4, 7, 2, 1, 3, 8, 6, 0, 5] + # Inverse of the above + expected_inverse = [8, 4, 3, 5, 1, 9, 7, 2, 6, 0] + + model = cp.Model(cp.Inverse(fwd, rev), fwd == fixed_fwd) + + solver = cp.SolverLookup.get("choco", model) + self.assertTrue(solver.solve()) + self.assertEqual(list(rev.value()), expected_inverse) + + @pytest.mark.skipif(not CPM_choco.supported(), + reason="pychoco not installed") + def test_choco_objective(self): + iv = cp.intvar(0,10, shape=2) + m = cp.Model(iv >= 1, iv <= 5, maximize=sum(iv)) + s = cp.SolverLookup.get("choco", m) + self.assertTrue( s.solve() ) + self.assertEqual( s.objective_value(), 10) + + m = cp.Model(iv >= 1, iv <= 5, minimize=sum(iv)) + s = cp.SolverLookup.get("choco", m) + self.assertTrue( s.solve() ) + self.assertEqual(s.objective_value(), 2) + @pytest.mark.skipif(not CPM_gurobi.supported(), reason="Gurobi not installed") def test_gurobi_element(self): @@ -555,3 +669,41 @@ def test_gurobi_element(self): s = cp.SolverLookup.get("gurobi", model) self.assertTrue(s.solve()) self.assertTrue(iv.value()[idx.value(), idx2.value()] == 8) + + + def test_vars_not_removed(self): + bvs = cp.boolvar(shape=3) + m = cp.Model([cp.any(bvs) <= 2]) + for name, cls in cp.SolverLookup.base_solvers(): + print(f"Testing with {name}") + if cls.supported(): + # reset value for vars + bvs.clear() + self.assertTrue(m.solve(solver=name)) + for v in bvs: + self.assertIsNotNone(v.value()) + #test solve_all + sols = set() + if name == 'gurobi': + self.assertEqual(m.solveAll(solver=name,solution_limit=20, display=lambda: sols.add(tuple([x.value() for x in bvs]))), 8) #test number of solutions is valid + self.assertEqual(m.solveAll(solver=name,solution_limit=20), 8) #test number of solutions is valid, no display + else: + self.assertEqual(m.solveAll(solver=name, display=lambda: sols.add(tuple([x.value() for x in bvs]))), 8) #test number of solutions is valid + self.assertEqual(m.solveAll(solver=name), 8) #test number of solutions is valid, no display + #test unique sols, should be same number + self.assertEqual(len(sols),8) + + + @pytest.mark.skipif(not CPM_minizinc.supported(), + reason="Minizinc not installed") + def test_count_mzn(self): + # bug #461 + from cpmpy.expressions.core import Operator + + iv = cp.intvar(0,10, shape=3) + x = cp.intvar(0,1) + y = cp.intvar(0,1) + wsum = Operator("wsum", [[1,2,3],[x,y,cp.Count(iv,3)]]) + + m = cp.Model([x + y == 2, wsum == 9]) + self.assertTrue(m.solve(solver="minizinc")) diff --git a/tests/test_tool_dimacs.py b/tests/test_tool_dimacs.py new file mode 100644 index 000000000..2b9b877a4 --- /dev/null +++ b/tests/test_tool_dimacs.py @@ -0,0 +1,88 @@ +import unittest +import tempfile + +import cpmpy as cp +from cpmpy.tools.dimacs import read_dimacs, write_dimacs +from cpmpy.transformations.get_variables import get_variables_model +class CNFTool(unittest.TestCase): + + + def setUp(self) -> None: + self.tmpfile = tempfile.NamedTemporaryFile() + + def test_read_cnf(self): + + """ + a | b | c, + ~b | ~c, + ~a + """ + cnf_txt = "p cnf \n-2 -3 0\n3 2 1 0\n-1 0\n" + with open(self.tmpfile.name, "w") as f: + f.write(cnf_txt) + + model = read_dimacs(self.tmpfile.name) + vars = sorted(get_variables_model(model), key=str) + + sols = set() + addsol = lambda : sols.add(tuple([v.value() for v in vars])) + + self.assertEqual(model.solveAll(display=addsol), 2) + self.assertSetEqual(sols, {(False, False, True), (False, True, False)}) + + + def test_badly_formatted(self): + + cases = [ + "p cnf 2 1\n1 \n2 \n0", "p cnf 2 1\n \n1 2 0" + ] + + for cnf_txt in cases: + with open(self.tmpfile.name, "w") as f: + f.write(cnf_txt) + + m = read_dimacs(self.tmpfile.name) + self.assertEqual(len(m.constraints), 1) + self.assertEqual(m.solveAll(), 3) + + def test_read_bigint(self): + + cnf_txt = "p cnf \n-2 -300 0\n300 2 1 0\n-1 0\n" + with open(self.tmpfile.name, "w") as f: + f.write(cnf_txt) + + model = read_dimacs(self.tmpfile.name) + vars = sorted(get_variables_model(model), key=str) + + self.assertEqual(model.solveAll(), 2) + + def test_with_comments(self): + cnf_txt = "c this file starts with some comments\nc\np cnf \n-2 -3 0\n3 2 1 0\n-1 0\n" + + with open(self.tmpfile.name, "w") as f: + f.write(cnf_txt) + + model = read_dimacs(self.tmpfile.name) + vars = sorted(get_variables_model(model), key=str) + + sols = set() + addsol = lambda : sols.add(tuple([v.value() for v in vars])) + + self.assertEqual(model.solveAll(display=addsol), 2) + self.assertSetEqual(sols, {(False, False, True), (False, True, False)}) + + def test_write_cnf(self): + + a,b,c = [cp.boolvar(name=n) for n in "abc"] + + m = cp.Model() + m += cp.any([a,b,c]) + m += b.implies(~c) + m += a <= 0 + + cnf_txt = write_dimacs(m) + gt_cnf = "p cnf 3 3\n1 2 3 0\n-2 -3 0\n-1 0\n" + + self.assertEqual(cnf_txt, gt_cnf) + + diff --git a/tests/test_trans_simplify.py b/tests/test_trans_simplify.py index b44c980ef..454e0c0f7 100644 --- a/tests/test_trans_simplify.py +++ b/tests/test_trans_simplify.py @@ -100,3 +100,11 @@ def test_with_floats(self): self.assertEqual(str(self.transform(expr)), '[(iv[0]) == (iv[1])]') + def test_nested_boolval(self): + + bv = cp.boolvar(name="bv") + x = cp.intvar(0, 3, name="x") + cons = (x == 2) == (bv == 4) + self.assertEquals(str(self.transform(cons)), "[x != 2]") + self.assertTrue(cp.Model(cons).solve()) + diff --git a/tests/test_transf_decompose.py b/tests/test_transf_decompose.py index e7649b166..afe9c2b83 100644 --- a/tests/test_transf_decompose.py +++ b/tests/test_transf_decompose.py @@ -89,15 +89,15 @@ def test_decompose_nested(self): cons = [min(ivs) == max(ivs)] self.assertEqual(str(decompose_in_tree(cons, supported={"min"})), - "[or([(x) >= (min(x,y,z)), (y) >= (min(x,y,z)), (z) >= (min(x,y,z))]), (x) <= (min(x,y,z)), (y) <= (min(x,y,z)), (z) <= (min(x,y,z))]") + '[(IV0) == (min(x,y,z)), or([(x) >= (IV0), (y) >= (IV0), (z) >= (IV0)]), (x) <= (IV0), (y) <= (IV0), (z) <= (IV0)]') self.assertEqual(str(decompose_in_tree(cons, supported={"max"})), - "[or([(x) <= (max(x,y,z)), (y) <= (max(x,y,z)), (z) <= (max(x,y,z))]), (x) >= (max(x,y,z)), (y) >= (max(x,y,z)), (z) >= (max(x,y,z))]") + "[(IV1) == (max(x,y,z)), or([(x) <= (IV1), (y) <= (IV1), (z) <= (IV1)]), (x) >= (IV1), (y) >= (IV1), (z) >= (IV1)]") # numerical in non-comparison context cons = [AllEqual([min(ivs[:-1]),ivs[-1]])] self.assertEqual(str(decompose_in_tree(cons, supported={"allequal"})), - "[allequal(IV0,z), ((x) <= (IV0)) or ((y) <= (IV0)), (x) >= (IV0), (y) >= (IV0)]") + "[allequal(IV2,z), (IV3) == (IV2), ((x) <= (IV3)) or ((y) <= (IV3)), (x) >= (IV3), (y) >= (IV3)]") self.assertEqual(str(decompose_in_tree(cons, supported={"min"})), "[(min(x,y)) == (z)]")