Repository: CPMpy/cpmpy

Issue:
choco solver crash using gcc constraint
crashed solver: choco
correctly working solvers: ortools, minizinc, z3, gurobi

In the following program:
- the error `'int' object has no attribute 'handle'` occurs in `make_intvar_array` of pychoco, so not sure if this is a cpmpy bug.
- the error disappears when the 3 gets changed into a variable (e.g. x)

```python
import cpmpy as cp

m = cp.Model()

x = cp.intvar(-3, 5, shape=1)
y = cp.intvar(-3, 5, shape=1)
z = cp.intvar(-3, 5, shape=1)

m += (cp.GlobalCardinalityCount([3, y, z],[1, 2, 3],[0,0,3])) >= 1

# following works correctly
# m += x == 3
# m += (cp.GlobalCardinalityCount([x, y, z],[1, 2, 3],[0,0,3])) >= 1

print(m.solve(solver='ortools'))   # returns no error
print(m.solve(solver='minizinc'))  # returns no error
print(m.solve(solver='z3'))        # returns no error
print(m.solve(solver='choco'))     # returns error: 'int' object has no attribute 'handle'
print(m.solve(solver='gurobi'))    # returns no error
```

Discussion:

Solution patch:
diff --git a/.github/workflows/python-test.yml b/.github/workflows/python-test.yml
index 9959d167..8653e773 100644
--- a/.github/workflows/python-test.yml
+++ b/.github/workflows/python-test.yml
@@ -14,6 +14,7 @@ jobs:
           python -m pip install --upgrade pip
           if [ -f requirements.txt ]; then pip install -r requirements.txt; fi
           pip install ".[test, z3, choco, exact, pysat, pysdd, choco, minizinc]"
+          pip install pypblib  # dependency of pysat for PB constraints
           pip install pytest-xdist
           sudo snap install minizinc --classic
       - name: Test with pytest
diff --git a/cpmpy/solvers/choco.py b/cpmpy/solvers/choco.py
index 10cfe9c6..85ccac66 100644
--- a/cpmpy/solvers/choco.py
+++ b/cpmpy/solvers/choco.py
@@ -592,7 +592,7 @@ class CPM_choco(SolverInterface):
                 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)
+                return self.chc_model.global_cardinality(self._to_vars(vars), self.solver_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.")
 
diff --git a/cpmpy/solvers/exact.py b/cpmpy/solvers/exact.py
index 2323dd45..33a5200f 100644
--- a/cpmpy/solvers/exact.py
+++ b/cpmpy/solvers/exact.py
@@ -54,12 +54,12 @@ from pkg_resources import VersionConflict
 
 from .solver_interface import SolverInterface, SolverStatus, ExitStatus
 from ..expressions.core import *
-from ..expressions.variables import _BoolVarImpl, NegBoolView, _IntVarImpl, _NumVarImpl, intvar
+from ..expressions.variables import _BoolVarImpl, NegBoolView, _IntVarImpl, _NumVarImpl
 from ..transformations.comparison import only_numexpr_equality
 from ..transformations.flatten_model import flatten_constraint, flatten_objective
 from ..transformations.get_variables import get_variables
 from ..transformations.decompose_global import decompose_in_tree
-from ..transformations.linearize import linearize_constraint, only_positive_bv
+from ..transformations.linearize import linearize_constraint, only_positive_bv, only_positive_bv_wsum
 from ..transformations.reification import only_implies, reify_rewrite, only_bv_reifies
 from ..transformations.normalize import toplevel_list
 from ..transformations.safening import no_partial_functions
@@ -361,10 +361,11 @@ class CPM_exact(SolverInterface):
         self.objective_ = expr
         self.objective_is_min_ = minimize
 
-        # make objective function non-nested
+        # make objective function non-nested and with positive BoolVars only
         (flat_obj, flat_cons) = flatten_objective(expr)
-        self += flat_cons  # add potentially created constraints
+        flat_obj = only_positive_bv_wsum(flat_obj)  # remove negboolviews
         self.user_vars.update(get_variables(flat_obj))  # add objvars to vars
+        self += flat_cons  # add potentially created constraints
 
         # make objective function or variable and post
         xct_cfvars,xct_rhs = self._make_numexpr(flat_obj,0)
diff --git a/cpmpy/solvers/gurobi.py b/cpmpy/solvers/gurobi.py
index 0b8afdde..3aaad500 100644
--- a/cpmpy/solvers/gurobi.py
+++ b/cpmpy/solvers/gurobi.py
@@ -52,7 +52,7 @@ from ..transformations.comparison import only_numexpr_equality
 from ..transformations.decompose_global import decompose_in_tree
 from ..transformations.flatten_model import flatten_constraint, flatten_objective
 from ..transformations.get_variables import get_variables
-from ..transformations.linearize import linearize_constraint, only_positive_bv
+from ..transformations.linearize import linearize_constraint, only_positive_bv, only_positive_bv_wsum
 from ..transformations.normalize import toplevel_list
 from ..transformations.reification import only_implies, reify_rewrite, only_bv_reifies
 from ..transformations.safening import no_partial_functions
@@ -102,8 +102,7 @@ class CPM_gurobi(SolverInterface):
             global GRB_ENV
             if GRB_ENV is None:
                 # initialise the native gurobi model object
-                GRB_ENV = gp.Env()
-                GRB_ENV.setParam("OutputFlag", 0)
+                GRB_ENV = gp.Env(params={"OutputFlag": 0})
                 GRB_ENV.start()
             return True
         except Exception as e:
@@ -268,9 +267,10 @@ class CPM_gurobi(SolverInterface):
         from gurobipy import GRB
 
         # make objective function non-nested
-        (flat_obj, flat_cons) = (flatten_objective(expr))
+        (flat_obj, flat_cons) = flatten_objective(expr)
+        flat_obj = only_positive_bv_wsum(flat_obj)  # remove negboolviews
+        get_variables(flat_obj, collect=self.user_vars)  # add potentially created variables
         self += flat_cons
-        get_variables(flat_obj, collect=self.user_vars)  # add potentially created constraints
 
         # make objective function or variable and post
         obj = self._make_numexpr(flat_obj)
diff --git a/cpmpy/solvers/pysat.py b/cpmpy/solvers/pysat.py
index 6e8e614e..e7f6557f 100644
--- a/cpmpy/solvers/pysat.py
+++ b/cpmpy/solvers/pysat.py
@@ -21,12 +21,17 @@
     Installation
     ============
 
-    Requires that the 'python-sat' package is installed. If you want to also solve
-    pseudo-Boolean constraints, you should also install its optional dependency 'pblib', as follows:
+    Requires that the 'python-sat' package is installed:
 
     .. code-block:: console
 
-        $ pip install python-sat[pblib]
+        $ pip install pysat
+
+    If you want to also solve pseudo-Boolean constraints, you should also install its optional dependency 'pypblib', as follows:
+
+    .. code-block:: console
+
+        $ pip install pypblib
 
     See detailed installation instructions at:
     https://pysathq.github.io/installation
@@ -82,8 +87,6 @@ class CPM_pysat(SolverInterface):
 
     """
 
-    IMPORT_PYSAT_PBLIB_ERROR = ImportError("PySAT was installed without optional dependency `pblib`, which is required to encode (weighted) sums. To install PySAT with recommended optional dependencies, run `pip install cpmpy[pysat]`.")
-
     @staticmethod
     def supported():
         # try to import the package
@@ -93,6 +96,19 @@ class CPM_pysat(SolverInterface):
             # while we need the 'python-sat' package, some more checks:
             from pysat.formula import IDPool
             from pysat.solvers import Solver
+
+            from pysat import card
+            CPM_pysat._card = card  # native
+
+            # try to import pypblib and avoid ever re-import by setting `_pb`
+            if not hasattr(CPM_pysat, ("_pb")):
+                try:
+                    from pysat import pb  # require pypblib
+                    """The `pysat.pb` module if its dependency `pypblib` installed, `None` if we have not checked it yet, or `False` if we checked and it is *not* installed"""
+                    CPM_pysat._pb = pb
+                except (ModuleNotFoundError, NameError):  # pysat returns the wrong error type (latter i/o former)
+                    CPM_pysat._pb = None  # not installed, avoid reimporting
+
             return True
         except ModuleNotFoundError:
             return False
@@ -142,15 +158,6 @@ class CPM_pysat(SolverInterface):
         elif subsolver.startswith('pysat:'):
             subsolver = subsolver[6:] # strip 'pysat:'
 
-        # try to import the pblib module once
-        try:
-            from pysat import card, pb
-            self._card = card
-            self._pb = pb
-        except NameError:
-            self._card = None
-            self._pb = None
-
         # initialise the native solver object
         self.pysat_vpool = IDPool()
         self.pysat_solver = Solver(use_timer=True, name=subsolver)
@@ -441,8 +448,6 @@ class CPM_pysat(SolverInterface):
 
     def _pysat_cardinality(self, cpm_expr, reified=False):
         """ Convert CPMpy comparison of `sum` (over Boolean variables) into PySAT list of clauses """
-        if self._card is None:
-            raise self.IMPORT_PYSAT_PBLIB_ERROR
 
         # unpack and transform to PySAT argument
         lhs, rhs = cpm_expr.args
@@ -468,9 +473,9 @@ class CPM_pysat(SolverInterface):
             raise ValueError(f"PySAT: Expected Comparison to be either <=, ==, or >=, but was {cpm_expr.name}")
 
     def _pysat_pseudoboolean(self, cpm_expr):
-        """ Convert CPMpy comparison of `wsum` (over Boolean variables) into PySAT list of clauses """
+        """Convert CPMpy comparison of `wsum` (over Boolean variables) into PySAT list of clauses."""
         if self._pb is None:
-            raise self.IMPORT_PYSAT_PBLIB_ERROR
+            raise ImportError("The model contains a PB constraint, for which PySAT needs an additional dependency (PBLib). To install it, run `pip install pypblib`.")
 
         if cpm_expr.args[0].name != "wsum":
             raise NotSupportedError(
@@ -482,6 +487,7 @@ class CPM_pysat(SolverInterface):
         lits = self.solver_vars(lhs.args[1])
         pysat_args = {"weights": lhs.args[0], "lits": lits, "bound": rhs, "vpool":self.pysat_vpool }
 
+
         if cpm_expr.name == "<=":
             return self._pb.PBEnc.atmost(**pysat_args).clauses
         elif cpm_expr.name == ">=":
diff --git a/cpmpy/solvers/utils.py b/cpmpy/solvers/utils.py
index 45c35c9f..e4ea011f 100644
--- a/cpmpy/solvers/utils.py
+++ b/cpmpy/solvers/utils.py
@@ -82,7 +82,27 @@ class SolverLookup():
                ]
 
     @classmethod
-    def solvernames(cls):
+    def print_status(cls):
+        """
+            Print all CPMpy solvers and their installation status on this system.
+        """
+        for (basename, CPM_slv) in cls.base_solvers():
+            if CPM_slv.supported():
+                print(f"{basename}: Supported, ready to use.")
+            else:
+                print(f"{basename}: Not supported (missing Python package, binary or license).")
+
+    @classmethod
+    def supported(cls):
+        """
+            Return the list of names of all solvers (and subsolvers) supported on this system.
+
+            If a solver name is returned, it means that the solver's `.supported()` function returns True
+            and it is hence ready for immediate use
+            (e.g. any separate binaries are also installed if necessary, and licenses are active if needed).
+
+            Typical use case is to use these names in `SolverLookup.get(name)`.
+        """
         names = []
         for (basename, CPM_slv) in cls.base_solvers():
             if CPM_slv.supported():
@@ -93,6 +113,12 @@ class SolverLookup():
                         names.append(basename+":"+subn)
         return names
 
+    @classmethod
+    def solvernames(cls):
+        # The older (more indirectly named) way to get the list of names of *supported* solvers.
+        # Will be deprecated at some point.
+        return cls.supported()
+
     @classmethod
     def get(cls, name=None, model=None, **init_kwargs):
         """
diff --git a/cpmpy/transformations/linearize.py b/cpmpy/transformations/linearize.py
index 13ab44fc..be560aa6 100644
--- a/cpmpy/transformations/linearize.py
+++ b/cpmpy/transformations/linearize.py
@@ -389,43 +389,43 @@ def linearize_constraint(lst_of_expr, supported={"sum","wsum"}, reified=False):
 
 def only_positive_bv(lst_of_expr):
     """
-        Replaces constraints containing :class:`~cpmpy.expressions.variables.NegBoolView` with equivalent expression using only :class:`~cpmpy.expressions.variables.BoolVar`.
-        `lst_of_expr` is expected to be linearized. Only apply after applying :func:`linearize_constraint(cpm_expr) <linearize_constraint>`.
+        Replaces :class:`~cpmpy.expressions.comparison.Comparison` containing :class:`~cpmpy.expressions.variables.NegBoolView` with equivalent expression using only :class:`~cpmpy.expressions.variables.BoolVar`.
+        Comparisons are expected to be linearized. Only apply after applying :func:`linearize_constraint(cpm_expr) <linearize_constraint>`.
 
-        Resulting expression is linear.
+        Resulting expression is linear if the original expression was linear.
     """
     newlist = []
     for cpm_expr in lst_of_expr:
 
         if isinstance(cpm_expr, Comparison):
             lhs, rhs = cpm_expr.args
+            new_lhs = lhs
             new_cons = []
 
-            if isinstance(lhs, _NumVarImpl):
-                if isinstance(lhs,NegBoolView):
-                    lhs, rhs = Operator("wsum",[[-1], [lhs._bv]]), 1 - rhs
-
-            if lhs.name == "sum" and any(isinstance(a, NegBoolView) for a in lhs.args):
-                lhs = Operator("wsum",[[1]*len(lhs.args), lhs.args])
+            if isinstance(lhs, _NumVarImpl) or lhs.name in {"sum","wsum"}:
+                new_lhs, const = only_positive_bv_wsum_const(lhs)
+                rhs -= const
+            else:
+                # other operators in comparison such as "min", "max"
+                nbv_sel = [isinstance(a, NegBoolView) for a in lhs.args]
+                if any(nbv_sel):
+                    new_args = []
+                    for i, nbv in enumerate(nbv_sel):
+                        if nbv:
+                            aux = cp.boolvar()
+                            new_args.append(aux)
+                            new_cons += [aux + lhs.args[i]._bv == 1]  # aux == 1 - arg._bv
+                        else:
+                            new_args.append(lhs.args[i])
 
-            if lhs.name == "wsum":
-                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", [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"}:
-            # other operators in comparison such as "min", "max"
-                lhs = copy.copy(lhs)
-                for i,arg in enumerate(list(lhs.args)):
-                    if isinstance(arg, NegBoolView):
-                        new_arg, cons = get_or_make_var(1 - arg)
-                        lhs.args[i] = new_arg
-                        new_cons += cons
+                    new_lhs = copy.copy(lhs)
+                    new_lhs.update_args(new_args)
 
-            newlist.append(eval_comparison(cpm_expr.name, lhs, rhs))
-            newlist += linearize_constraint(new_cons)
+            if new_lhs != lhs:
+                newlist.append(eval_comparison(cpm_expr.name, new_lhs, rhs))
+                newlist += new_cons  # already linear
+            else:
+                newlist.append(cpm_expr)
 
         # reification
         elif isinstance(cpm_expr, Operator) and cpm_expr.name == "->":
@@ -445,6 +445,91 @@ def only_positive_bv(lst_of_expr):
 
     return newlist
 
+def only_positive_bv_wsum(expr):
+    """
+        Replaces a var/sum/wsum expression containing :class:`~cpmpy.expressions.variables.NegBoolView` with an equivalent expression 
+        using only :class:`~cpmpy.expressions.variables.BoolVar`. 
+
+        It might add a constant term to the expression, if you want the constant separately, use :func:`only_positive_bv_wsum_const`.
+        
+        Arguments:
+        - `cpm_expr`: linear expression (sum, wsum, var)
+        
+        Returns tuple of:
+        - `pos_expr`: linear expression (sum, wsum, var) without NegBoolView
+    """
+    if isinstance(expr, _NumVarImpl) or expr.name in {"sum","wsum"}:
+        pos_expr, const = only_positive_bv_wsum_const(expr)
+        if const == 0:
+            return pos_expr
+        else:
+            assert isinstance(pos_expr, Operator) and pos_expr.name == "wsum", f"unexpected expression, should be wsum but got {pos_expr}"
+            # should we check if it already has a constant term?
+            return Operator("wsum", [pos_expr.args[0]+[1], pos_expr.args[1]+[const]])
+    else:
+        return expr
+
+def only_positive_bv_wsum_const(cpm_expr):
+    """
+        Replaces a var/sum/wsum expression containing :class:`~cpmpy.expressions.variables.NegBoolView` with an equivalent expression 
+        using only :class:`~cpmpy.expressions.variables.BoolVar` as well as a constant term that must be added to the new expression to be equivalent.
+
+        If you want the expression where the constant term is part of the wsum returned, use :func:`only_positive_bv_wsum`.
+        
+        Arguments:
+        - `cpm_expr`: linear expression (sum, wsum, var)
+        
+        Returns tuple of:
+        - `pos_expr`: linear expression (sum, wsum, var) without NegBoolView
+        - `const`: The difference between the original expression and the new expression, 
+                   i.e. a constant term that must be added to pos_expr to be an equivalent linear expression.
+    """
+    if isinstance(cpm_expr, _NumVarImpl):
+        if isinstance(cpm_expr,NegBoolView):
+            return Operator("wsum",[[-1], [cpm_expr._bv]]), 1
+        else:
+            return cpm_expr, 0
+
+    elif cpm_expr.name == "sum":
+        # indicator on arguments being negboolviews
+        nbv_sel = [isinstance(a, NegBoolView) for a in cpm_expr.args]
+        if any(nbv_sel):
+            const = 0
+            weights = []
+            variables = []
+            for i, nbv in enumerate(nbv_sel):
+                if nbv:
+                    const += 1
+                    weights.append(-1)
+                    variables.append(cpm_expr.args[i]._bv)
+                else:
+                    weights.append(1)
+                    variables.append(cpm_expr.args[i])
+            return Operator("wsum", [weights, variables]), const
+        else:
+            return cpm_expr, 0
+
+    elif cpm_expr.name == "wsum":
+        # indicator on arguments of the wsum variable being negboolviews
+        nbv_sel = [isinstance(a, NegBoolView) for a in cpm_expr.args[1]]
+        if any(nbv_sel):
+            # copy weights and variables
+            weights = [w for w in cpm_expr.args[0]]
+            variables = [v for v in cpm_expr.args[1]]
+            const = 0
+            for i, nbv in enumerate(nbv_sel):
+                if nbv:
+                    const += weights[i]
+                    weights[i] = -weights[i]
+                    variables[i] = variables[i]._bv
+            return Operator("wsum", [weights, variables]), const
+        else:
+            return cpm_expr, 0
+
+    else:
+        raise ValueError(f"unexpected expression, should be sum, wsum or var but got {cpm_expr}")
+
+
 def canonical_comparison(lst_of_expr):
 
     lst_of_expr = toplevel_list(lst_of_expr)               # ensure it is a list
@@ -473,6 +558,8 @@ def canonical_comparison(lst_of_expr):
                 lhs2 = []
                 if isinstance(rhs, _NumVarImpl):
                     lhs2, rhs = [-1 * rhs], 0
+                elif isinstance(rhs, BoolVal):
+                    lhs2, rhs = [-1] if rhs.value() else [], 0 
                 elif isinstance(rhs, Operator) and rhs.name == "-":
                     lhs2, rhs = [rhs.args[0]], 0
                 elif isinstance(rhs, Operator) and rhs.name == "sum":
@@ -484,6 +571,7 @@ def canonical_comparison(lst_of_expr):
                                     if isinstance(b, _NumVarImpl)], \
                                     sum(-a * b for a, b in zip(rhs.args[0], rhs.args[1])
                                     if not isinstance(b, _NumVarImpl))
+                
                 # 2) add collected variables to lhs
                 if isinstance(lhs, Operator) and lhs.name == "sum":
                     lhs, rhs = sum([1 * a for a in lhs.args] + lhs2), rhs
@@ -571,4 +659,3 @@ def only_positive_coefficients(lst_of_expr):
             newlist.append(cpm_expr)
 
     return newlist
-
diff --git a/setup.py b/setup.py
index 26f90def..98b9f074 100644
--- a/setup.py
+++ b/setup.py
@@ -44,7 +44,7 @@ setup(
         "choco": ["pychoco>=0.2.1"],
         "exact": ["exact>=2.1.0"],
         "minizinc": ["minizinc"],
-        "pysat": ["python-sat[pblib]"],
+        "pysat": ["python-sat"],
         "gurobi": ["gurobipy"],
         "pysdd": ["pysdd"],
         "gcs": ["gcspy"],


Test patch:
diff --git a/tests/test_pysat_wsum.py b/tests/test_pysat_wsum.py
index 3af3191d..bbd8451f 100644
--- a/tests/test_pysat_wsum.py
+++ b/tests/test_pysat_wsum.py
@@ -1,9 +1,24 @@
 import unittest
+import pytest
 import cpmpy as cp 
 from cpmpy import *
 from cpmpy.solvers.pysat import CPM_pysat
-from cpmpy.transformations.to_cnf import to_cnf
 
+import importlib # can check for modules *without* importing them
+pblib_available = importlib.util.find_spec("pypblib") is not None
+
+@pytest.mark.skipif(pblib_available, reason="`pypblib` is installed")
+def test_pypblib_error():
+    # NOTE if you want to run this but pypblib is already installed, run `pip uninstall pypblib && pip install -e .[pysat]`
+    unittest.TestCase().assertRaises(
+            ImportError, # just solve a pb constraint with pypblib not installed
+            lambda : CPM_pysat(cp.Model(2*cp.boolvar() + 3 * cp.boolvar() + 5 * cp.boolvar() <= 6)).solve()
+        )
+
+    # this one should still work without `pypblib`
+    assert CPM_pysat(cp.Model(1*cp.boolvar() + 1 * cp.boolvar() + 1 * cp.boolvar() <= 2)).solve()
+
+@pytest.mark.skipif(not pblib_available, reason="pypblib not installed")
 class TestEncodePseudoBooleanConstraint(unittest.TestCase):
     def setUp(self):
         self.bv = boolvar(shape=3)
diff --git a/tests/test_trans_linearize.py b/tests/test_trans_linearize.py
index 4a41d77f..5e8a45dd 100644
--- a/tests/test_trans_linearize.py
+++ b/tests/test_trans_linearize.py
@@ -4,7 +4,9 @@ import cpmpy as cp
 from cpmpy.expressions import boolvar, intvar
 from cpmpy.expressions.core import Operator
 from cpmpy.expressions.utils import argvals
-from cpmpy.transformations.linearize import linearize_constraint, canonical_comparison, only_positive_bv, only_positive_coefficients
+from cpmpy.transformations.decompose_global import decompose_global
+from cpmpy.transformations.flatten_model import flatten_objective
+from cpmpy.transformations.linearize import linearize_constraint, canonical_comparison, only_positive_bv, only_positive_coefficients, only_positive_bv_wsum_const, only_positive_bv_wsum
 from cpmpy.expressions.variables import _IntVarImpl, _BoolVarImpl
 
 
@@ -40,9 +42,20 @@ class TestTransLinearize(unittest.TestCase):
             iv = intvar(1, 9)
             e1 = (bv[0] * bv[1] == iv)
             s1 = cp.Model(e1).solve("gurobi")
-            s1 = cp.Model(e1).solve("ortools")
             self.assertTrue(s1)
             self.assertEqual([bv[0].value(), bv[1].value(), iv.value()],[True, True, 1])
+            
+    def test_bug_468(self):
+        from cpmpy.solvers import CPM_exact, CPM_gurobi
+        a, b, c = boolvar(shape=3)
+        m = cp.Model(cp.any([a, b, c]))
+        m.minimize(3*a + 4*~b + 3*~c)
+        if CPM_gurobi.supported():
+            m.solve("gurobi")
+            self.assertEqual([a.value(), b.value(), c.value()], [False, True, True])
+        if CPM_exact.supported():
+            m.solve("exact")
+            self.assertEqual([a.value(), b.value(), c.value()], [False, True, True])
 
     def test_constraint(self):
         x,y,z = [cp.intvar(0,5, name=n) for n in "xyz"]
@@ -385,6 +398,12 @@ class TestVarsLhs(unittest.TestCase):
         self.assertEqual("alldifferent(a,b,c)", str(cons))
 
 class testCanonical_comparison(unittest.TestCase):
+    def setUp(self):
+        _IntVarImpl.counter = 0
+        _BoolVarImpl.counter = 0
+        self.ivars = cp.intvar(1, 10, shape=(5,))
+        self.bvars = cp.boolvar((3,))
+    
     def test_sum(self):
         a,b,c = [cp.intvar(0,10,name=n) for n in "abc"]
         rhs = 5
@@ -465,3 +484,97 @@ class testCanonical_comparison(unittest.TestCase):
     def test_only_positive_bv_implied_by_negated_literal(self):
         p = cp.boolvar(name="p")
         self.assertEqual(str([p <= 0]), str(only_positive_bv(linearize_constraint([~p]))))
+        
+        
+class testOnlyPositiveBv(unittest.TestCase):
+    def setUp(self):
+        _IntVarImpl.counter = 0
+        _BoolVarImpl.counter = 0
+        self.ivars = cp.intvar(1, 10, shape=(5,))
+        self.bvars = cp.boolvar((3,))
+
+    def test_only_positive_bv_wsum_const_positive_literal(self):
+        p = cp.boolvar(name="p")
+        self.assertEqual(str((p, 0)), str(only_positive_bv_wsum_const(p)))
+        
+    def test_only_positive_bv_wsum_const_negated_litral(self):
+        p = cp.boolvar(name="p")
+        self.assertEqual(str((Operator("wsum",[[-1],[p]]), 1)), str(only_positive_bv_wsum_const(~p)))
+        
+    def test_only_positive_bv_wsum_const_sum_positive_input(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        self.assertEqual(str((Operator("sum",[a,b,c]), 0)), str(only_positive_bv_wsum_const(a+b+c)))
+        
+    def test_only_positive_bv_wsum_const_sum(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        self.assertEqual(str((Operator("wsum",[[-1, -1, 1],[a,b,c]]), 2)), str(only_positive_bv_wsum_const(~a+~b+c)))
+        
+    def test_only_positive_bv_wsum_const_wsum_positive_input(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        self.assertEqual(str((Operator("wsum",[[4, 5, 1],[a,b,c]]), 0)), str(only_positive_bv_wsum_const(4*a+5*b+c)))
+        
+    def test_only_positive_bv_wsum_const_wsum(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        self.assertEqual(str((Operator("wsum",[[-4, 5, 1],[a,b,c]]), 4)), str(only_positive_bv_wsum_const(4*~a+5*b+c)))
+        
+    def test_only_positive_bv_wsum_const_non_linear(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        with self.assertRaises(ValueError) as cm:
+            only_positive_bv_wsum_const(~a * b * c)
+            self.assertEqual(str(cm.exception), "unexpected expression, should be sum, wsum or var but got ((~a) * (b)) * (c)")
+            
+    def test_only_positive_bv_wsum_positive_literal(self):
+        p = cp.boolvar(name="p")
+        obj = only_positive_bv_wsum(p)
+        self.assertEqual(str(obj), str(p))
+            
+    def test_only_positive_bv_wsum_negated_literal(self):
+        p = cp.boolvar(name="p")
+        obj = only_positive_bv_wsum(~p)
+        self.assertEqual(str(obj), str(Operator("wsum",[[-1, 1],[p, 1]])))
+        
+    def test_only_positive_bv_wsum_sum_positive_input(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        obj = only_positive_bv_wsum(a + b + c)
+        self.assertEqual(str(obj), str(Operator("sum",[a,b,c])))
+        
+    def test_only_positive_bv_wsum_sum(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        obj = only_positive_bv_wsum(~a + ~b + c)
+        self.assertEqual(str(obj), str(Operator("wsum",[[-1, -1, 1, 1],[a,b,c,2]])))
+        
+    def test_only_positive_bv_wsum_wsum_positive_input(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        obj = only_positive_bv_wsum(4*a + 5*b + c)
+        self.assertEqual(str(obj), str(Operator("wsum",[[4, 5, 1],[a,b,c]])))
+        
+    def test_only_positive_bv_wsum_wsum(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        obj = only_positive_bv_wsum(4*~a + 5*~b + c)
+        self.assertEqual(str(obj), str(Operator("wsum",[[-4, -5, 1, 1],[a,b,c,9]])))
+        
+    def test_only_positive_bv_wsum_non_linear_positive_input(self): # TODO: make boolvars
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        flat_obj, flat_cons = flatten_objective(a * b * c)
+        flat_obj = only_positive_bv_wsum(flat_obj)
+        self.assertEqual(str((flat_obj, flat_cons)), "(IV6, [((IV5) * (c)) == (IV6), ((a) * (b)) == (IV5)])")
+        
+    def test_only_positive_bv_wsum_non_linear(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        flat_obj, flat_cons = flatten_objective(~a * b * c)
+        flat_obj = only_positive_bv_wsum(flat_obj)
+        self.assertEqual(str((flat_obj, flat_cons)), "(IV6, [((IV5) * (c)) == (IV6), ((~a) * (b)) == (IV5)])")
+        
+    def test_only_flat_positive_bv_wsum_max(self):
+        a, b = [cp.boolvar(name=n) for n in "ab"]
+        expr = cp.max((~a * b), (a * ~b ))
+        flat_obj, flat_cons = flatten_objective(expr)
+        # self.assertEqual(str(expr.args), "(max((~a) * (b), (a) * (~b)))")
+        obj = only_positive_bv_wsum(flat_obj)
+        self.assertEqual(str((obj, flat_cons)), "(IV7, [(max(IV5,IV6)) == (IV7), ((~a) * (b)) == (IV5), ((a) * (~b)) == (IV6)])")
+        
+    def test_only_flat_positive_bv_max(self):
+        a, b, c = [cp.boolvar(name=n) for n in "abc"]
+        obj = only_positive_bv(linearize_constraint([cp.max(~a,b) >= c], supported={"sum", "wsum", "max"}))
+        self.assertEqual(str(obj), "[(max(BV3,b)) >= (c), (BV3) + (a) == 1]")
+