diff --git a/doc/OnlineDocs/modeling_extensions/gdp/concepts.rst b/doc/OnlineDocs/modeling_extensions/gdp/concepts.rst index 3f571f988d2..95629bc48fd 100644 --- a/doc/OnlineDocs/modeling_extensions/gdp/concepts.rst +++ b/doc/OnlineDocs/modeling_extensions/gdp/concepts.rst @@ -115,7 +115,7 @@ These logical propositions can include: .. |equiv| replace:: :math:`Y_1 \Leftrightarrow Y_2` .. |land| replace:: :math:`Y_1 \land Y_2` .. |lor| replace:: :math:`Y_1 \lor Y_2` -.. |xor| replace:: :math:`Y_1 \underline{\lor} Y_2` +.. |xor| replace:: :math:`Y_1 \veebar Y_2` .. |impl| replace:: :math:`Y_1 \Rightarrow Y_2` +-----------------+---------+-------------+-------------+-------------+ diff --git a/doc/OnlineDocs/modeling_extensions/gdp/index.rst b/doc/OnlineDocs/modeling_extensions/gdp/index.rst index a5e47333377..0c8529c60cb 100644 --- a/doc/OnlineDocs/modeling_extensions/gdp/index.rst +++ b/doc/OnlineDocs/modeling_extensions/gdp/index.rst @@ -54,7 +54,7 @@ These can be expressed as a disjunction as follows: \text{constraints} \\ \text{for }\textit{else} \end{gathered}\right] \\ - Y_1 \underline{\vee} Y_2 + Y_1 \veebar Y_2 \end{gather*} Here, if the Boolean :math:`Y_1` is ``True``, then the constraints in the first disjunct are enforced; otherwise, the constraints in the second disjunct are enforced. diff --git a/doc/OnlineDocs/modeling_extensions/gdp/modeling.rst b/doc/OnlineDocs/modeling_extensions/gdp/modeling.rst index 9627698634e..3a666e63233 100644 --- a/doc/OnlineDocs/modeling_extensions/gdp/modeling.rst +++ b/doc/OnlineDocs/modeling_extensions/gdp/modeling.rst @@ -86,7 +86,7 @@ When the ``Disjunction`` object constructor is passed a list of lists, the outer By default, Pyomo.GDP ``Disjunction`` objects enforce an implicit "exactly one" relationship among the selection of the disjuncts (generalization of exclusive-OR). That is, exactly one of the ``Disjunct`` indicator variables should take a ``True`` value. - This can be seen as an implicit logical proposition, in our example, :math:`Y_1 \underline{\lor} Y_2`. + This can be seen as an implicit logical proposition, in our example, :math:`Y_1 \veebar Y_2`. Logical Propositions ==================== @@ -139,6 +139,11 @@ Pyomo.GDP logical expression system supported operators and their usage are list | Equivalence | | :code:`Y[1].equivalent_to(Y[2])` | :code:`equivalent(Y[1], Y[2])` | +--------------+------------------------+-----------------------------------+--------------------------------+ +.. note:: + + We omit support for most infix operators, e.g. :code:`Y[1] >> Y[2]`, due to concerns about non-intuitive Python operator precedence. + That is :code:`Y[1] | Y[2] >> Y[3]` would translate to :math:`Y_1 \lor (Y_2 \Rightarrow Y_3)` rather than :math:`(Y_1 \lor Y_2) \Rightarrow Y_3` + In addition, the following constraint-programming-inspired operators are provided: ``exactly``, ``atmost``, and ``atleast``. These predicates enforce, respectively, that exactly, at most, or at least N of their ``BooleanVar`` arguments are ``True``. @@ -148,26 +153,23 @@ Usage: - :code:`atmost(3, Y)` - :code:`exactly(3, Y)` -.. note:: - - We omit support for most infix operators, e.g. :code:`Y[1] >> Y[2]`, due to concerns about non-intuitive Python operator precedence. - That is :code:`Y[1] | Y[2] >> Y[3]` would translate to :math:`Y_1 \lor (Y_2 \Rightarrow Y_3)` rather than :math:`(Y_1 \lor Y_2) \Rightarrow Y_3` - .. doctest:: >>> m = ConcreteModel() >>> m.my_set = RangeSet(4) >>> m.Y = BooleanVar(m.my_set) >>> m.p = LogicalConstraint(expr=atleast(3, m.Y)) - >>> TransformationFactory('core.logical_to_linear').apply_to(m) - >>> m.logic_to_linear.transformed_constraints.pprint() # constraint auto-generated by transformation - transformed_constraints : Size=1, Index=logic_to_linear.transformed_constraints_index, Active=True - Key : Lower : Body : Upper : Active - 1 : 3.0 : Y_asbinary[1] + Y_asbinary[2] + Y_asbinary[3] + Y_asbinary[4] : +Inf : True >>> m.p.pprint() p : Size=1, Index=None, Active=False Key : Body : Active None : atleast(3: [Y[1], Y[2], Y[3], Y[4]]) : False + >>> TransformationFactory('core.logical_to_linear').apply_to(m) + ... + >>> # constraint auto-generated by transformation + >>> m.logic_to_linear.transformed_constraints.pprint() + transformed_constraints : Size=1, Index=logic_to_linear.transformed_constraints_index, Active=True + Key : Lower : Body : Upper : Active + 1 : 3.0 : Y_asbinary[1] + Y_asbinary[2] + Y_asbinary[3] + Y_asbinary[4] : +Inf : True We elaborate on the ``logical_to_linear`` transformation :ref:`on the next page `. @@ -223,8 +225,8 @@ Here, we demonstrate this capability with a toy example: \min~&x\\ \text{s.t.}~&\left[\begin{gathered}Y_1\\x \geq 2\end{gathered}\right] \vee \left[\begin{gathered}Y_2\\x \geq 3\end{gathered}\right]\\ &\left[\begin{gathered}Y_3\\x \leq 8\end{gathered}\right] \vee \left[\begin{gathered}Y_4\\x = 2.5\end{gathered}\right] \\ - &Y_1 \underline{\vee} Y_2\\ - &Y_3 \underline{\vee} Y_4\\ + &Y_1 \veebar Y_2\\ + &Y_3 \veebar Y_4\\ &Y_1 \Rightarrow Y_4 .. doctest:: @@ -359,7 +361,7 @@ In the ``logical_to_linear`` transformation, we automatically convert these spec Additional Examples =================== -The following models all work and are equivalent for :math:`\left[x = 0\right] \underline{\lor} \left[y = 0\right]`: +The following models all work and are equivalent for :math:`\left[x = 0\right] \veebar \left[y = 0\right]`: .. doctest:: diff --git a/doc/OnlineDocs/modeling_extensions/gdp/solving.rst b/doc/OnlineDocs/modeling_extensions/gdp/solving.rst index 131b4790cb9..2f3076862e6 100644 --- a/doc/OnlineDocs/modeling_extensions/gdp/solving.rst +++ b/doc/OnlineDocs/modeling_extensions/gdp/solving.rst @@ -29,15 +29,25 @@ Logical constraints .. note:: - Historically it was required to convert logical propositions to - algebraic form prior to use of the MI(N)LP reformulations and the - GDPopt solver. However, this is mathematically incorrect since these - reformulations convert logical formulations to algebraic formulations. - It is therefore recommended to use both the MI(N)LP reformulations - and GDPopt directly to transform or solve GDPs that include logical - propositions. + Historically users needed to explicitly convert logical propositions + to algebraic form prior to invoking the GDP MI(N)LP reformulations + or the GDPopt solver. However, this is mathematically incorrect + since the GDP MI(N)LP reformulations themselves convert logical + formulations to algebraic formulations. The current recommended + practice is to pass the entire (mixed logical / algebraic) model to + the MI(N)LP reformulations or GDPopt directly. + +There are several approaches to convert logical constraints into +algebraic form. + +Conjunctive Normal Form +^^^^^^^^^^^^^^^^^^^^^^^ -The following transforms logical propositions on the model to algebraic form: +The first transformation (`core.logical_to_linear`) leverages the +`sympy` package to generate the conjunctive normal form of the logical +constraints and then adds the equivalent as a list algebraic +constraints. The following transforms logical propositions on the model +to algebraic form: .. code:: @@ -61,6 +71,29 @@ Following solution of the GDP model, values of the Boolean variables may be upda .. autofunction:: pyomo.core.plugins.transform.logical_to_linear.update_boolean_vars_from_binary +Factorable Programming +^^^^^^^^^^^^^^^^^^^^^^ + +The second transformation (`contrib.logical_to_disjunctive`) leverages +ideas from factorable programming to first generate an equivalent set of +"factored" logical constraints form by traversing each logical +proposition and replacing each logical operator with an additional +Boolean variable and then adding the "simple" logical constraint that +equates the new Boolean variable with the single logical operator. + +The resulting "simple" logical constraints are converted to either MIP +or GDP form: if the constraint contains only Boolean variables, then +then MIP representation is emitted. Logical constraints with mixed +integer-Boolean arguments (e.g., `atmost`, `atleast`, `exactly`, etc.) +are converted to a disjunctive representation. + +As this transformation both avoids the conversion into `sympy` and only +requires a single traversal of each logical constraint, +`contrib.logical_to_disjunctive` is significantly faster than +`core.logical_to_linear` at the cost of a larger model. In practice, +the cost of the larger model is negated by the effectiveness of the MIP +presolve in most solvers. + Reformulation to MI(N)LP ------------------------ diff --git a/pyomo/common/errors.py b/pyomo/common/errors.py index cc5a0231876..730bcc0b115 100644 --- a/pyomo/common/errors.py +++ b/pyomo/common/errors.py @@ -70,3 +70,23 @@ class TempfileContextError(PyomoException, IndexError): """ pass + + +class MouseTrap(PyomoException, NotImplementedError): + """ + Exception class used to throw errors for not-implemented functionality + that might be rational to support (i.e., we already gave you a cookie) + but risks taking Pyomo's flexibility a step beyond what is sane, + or solvable, or communicable to a solver, etc. (i.e., Really? Now you + want a glass of milk too?) + """ + def __init__(self, val): + self.parameter = val + + def __str__(self): + return ("Sorry, mouse, no cookies here!\n\t%s\n" + "\tThis is functionality we think may be rational to " + "support, but is not yet implemented since it might go " + "beyond what can practically be solved. However, please " + "feed the mice: pull requests are always welcome!" + % (repr(self.parameter), )) diff --git a/pyomo/contrib/cp/__init__.py b/pyomo/contrib/cp/__init__.py index bea5d7b24f4..1cc3a24f0e5 100644 --- a/pyomo/contrib/cp/__init__.py +++ b/pyomo/contrib/cp/__init__.py @@ -5,3 +5,5 @@ DocplexWriter, CPOptimizerSolver) from pyomo.contrib.cp.scheduling_expr.step_function_expressions import ( AlwaysIn, Step, Pulse) +# register logical_to_disjunctive transformation +import pyomo.contrib.cp.transform.logical_to_disjunctive_program diff --git a/pyomo/contrib/cp/plugins.py b/pyomo/contrib/cp/plugins.py new file mode 100644 index 00000000000..cc5aa60e89e --- /dev/null +++ b/pyomo/contrib/cp/plugins.py @@ -0,0 +1,15 @@ +# ___________________________________________________________________________ +# +# Pyomo: Python Optimization Modeling Objects +# Copyright (c) 2008-2022 +# National Technology and Engineering Solutions of Sandia, LLC +# Under the terms of Contract DE-NA0003525 with National Technology and +# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain +# rights in this software. +# This software is distributed under the 3-clause BSD License. +# ___________________________________________________________________________ + +def load(): + from . import interval_var + from .repn import docplex_writer + from .transform import logical_to_disjunctive_program diff --git a/pyomo/contrib/cp/repn/docplex_writer.py b/pyomo/contrib/cp/repn/docplex_writer.py index 4d7539810b1..c47d0915add 100644 --- a/pyomo/contrib/cp/repn/docplex_writer.py +++ b/pyomo/contrib/cp/repn/docplex_writer.py @@ -10,8 +10,6 @@ # ___________________________________________________________________________ from pyomo.common.dependencies import attempt_import -cp, docplex_available = attempt_import('docplex.cp.model') -cp_solver, docplex_available = attempt_import('docplex.cp.solver') import itertools import logging @@ -68,6 +66,33 @@ from pyomo.network import Port ### +def _finalize_docplex(module, available): + if not available: + return + _deferred_element_getattr_dispatcher['start_time'] = module.start_of + _deferred_element_getattr_dispatcher['end_time'] = module.end_of + _deferred_element_getattr_dispatcher['length'] = module.length_of + _deferred_element_getattr_dispatcher['is_present'] = module.presence_of + + # Scheduling dispatchers + _before_dispatchers[_START_TIME, _START_TIME] = module.start_before_start + _before_dispatchers[_START_TIME, _END_TIME] = module.start_before_end + _before_dispatchers[_END_TIME, _START_TIME] = module.end_before_start + _before_dispatchers[_END_TIME, _END_TIME] = module.end_before_end + + _at_dispatchers[_START_TIME, _START_TIME] = module.start_at_start + _at_dispatchers[_START_TIME, _END_TIME] = module.start_at_end + _at_dispatchers[_END_TIME, _START_TIME] = module.end_at_start + _at_dispatchers[_END_TIME, _END_TIME] = module.end_at_end + + _time_point_dispatchers[_START_TIME] = module.start_of + _time_point_dispatchers[_END_TIME] = module.end_of + + +cp, docplex_available = attempt_import( + 'docplex.cp.model', callback=_finalize_docplex) +cp_solver, docplex_available = attempt_import('docplex.cp.solver') + logger = logging.getLogger('pyomo.contrib.cp') # These are things that don't need special handling: @@ -217,13 +242,8 @@ def _handle_getitem(visitor, node, *data): 'xor': _XOR, 'equivalent_to': _EQUIVALENT_TO, } -if docplex_available: - _deferred_element_getattr_dispatcher = { - 'start_time': cp.start_of, - 'end_time': cp.end_of, - 'length': cp.length_of, - 'is_present': cp.presence_of, - } +# This will get populated when cp is finally imported +_deferred_element_getattr_dispatcher = {} def _handle_getattr(visitor, node, obj, attr): # We either end up here because we do not yet know the list of variables to @@ -676,25 +696,13 @@ def _handle_call(visitor, node, *args): # Scheduling ## -if docplex_available: - _before_dispatchers = { - (_START_TIME, _START_TIME): cp.start_before_start, - (_START_TIME, _END_TIME): cp.start_before_end, - (_END_TIME, _START_TIME): cp.end_before_start, - (_END_TIME, _END_TIME): cp.end_before_end, - } - _at_dispatchers = { - (_START_TIME, _START_TIME): cp.start_at_start, - (_START_TIME, _END_TIME): cp.start_at_end, - (_END_TIME, _START_TIME): cp.end_at_start, - (_END_TIME, _END_TIME): cp.end_at_end - } - _time_point_dispatchers = { - _START_TIME: cp.start_of, - _END_TIME: cp.end_of, - _GENERAL: lambda x: x, - _ELEMENT_CONSTRAINT: lambda x: x, - } +# This will get populated when cp is finally imported +_before_dispatchers = {} +_at_dispatchers = {} +_time_point_dispatchers = { + _GENERAL: lambda x: x, + _ELEMENT_CONSTRAINT: lambda x: x, +} _non_precedence_types = {_GENERAL, _ELEMENT_CONSTRAINT} diff --git a/pyomo/contrib/cp/tests/test_logical_to_disjunctive.py b/pyomo/contrib/cp/tests/test_logical_to_disjunctive.py new file mode 100755 index 00000000000..fe9ca9ec082 --- /dev/null +++ b/pyomo/contrib/cp/tests/test_logical_to_disjunctive.py @@ -0,0 +1,558 @@ +# ___________________________________________________________________________ +# +# Pyomo: Python Optimization Modeling Objects +# Copyright (c) 2008-2022 +# National Technology and Engineering Solutions of Sandia, LLC +# Under the terms of Contract DE-NA0003525 with National Technology and +# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain +# rights in this software. +# This software is distributed under the 3-clause BSD License. +# ___________________________________________________________________________ + +from pyomo.common.errors import MouseTrap +import pyomo.common.unittest as unittest +from pyomo.contrib.cp.transform.logical_to_disjunctive_program import ( + LogicalToDisjunctive) +from pyomo.contrib.cp.transform.logical_to_disjunctive_walker import ( + LogicalToDisjunctiveVisitor +) +from pyomo.core.expr.compare import assertExpressionsEqual +from pyomo.environ import ( + atmost, atleast, exactly, Block, BooleanVar, Binary, ConcreteModel, + Expression, Integers, land, lnot, lor, LogicalConstraint, Param, value, + Var, TransformationFactory) + +class TestLogicalToDisjunctiveVisitor(unittest.TestCase): + def make_model(self): + m = ConcreteModel() + + m.a = BooleanVar() + m.b = BooleanVar() + m.c = BooleanVar() + + return m + + def test_logical_or(self): + m = self.make_model() + e = lor(m.a, m.b, m.c) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + + visitor.walk_expression(e) + + self.assertIs(m.a.get_associated_binary(), m.z[1]) + self.assertIs(m.b.get_associated_binary(), m.z[2]) + self.assertIs(m.c.get_associated_binary(), m.z[3]) + + self.assertEqual(len(m.cons), 5) + self.assertEqual(len(m.z), 4) + # !z4 v a v b v c + assertExpressionsEqual( + self, m.cons[1].expr, 1 - m.z[4] + m.z[1] + m.z[2] + m.z[3] >= 1) + # z4 v !a + assertExpressionsEqual( + self, m.cons[2].expr, m.z[4] + (1 - m.z[1]) >= 1) + # z4 v !b + assertExpressionsEqual( + self, m.cons[3].expr, m.z[4] + (1 - m.z[2]) >= 1) + # z4 v !c + assertExpressionsEqual( + self, m.cons[4].expr, m.z[4] + (1 - m.z[3]) >= 1) + + # z4 is constrained to be 'True' + assertExpressionsEqual( + self, m.cons[5].expr, m.z[4] >= 1) + + def test_logical_and(self): + m = self.make_model() + e = land(m.a, m.b, m.c) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + + visitor.walk_expression(e) + + self.assertIs(m.a.get_associated_binary(), m.z[1]) + self.assertIs(m.b.get_associated_binary(), m.z[2]) + self.assertIs(m.c.get_associated_binary(), m.z[3]) + + self.assertEqual(len(m.cons), 4) + self.assertEqual(len(m.z), 4) + assertExpressionsEqual( + self, m.cons[1].expr, m.z[4] <= m.z[1]) + assertExpressionsEqual( + self, m.cons[2].expr, m.z[4] <= m.z[2]) + assertExpressionsEqual( + self, m.cons[3].expr, m.z[4] <= m.z[3]) + assertExpressionsEqual( + self, m.cons[4].expr, m.z[4] >= 1) + + def test_logical_not(self): + m = self.make_model() + e = lnot(m.a) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + + visitor.walk_expression(e) + self.assertIs(m.a.get_associated_binary(), m.z[1]) + self.assertEqual(len(m.cons), 2) + self.assertEqual(len(m.z), 2) + assertExpressionsEqual( + self, m.cons[1].expr, m.z[2] == 1 - m.z[1]) + assertExpressionsEqual( + self, m.cons[2].expr, m.z[2] >= 1) + + def test_implication(self): + m = self.make_model() + e = m.a.implies(m.b.land(m.c)) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + + visitor.walk_expression(e) + self.assertIs(m.a.get_associated_binary(), m.z[1]) + self.assertIs(m.b.get_associated_binary(), m.z[2]) + self.assertIs(m.c.get_associated_binary(), m.z[3]) + + self.assertEqual(len(m.cons), 6) + # z4 = b ^ c + assertExpressionsEqual(self, m.cons[1].expr, m.z[4] <= m.z[2]) + assertExpressionsEqual(self, m.cons[2].expr, m.z[4] <= m.z[3]) + # z5 = a -> z4 + # which means z5 = !a v z4 + assertExpressionsEqual(self, m.cons[3].expr, + (1 - m.z[5]) + (1 - m.z[1]) + m.z[4] >= 1) + # z5 >= 1 - z1 + assertExpressionsEqual(self, m.cons[4].expr, + m.z[5] + (1 - (1 - m.z[1])) >= 1) + # z5 >= z4 + assertExpressionsEqual(self, m.cons[5].expr, m.z[5] + (1 - m.z[4]) >= 1) + + # z5 is constrained to be 'True' + assertExpressionsEqual(self, m.cons[6].expr, m.z[5] >= 1) + + def test_equivalence(self): + m = self.make_model() + e = m.a.equivalent_to(m.c) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + + visitor.walk_expression(e) + self.assertIs(m.a.get_associated_binary(), m.z[1]) + self.assertIs(m.c.get_associated_binary(), m.z[2]) + self.assertEqual(len(m.z), 5) + self.assertEqual(len(m.cons), 9) + + assertExpressionsEqual( + self, m.cons[1].expr, (1 - m.z[3]) + (1 - m.z[1]) + m.z[2] >= 1) + assertExpressionsEqual( + self, m.cons[2].expr, 1 - (1 - m.z[1]) + m.z[3] >= 1) + assertExpressionsEqual( + self, m.cons[3].expr, m.z[3] + (1 - m.z[2]) >= 1) + + assertExpressionsEqual( + self, m.cons[4].expr, (1 - m.z[4]) + (1 - m.z[2]) + m.z[1] >= 1) + assertExpressionsEqual( + self, m.cons[5].expr, m.z[4] + (1 - m.z[1]) >= 1) + assertExpressionsEqual( + self, m.cons[6].expr, 1 - (1 - m.z[2]) + m.z[4] >= 1) + + assertExpressionsEqual( + self, m.cons[7].expr, m.z[5] <= m.z[3]) + assertExpressionsEqual( + self, m.cons[8].expr, m.z[5] <= m.z[4]) + + assertExpressionsEqual(self, m.cons[9].expr, m.z[5] >= 1) + + def test_xor(self): + m = self.make_model() + e = m.a.xor(m.b) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + m.disjuncts = visitor.disjuncts + m.disjunctions = visitor.disjunctions + + visitor.walk_expression(e) + self.assertIs(m.a.get_associated_binary(), m.z[1]) + self.assertIs(m.b.get_associated_binary(), m.z[2]) + + self.assertEqual(len(m.z), 2) + self.assertEqual(len(m.cons), 1) + self.assertEqual(len(m.disjuncts), 2) + self.assertEqual(len(m.disjunctions), 1) + + assertExpressionsEqual( + self, m.disjuncts[0].constraint.expr, m.z[1] + m.z[2] == 1) + assertExpressionsEqual( + self, + m.disjuncts[1].disjunction.disjuncts[0].constraint[1].expr, + m.z[1] + m.z[2] <= 0) + assertExpressionsEqual( + self, + m.disjuncts[1].disjunction.disjuncts[1].constraint[1].expr, + m.z[1] + m.z[2] >= 2) + + assertExpressionsEqual(self, m.cons[1].expr, + m.disjuncts[0].binary_indicator_var >= 1) + + def test_at_most(self): + m = self.make_model() + e = atmost(2, m.a, (m.a.land(m.b)), m.c) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + m.disjuncts = visitor.disjuncts + m.disjunctions = visitor.disjunctions + + visitor.walk_expression(e) + self.assertIs(m.a.get_associated_binary(), m.z[1]) + a = m.z[1] + self.assertIs(m.b.get_associated_binary(), m.z[2]) + b = m.z[2] + self.assertIs(m.c.get_associated_binary(), m.z[4]) + c = m.z[4] + + self.assertEqual(len(m.z), 4) + self.assertEqual(len(m.cons), 3) + self.assertEqual(len(m.disjuncts), 2) + self.assertEqual(len(m.disjunctions), 1) + + # z3 = a ^ b + assertExpressionsEqual( + self, m.cons[1].expr, m.z[3] <= a) + assertExpressionsEqual( + self, m.cons[2].expr, m.z[3] <= b) + + # atmost in disjunctive form + assertExpressionsEqual( + self, m.disjuncts[0].constraint.expr, m.z[1] + m.z[3] + m.z[4] <= 2) + assertExpressionsEqual( + self, + m.disjuncts[1].constraint.expr, + m.z[1] + m.z[3] + m.z[4] >= 3) + assertExpressionsEqual(self, m.cons[3].expr, + m.disjuncts[0].binary_indicator_var >= 1) + + def test_at_least(self): + m = self.make_model() + e = atleast(2, m.a, m.b, m.c) + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + m.disjuncts = visitor.disjuncts + m.disjunctions = visitor.disjunctions + + visitor.walk_expression(e) + self.assertIs(m.a.get_associated_binary(), m.z[1]) + a = m.z[1] + self.assertIs(m.b.get_associated_binary(), m.z[2]) + b = m.z[2] + self.assertIs(m.c.get_associated_binary(), m.z[3]) + c = m.z[3] + + self.assertEqual(len(m.z), 3) + self.assertEqual(len(m.cons), 1) + + # atleast in disjunctive form + assertExpressionsEqual( + self, m.disjuncts[0].constraint.expr, m.z[1] + m.z[2] + m.z[3] >= 2) + assertExpressionsEqual( + self, + m.disjuncts[1].constraint.expr, + m.z[1] + m.z[2] + m.z[3] <= 1) + + assertExpressionsEqual(self, m.cons[1].expr, + m.disjuncts[0].binary_indicator_var >= 1) + + def test_no_need_to_walk(self): + m = self.make_model() + e = m.a + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + + visitor.walk_expression(e) + self.assertEqual(len(m.z), 1) + self.assertIs(m.a.get_associated_binary(), m.z[1]) + self.assertEqual(len(m.cons), 1) + assertExpressionsEqual(self, m.cons[1].expr, m.z[1] >= 1) + + def test_binary_already_associated(self): + m = self.make_model() + m.mine = Var(domain=Binary) + m.a.associate_binary_var(m.mine) + + e = m.a.land(m.b) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + + visitor.walk_expression(e) + + self.assertEqual(len(m.z), 2) + self.assertIs(m.b.get_associated_binary(), m.z[1]) + self.assertEqual(len(m.cons), 3) + assertExpressionsEqual(self, m.cons[1].expr, m.z[2] <= m.mine) + assertExpressionsEqual(self, m.cons[2].expr, m.z[2] <= m.z[1]) + assertExpressionsEqual(self, m.cons[3].expr, m.z[2] >= 1) + + # [ESJ 11/22]: We'll probably eventually support all of these examples, but + # for now test that we handle them gracefully: + def test_integer_var_in_at_least(self): + m = self.make_model() + m.x = Var(bounds=(0, 10), domain=Integers) + e = atleast(m.x, m.a, m.b, m.c) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + + with self.assertRaisesRegex( + MouseTrap, + r"The first argument 'x' to " + r"'atleast\(x: \[a, b, c\]\)' is potentially variable. " + r"This may be a mathematically coherent expression; However " + r"it is not yet supported to convert it to a disjunctive " + r"program."): + visitor.walk_expression(e) + + def test_numeric_expression_in_at_most(self): + m = self.make_model() + m.x = Var([1, 2], bounds=(0, 10), domain=Integers) + m.y = Var(domain=Integers) + m.e = Expression(expr=m.x[1]*m.x[2]) + e = atmost(m.e + m.y, m.a, m.b, m.c) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + + with self.assertRaisesRegex( + MouseTrap, + r"The first argument '\(x\[1\]\*x\[2\]\) \+ y' to " + r"'atmost\(\(x\[1\]\*x\[2\]\) \+ y: \[a, b, c\]\)' is " + r"potentially variable. " + r"This may be a mathematically coherent expression; However " + r"it is not yet supported to convert it to a disjunctive " + r"program"): + visitor.walk_expression(e) + + def test_named_expression_in_at_most(self): + m = self.make_model() + m.x = Var([1, 2], bounds=(0, 10), domain=Integers) + m.y = Var(domain=Integers) + m.e = Expression(expr=m.x[1]*m.x[2]) + e = atmost(m.e, m.a, m.b, m.c) + + visitor = LogicalToDisjunctiveVisitor() + m.cons = visitor.constraints + m.z = visitor.z_vars + + with self.assertRaisesRegex( + MouseTrap, + r"The first argument 'e' to " + r"'atmost\(\(x\[1\]\*x\[2\]\): \[a, b, c\]\)' is " + r"potentially variable. " + r"This may be a mathematically coherent expression; However " + r"it is not yet supported to convert it to a disjunctive " + r"program"): + visitor.walk_expression(e) + + def test_relational_expr_as_boolean_atom(self): + m = self.make_model() + m.x = Var() + e = m.a.land(m.x >= 3) + visitor = LogicalToDisjunctiveVisitor() + + with self.assertRaisesRegex( + MouseTrap, + "The RelationalExpression '3 <= x' was used as a Boolean " + "term in a logical proposition. This is not yet supported " + "when transforming to disjunctive form."): + visitor.walk_expression(e) + +class TestLogicalToDisjunctiveTransformation(unittest.TestCase): + def make_model(self): + m = ConcreteModel() + + m.a = BooleanVar() + m.b = BooleanVar([1, 2]) + m.p = Param(initialize=1) + m.p2 = Param([1, 2], mutable=True) + m.p2[1] = 1 + m.p2[2] = 2 + + m.block = Block() + m.block.c1 = LogicalConstraint(expr=m.a.land(m.b[1])) + m.block.c2 = LogicalConstraint( + expr=exactly(m.p2[2], m.a, m.b[1], m.b[2].lor(m.b[1]))) + + m.c1 = LogicalConstraint( + expr=atmost(m.p + m.p2[1], m.a, m.b[1], m.b[2])) + + return m + + def check_and_constraints(self, a, b1, z, transBlock): + assertExpressionsEqual( + self, + transBlock.transformed_constraints[1].expr, + z <= a) + assertExpressionsEqual( + self, + transBlock.transformed_constraints[2].expr, + z <= b1) + assertExpressionsEqual( + self, transBlock.transformed_constraints[3].expr, z >= 1) + + def check_block_c1_transformed(self, m, transBlock): + self.assertFalse(m.block.c1.active) + self.assertIs(m.a.get_associated_binary(), transBlock.auxiliary_vars[1]) + self.assertIs(m.b[1].get_associated_binary(), + transBlock.auxiliary_vars[2]) + self.check_and_constraints(transBlock.auxiliary_vars[1], + transBlock.auxiliary_vars[2], + transBlock.auxiliary_vars[3], transBlock) + + def check_block_exactly(self, a, b1, b2, z4, transBlock): + m = transBlock.model() + + # z[4] = b[2] v b[1] + assertExpressionsEqual( + self, + transBlock.transformed_constraints[4].expr, + (1 - z4) + b2 + b1 >= 1) + assertExpressionsEqual( + self, + transBlock.transformed_constraints[5].expr, + z4 + (1 - b2) >= 1) + assertExpressionsEqual( + self, + transBlock.transformed_constraints[6].expr, + z4 + (1 - b1) >= 1) + + # exactly in disjunctive form + assertExpressionsEqual( + self, + transBlock.auxiliary_disjuncts[0].constraint.expr, + a + b1 + z4 == m.p2[2]) + assertExpressionsEqual( + self, + transBlock.auxiliary_disjuncts[1].disjunction.disjuncts[0].\ + constraint[1].expr, + a + b1 + z4 <= m.p2[2] - 1) + assertExpressionsEqual( + self, + transBlock.auxiliary_disjuncts[1].disjunction.disjuncts[1].\ + constraint[1].expr, + a + b1 + z4 >= m.p2[2] + 1) + + assertExpressionsEqual( + self, + transBlock.transformed_constraints[7].expr, + transBlock.auxiliary_disjuncts[0].binary_indicator_var >= 1) + + def check_block_transformed(self, m): + self.assertFalse(m.block.c2.active) + transBlock = m.block._logical_to_disjunctive + self.assertEqual(len(transBlock.auxiliary_vars), 5) + self.assertEqual(len(transBlock.transformed_constraints), 7) + self.assertEqual(len(transBlock.auxiliary_disjuncts), 2) + self.assertEqual(len(transBlock.auxiliary_disjunctions), 1) + + self.check_block_c1_transformed(m, transBlock) + + self.assertIs(m.b[2].get_associated_binary(), + transBlock.auxiliary_vars[4]) + + z4 = transBlock.auxiliary_vars[5] + a = transBlock.auxiliary_vars[1] + b1 = transBlock.auxiliary_vars[2] + b2 = transBlock.auxiliary_vars[4] + self.check_block_exactly(a, b1, b2, z4, transBlock) + + def test_constraint_target(self): + m = self.make_model() + TransformationFactory('contrib.logical_to_disjunctive').apply_to( + m, + targets=[m.block.c1]) + + transBlock = m.block._logical_to_disjunctive + self.assertEqual(len(transBlock.auxiliary_vars), 3) + self.assertEqual(len(transBlock.transformed_constraints), 3) + self.assertEqual(len(transBlock.auxiliary_disjuncts), 0) + self.assertEqual(len(transBlock.auxiliary_disjunctions), 0) + self.check_block_c1_transformed(m, transBlock) + self.assertTrue(m.block.c2.active) + self.assertTrue(m.c1.active) + + def test_block_target(self): + m = self.make_model() + TransformationFactory('contrib.logical_to_disjunctive').apply_to( + m, + targets=[m.block]) + + self.check_block_transformed(m) + self.assertTrue(m.c1.active) + + def test_transform_block(self): + m = self.make_model() + TransformationFactory('contrib.logical_to_disjunctive').apply_to( + m.block) + + self.check_block_transformed(m) + self.assertTrue(m.c1.active) + + def test_transform_model(self): + m = self.make_model() + TransformationFactory('contrib.logical_to_disjunctive').apply_to(m) + + # c1 got transformed first + self.assertFalse(m.c1.active) + transBlock = m._logical_to_disjunctive + self.assertEqual(len(transBlock.auxiliary_vars), 3) + self.assertEqual(len(transBlock.transformed_constraints), 1) + self.assertEqual(len(transBlock.auxiliary_disjuncts), 2) + self.assertEqual(len(transBlock.auxiliary_disjunctions), 1) + + a = m._logical_to_disjunctive.auxiliary_vars[1] + b1 = m._logical_to_disjunctive.auxiliary_vars[2] + b2 = m._logical_to_disjunctive.auxiliary_vars[3] + + # atmost in disjunctive form + assertExpressionsEqual( + self, + transBlock.auxiliary_disjuncts[0].constraint.expr, + a + b1 + b2 <= 1 + m.p2[1]) + assertExpressionsEqual( + self, + transBlock.auxiliary_disjuncts[1].constraint.expr, + a + b1 + b2 >= 1 + m.p2[1] + 1) + + assertExpressionsEqual( + self, + transBlock.transformed_constraints[1].expr, + transBlock.auxiliary_disjuncts[0].binary_indicator_var >= 1) + + # and everything on the block is transformed too + transBlock = m.block._logical_to_disjunctive + self.assertEqual(len(transBlock.auxiliary_vars), 2) + self.assertEqual(len(transBlock.transformed_constraints), 7) + self.assertEqual(len(transBlock.auxiliary_disjuncts), 2) + self.assertEqual(len(transBlock.auxiliary_disjunctions), 1) + self.check_and_constraints(a, b1, transBlock.auxiliary_vars[1], + transBlock) + self.check_block_exactly(a, b1, b2, transBlock.auxiliary_vars[2], + transBlock) diff --git a/pyomo/contrib/cp/transform/__init__.py b/pyomo/contrib/cp/transform/__init__.py new file mode 100644 index 00000000000..e69de29bb2d diff --git a/pyomo/contrib/cp/transform/logical_to_disjunctive_program.py b/pyomo/contrib/cp/transform/logical_to_disjunctive_program.py new file mode 100644 index 00000000000..20e0455cfc0 --- /dev/null +++ b/pyomo/contrib/cp/transform/logical_to_disjunctive_program.py @@ -0,0 +1,128 @@ +# ___________________________________________________________________________ +# +# Pyomo: Python Optimization Modeling Objects +# Copyright (c) 2008-2022 +# National Technology and Engineering Solutions of Sandia, LLC +# Under the terms of Contract DE-NA0003525 with National Technology and +# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain +# rights in this software. +# This software is distributed under the 3-clause BSD License. +# ___________________________________________________________________________ + +from pyomo.contrib.cp.transform.logical_to_disjunctive_walker import ( + LogicalToDisjunctiveVisitor +) +from pyomo.common.collections import ComponentMap +from pyomo.common.modeling import unique_component_name +from pyomo.common.config import ConfigDict, ConfigValue + +from pyomo.core import ( + TransformationFactory, VarList, Binary, LogicalConstraint, Block, + ConstraintList, Transformation, NonNegativeIntegers) +from pyomo.core.base.block import _BlockData +from pyomo.core.util import target_list +from pyomo.gdp import Disjunct, Disjunction + +@TransformationFactory.register( + "contrib.logical_to_disjunctive", + doc="Convert logical propositions with only Boolean arguments to MIP " + "representation and convert logical expressions with mixed " + "integer-Boolean arguments (such as atleast, atmost, and exactly) to " + "disjunctive representation") +class LogicalToDisjunctive(Transformation): + """ + Re-encode logical constraints as linear constraints, + converting Boolean variables to binary. + """ + CONFIG = ConfigDict('core.logical_to_disjunctive') + CONFIG.declare('targets', ConfigValue( + default=None, + domain=target_list, + description="target or list of targets that will be relaxed", + doc=""" + This specifies the list of LogicalConstraints to transform, or the + list of Blocks or Disjuncts on which to transform all of the + LogicalConstraints. Note that if the transformation is done out + of place, the list of targets should be attached to the model before it + is cloned, and the list will specify the targets on the cloned + instance. + """ + )) + + def _apply_to(self, model, **kwds): + config = self.CONFIG(kwds.pop('options', {})) + config.set_value(kwds) + targets = config.targets + if targets is None: + targets = (model, ) + + transBlocks = {} + visitor = LogicalToDisjunctiveVisitor() + for t in targets: + if t.ctype is Block or isinstance(t, _BlockData): + self._transform_block(t, model, visitor, transBlocks) + elif t.ctype is LogicalConstraint: + if t.is_indexed(): + self._transform_constraint(t, visitor, transBlocks) + else: + self._transform_constraintData(t, visitor, transBlocks) + else: + raise RuntimeError("Target '%s' was not a Block, Disjunct, or" + " LogicalConstraint. It was of type %s " + "and can't be transformed." % (t.name, + type(t))) + + def _transform_constraint(self, constraint, visitor, transBlocks): + for i in constraint.keys(ordered=True): + self._transform_constraintData(constraint[i], visitor, transBlocks) + constraint.deactivate() + + def _transform_block(self, target_block, model, new_varlists, transBlocks): + _blocks = target_block.values() if target_block.is_indexed() else \ + (target_block,) + for block in _blocks: + # Note that this changes the current (though not the original) + # behavior of logical-to-linear because we descend into Disjuncts in + # order to find logical constraints. In the context of creating a + # traditional disjunctive program, this makes sense--we cannot have + # logical constraints *anywhere* in the active tree after this + # transformation. + for logical_constraint in block.component_objects( + ctype=LogicalConstraint, + active=True, + descend_into=(Block, Disjunct)): + self._transform_constraint(logical_constraint, new_varlists, + transBlocks) + + def _transform_constraintData(self, logical_constraint, visitor, + transBlocks): + # now create a transformation block on the constraint's parent block (if + # we don't have one already) + parent_block = logical_constraint.parent_block() + xfrm_block = transBlocks.get(parent_block) + if xfrm_block is None: + xfrm_block = self._create_transformation_block(parent_block) + transBlocks[parent_block] = xfrm_block + + # This is may be too cute, but just deceive the walker so it puts stuff + # in the right place. + visitor.constraints = xfrm_block.transformed_constraints + visitor.z_vars = xfrm_block.auxiliary_vars + visitor.disjuncts = xfrm_block.auxiliary_disjuncts + visitor.disjunctions = xfrm_block.auxiliary_disjunctions + visitor.walk_expression(logical_constraint.expr) + logical_constraint.deactivate() + + def _create_transformation_block(self, context): + new_xfrm_block_name = unique_component_name( + context, '_logical_to_disjunctive') + new_xfrm_block = Block( + doc="Transformation objects for logical_to_disjunctive") + context.add_component(new_xfrm_block_name, new_xfrm_block) + + new_xfrm_block.transformed_constraints = ConstraintList() + new_xfrm_block.auxiliary_vars = VarList(domain=Binary) + new_xfrm_block.auxiliary_disjuncts = Disjunct(NonNegativeIntegers) + new_xfrm_block.auxiliary_disjunctions = Disjunction(NonNegativeIntegers) + + return new_xfrm_block diff --git a/pyomo/contrib/cp/transform/logical_to_disjunctive_walker.py b/pyomo/contrib/cp/transform/logical_to_disjunctive_walker.py new file mode 100644 index 00000000000..592adf3dbcb --- /dev/null +++ b/pyomo/contrib/cp/transform/logical_to_disjunctive_walker.py @@ -0,0 +1,238 @@ +# ___________________________________________________________________________ +# +# Pyomo: Python Optimization Modeling Objects +# Copyright (c) 2008-2022 +# National Technology and Engineering Solutions of Sandia, LLC +# Under the terms of Contract DE-NA0003525 with National Technology and +# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain +# rights in this software. +# This software is distributed under the 3-clause BSD License. +# ___________________________________________________________________________ + +import collections + +from pyomo.common.collections import ComponentMap +from pyomo.common.errors import MouseTrap +from pyomo.core.expr.expr_common import ExpressionType +from pyomo.core.expr.visitor import StreamBasedExpressionVisitor +from pyomo.core.expr.numeric_expr import NumericExpression +from pyomo.core.expr.relational_expr import RelationalExpression +import pyomo.core.expr.current as EXPR +import pyomo.core.expr.logical_expr as LE +from pyomo.core.base import ( + Binary, Constraint, ConstraintList, NonNegativeIntegers, VarList, value) +import pyomo.core.base.boolean_var as BV +from pyomo.core.base.expression import ScalarExpression, _GeneralExpressionData +from pyomo.core.base.param import ScalarParam, _ParamData +from pyomo.core.base.var import ScalarVar, _GeneralVarData +from pyomo.gdp.disjunct import AutoLinkedBooleanVar, Disjunct, Disjunction + +def _dispatch_boolean_var(visitor, node): + if node not in visitor.boolean_to_binary_map: + binary = node.get_associated_binary() + if binary is not None: + visitor.boolean_to_binary_map[node] = binary + else: + z = visitor.z_vars.add() + visitor.boolean_to_binary_map[node] = z + node.associate_binary_var(z) + return False, visitor.boolean_to_binary_map[node] + +def _dispatch_var(visitor, node): + return False, node + +def _dispatch_param(visitor, node): + if int(value(node)) == value(node): + return False, node + else: + raise ValueError("Found non-integer valued Param '%s' in a logical " + "expression. This cannot be written to a disjunctive " + "form." % node.name) + +def _dispatch_expression(visitor, node): + return False, node.expr + +def _before_relational_expr(visitor, node): + raise MouseTrap("The RelationalExpression '%s' was used as a Boolean term " + "in a logical proposition. This is not yet supported " + "when transforming to disjunctive form." % node) + +def _dispatch_not(visitor, node, a): + # z == !a + if a not in visitor.expansions: + z = visitor.z_vars.add() + visitor.constraints.add(z == 1 - a) + visitor.expansions[a] = z + return visitor.expansions[a] + +def _dispatch_implication(visitor, node, a, b): + # z == !a v b + return _dispatch_or(visitor, node, 1 - a, b) + +def _dispatch_equivalence(visitor, node, a, b): + # z == (!a v b) ^ (a v !b) + return _dispatch_and( + visitor, node, + _dispatch_or(visitor, node, 1 - a, b), + _dispatch_or(visitor, node, a, 1 - b), + ) + +def _dispatch_and(visitor, node, *args): + # z == a ^ b ^ ... + z = visitor.z_vars.add() + for arg in args: + visitor.constraints.add(arg >= z) + return z + +def _dispatch_or(visitor, node, *args): + # z == a v b v ... + # (!z v a v b v ...) ^ (z v !a) ^ (z v !b) ^ ... + z = visitor.z_vars.add() + visitor.constraints.add((1 - z) + sum(args) >= 1) + for arg in args: + visitor.constraints.add(z + (1 - arg) >= 1) + return z + +def _dispatch_xor(visitor, node, a, b): + # z == a XOR b + # This is a special case of exactly + return _dispatch_exactly(visitor, node, 1, a, b) + +def _get_integer_value(n, node): + if n.__class__ in EXPR.native_numeric_types and int(n) == n: + return n + if n.__class__ not in EXPR.native_types: + if n.is_potentially_variable(): + # [ESJ 11/22]: This is probably worth supporting sometime, but right + # now we are abiding by what docplex allows in their 'count' + # function. Part of supporting this will be making sure we catch + # strict inequalities in the GDP transformations. Because if we + # don't know that n is integer-valued we will be forced to write + # strict inequalities instead of incrememting or decrementing by 1 + # in the disjunctions. + raise MouseTrap( + "The first argument '%s' to '%s' is potentially variable. " + "This may be a mathematically coherent expression; However " + "it is not yet supported to convert it to a disjunctive " + "program." % (n, node)) + else: + return n + raise ValueError("The first argument to '%s' must be an integer.\n\t" + "Recieved: %s" % (node, n)) + +def _dispatch_exactly(visitor, node, *args): + # z = sum(args[1:] == args[0] + # This is currently implemented as: + # [sum(args[1:] = n] v [[sum(args[1:]) < n] v [sum(args[1:]) > n]] + M = len(args) - 1 + n = _get_integer_value(args[0], node) + sum_expr = sum(args[1:]) + equality_disj = visitor.disjuncts[len(visitor.disjuncts)] + equality_disj.constraint = Constraint(expr=sum_expr == n) + inequality_disj = visitor.disjuncts[len(visitor.disjuncts)] + inequality_disj.disjunction = Disjunction( + expr=[[sum_expr <= n - 1], [sum_expr >= n + 1]]) + visitor.disjunctions[len(visitor.disjunctions)] = [equality_disj, + inequality_disj] + return equality_disj.indicator_var.get_associated_binary() + +def _dispatch_atleast(visitor, node, *args): + # z = sum[args[1:] >= n + # This is implemented as: + # [sum(args[1:] >= n] v [sum(args[1:] < n] + n = _get_integer_value(args[0], node) + sum_expr = sum(args[1:]) + atleast_disj = visitor.disjuncts[len(visitor.disjuncts)] + less_disj = visitor.disjuncts[len(visitor.disjuncts)] + atleast_disj.constraint = Constraint(expr=sum_expr >= n) + less_disj.constraint = Constraint(expr=sum_expr <= n - 1) + visitor.disjunctions[len(visitor.disjunctions)] = [atleast_disj, less_disj] + return atleast_disj.indicator_var.get_associated_binary() + +def _dispatch_atmost(visitor, node, *args): + # z = sum[args[1:] <= n + # This is implemented as: + # [sum(args[1:] <= n] v [sum(args[1:] > n] + n = _get_integer_value(args[0], node) + sum_expr = sum(args[1:]) + atmost_disj = visitor.disjuncts[len(visitor.disjuncts)] + more_disj = visitor.disjuncts[len(visitor.disjuncts)] + atmost_disj.constraint = Constraint(expr=sum_expr <= n) + more_disj.constraint = Constraint(expr=sum_expr >= n + 1) + visitor.disjunctions[len(visitor.disjunctions)] = [atmost_disj, more_disj] + return atmost_disj.indicator_var.get_associated_binary() + +_operator_dispatcher = {} +_operator_dispatcher[LE.ImplicationExpression] = _dispatch_implication +_operator_dispatcher[LE.EquivalenceExpression] = _dispatch_equivalence +_operator_dispatcher[LE.NotExpression] = _dispatch_not +_operator_dispatcher[LE.AndExpression] = _dispatch_and +_operator_dispatcher[LE.OrExpression] = _dispatch_or +_operator_dispatcher[LE.XorExpression] = _dispatch_xor +_operator_dispatcher[LE.ExactlyExpression] = _dispatch_exactly +_operator_dispatcher[LE.AtLeastExpression] = _dispatch_atleast +_operator_dispatcher[LE.AtMostExpression] = _dispatch_atmost + +_before_child_dispatcher = {} +_before_child_dispatcher[BV.ScalarBooleanVar] = _dispatch_boolean_var +_before_child_dispatcher[BV._GeneralBooleanVarData] = _dispatch_boolean_var +_before_child_dispatcher[AutoLinkedBooleanVar] = _dispatch_boolean_var +_before_child_dispatcher[_ParamData] = _dispatch_param +_before_child_dispatcher[ScalarParam] = _dispatch_param +# for the moment, these are all just so we can get good error messages when we +# don't handle them: +_before_child_dispatcher[ScalarVar] = _dispatch_var +_before_child_dispatcher[_GeneralVarData] = _dispatch_var +_before_child_dispatcher[_GeneralExpressionData] = _dispatch_expression +_before_child_dispatcher[ScalarExpression] = _dispatch_expression + +class LogicalToDisjunctiveVisitor(StreamBasedExpressionVisitor): + """Converts BooleanExpressions to Linear (MIP) representation + + This converter eschews conjunctive normal form, and instead follows + the well-trodden MINLP path of factorable programming. + + """ + + def __init__(self): + super().__init__() + self.z_vars = VarList(domain=Binary) + self.z_vars.construct() + self.constraints = ConstraintList() + self.disjuncts = Disjunct(NonNegativeIntegers, concrete=True) + self.disjunctions = Disjunction(NonNegativeIntegers) + self.disjunctions.construct() + self.expansions = ComponentMap() + self.boolean_to_binary_map = ComponentMap() + + def initializeWalker(self, expr): + walk, result = self.beforeChild(None, expr, 0) + if not walk: + return False, self.finalizeResult(result) + return True, expr + + def beforeChild(self, node, child, child_idx): + if child.__class__ in EXPR.native_types: + return False, child + + if child.is_numeric_type(): + # Just pass it through, we'll figure it out later + return False, child + if child.is_expression_type(ExpressionType.RELATIONAL): + # Eventually we'll handle these. Right now we set a MouseTrap + return _before_relational_expr(self, child) + + if not child.is_expression_type() or child.is_named_expression_type(): + return _before_child_dispatcher[child.__class__](self, child) + + return True, None + + def exitNode(self, node, data): + return _operator_dispatcher[node.__class__](self, node, *data) + + def finalizeResult(self, result): + # This LogicalExpression must evaluate to True (but note that we cannot + # fix this variable to 1 since this logical expression could be living + # on a Disjunct and later need to be relaxed.) + self.constraints.add(result >= 1) + return result diff --git a/pyomo/contrib/gdpopt/branch_and_bound.py b/pyomo/contrib/gdpopt/branch_and_bound.py index f5180027667..358e54e7b51 100644 --- a/pyomo/contrib/gdpopt/branch_and_bound.py +++ b/pyomo/contrib/gdpopt/branch_and_bound.py @@ -88,12 +88,12 @@ def _solve_gdp(self, model, config): add_boolean_variable_lists(util_block) root_node = TransformationFactory( - 'core.logical_to_linear').create_using(model) + 'contrib.logical_to_disjunctive').create_using(model) root_util_blk = root_node.component( self.original_util_block.name) # Add to root utility block what we will need during the algorithm add_disjunction_list(root_util_blk) - # Now that logical_to_linear has been called. + # Now that logical_to_disjunctive has been called. add_transformed_boolean_variable_list(root_util_blk) # Map unfixed disjunct -> list of deactivated constraints diff --git a/pyomo/contrib/gdpopt/create_oa_subproblems.py b/pyomo/contrib/gdpopt/create_oa_subproblems.py index 0816d92213f..110ad99b1f7 100644 --- a/pyomo/contrib/gdpopt/create_oa_subproblems.py +++ b/pyomo/contrib/gdpopt/create_oa_subproblems.py @@ -165,16 +165,17 @@ def add_boolean_variable_lists(util_block): util_block.non_indicator_boolean_variable_list.append(v) # For the discrete problem, we want the corresponding binaries for all of the -# BooleanVars. This must be called after logical_to_linear has been called. +# BooleanVars. This must be called after logical_to_disjunctive has been called. def add_transformed_boolean_variable_list(util_block): util_block.transformed_boolean_variable_list = [ v.get_associated_binary() for v in util_block.boolean_variable_list] def get_subproblem(original_model, util_block): """Clone the original, and reclassify all the Disjuncts to Blocks. - We'll also call logical_to_linear in case any of the indicator_vars are - used in logical constraints and to make sure that the rest of the model is - algebraic (assuming it was a proper GDP to begin with). + We'll also call logical_to_disjunctive and bigm the disjunctive parts in + case any of the indicator_vars are used in logical constraints and to make + sure that the rest of the model is algebraic (assuming it was a proper + GDP to begin with). """ subproblem = original_model.clone() subproblem.name = subproblem.name + ": subproblem" @@ -203,7 +204,9 @@ def get_subproblem(original_model, util_block): disjunction.deactivate() - TransformationFactory('core.logical_to_linear').apply_to(subproblem) + TransformationFactory('contrib.logical_to_disjunctive').apply_to(subproblem) + # transform any of the Disjuncts we created above with bigm. + TransformationFactory('gdp.bigm').apply_to(subproblem) subproblem_util_block = subproblem.component(util_block.local_name) save_initial_values(subproblem_util_block) diff --git a/pyomo/contrib/gdpopt/solve_discrete_problem.py b/pyomo/contrib/gdpopt/solve_discrete_problem.py index 697379518be..c4f6288dca4 100644 --- a/pyomo/contrib/gdpopt/solve_discrete_problem.py +++ b/pyomo/contrib/gdpopt/solve_discrete_problem.py @@ -43,8 +43,9 @@ def solve_MILP_discrete_problem(util_block, solver, config): # are okay to leave in because if you tighten the bounds now, they # could only get tighter in later iterations, since you are # tightening this relaxation - except InfeasibleConstraintException: - config.logger.debug("MIP preprocessing detected infeasibility.") + except InfeasibleConstraintException as e: + config.logger.debug( + "MIP preprocessing detected infeasibility:\n\t%s" % str(e)) return tc.infeasible # Deactivate extraneous IMPORT/EXPORT suffixes diff --git a/pyomo/contrib/gdpopt/tests/test_gdpopt.py b/pyomo/contrib/gdpopt/tests/test_gdpopt.py index b21b25955cf..feebe94a4b8 100644 --- a/pyomo/contrib/gdpopt/tests/test_gdpopt.py +++ b/pyomo/contrib/gdpopt/tests/test_gdpopt.py @@ -31,7 +31,6 @@ from pyomo.contrib.mcpp.pyomo_mcpp import mcpp_available from pyomo.contrib.gdpopt.solve_discrete_problem import ( solve_MILP_discrete_problem, distinguish_mip_infeasible_or_unbounded) -from pyomo.core.expr.sympy_tools import sympy_available from pyomo.environ import ( Block, ConcreteModel, Constraint, Integers, LogicalConstraint, maximize, Objective, RangeSet, TransformationFactory, SolverFactory, sqrt, value, Var) @@ -605,18 +604,18 @@ def test_reverse_numeric_differentiation_in_LOA(self): self.assertAlmostEqual(results.problem.upper_bound, 1300*x_val, places=6) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_logical_constraints_on_disjuncts(self): m = models.makeLogicalConstraintsOnDisjuncts() SolverFactory('gdpopt.loa').solve(m, mip_solver=mip_solver, nlp_solver=nlp_solver) self.assertAlmostEqual(value(m.x), 8) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_logical_constraints_on_disjuncts_nonlinear_convex(self): m = models.makeLogicalConstraintsOnDisjuncts_NonlinearConvex() - SolverFactory('gdpopt.loa').solve(m, mip_solver=mip_solver, - nlp_solver=nlp_solver, tee=True) + results = SolverFactory('gdpopt.loa').solve(m, mip_solver=mip_solver, + nlp_solver=nlp_solver) + self.assertEqual(results.solver.termination_condition, + TerminationCondition.optimal) self.assertAlmostEqual(value(m.x), 4) def test_nested_disjunctions_no_init(self): @@ -635,7 +634,6 @@ def test_nested_disjunctions_max_binary(self): self.assertAlmostEqual(value(m.x), sqrt(2)/2) self.assertAlmostEqual(value(m.y), sqrt(2)/2) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_boolean_vars_on_disjuncts(self): m = models.makeBooleanVarsOnDisjuncts() SolverFactory('gdpopt.loa').solve(m, mip_solver=mip_solver, @@ -682,7 +680,6 @@ def test_time_limit(self): self.assertEqual(results.solver.termination_condition, TerminationCondition.maxTimeLimit) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_LOA_8PP_logical_default_init(self): """Test logic-based outer approximation with 8PP.""" exfile = import_file( @@ -728,7 +725,6 @@ def test_LOA_strip_pack_default_init(self): self.assertTrue( fabs(value(strip_pack.total_length.expr) - 11) <= 1E-2) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_LOA_strip_pack_logical_constraints(self): """Test logic-based outer approximation with variation of strip packing with some logical constraints.""" @@ -780,7 +776,6 @@ def test_LOA_8PP_maxBinary(self): nlp_solver=nlp_solver) ct.check_8PP_solution(self, eight_process, results) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_LOA_8PP_logical_maxBinary(self): """Test logic-based OA with max_binary initialization.""" exfile = import_file( @@ -804,7 +799,6 @@ def test_LOA_strip_pack_maxBinary(self): self.assertTrue( fabs(value(strip_pack.total_length.expr) - 11) <= 1E-2) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_LOA_strip_pack_maxBinary_logical_constraints(self): """Test LOA with strip packing using max_binary initialization and logical constraints.""" @@ -999,14 +993,12 @@ def test_GDP_nonlinear_objective(self): nlp_solver=nlp_solver) self.assertAlmostEqual(value(m.o), 0) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_logical_constraints_on_disjuncts(self): m = models.makeLogicalConstraintsOnDisjuncts() SolverFactory('gdpopt.ric').solve(m, mip_solver=mip_solver, nlp_solver=nlp_solver) self.assertAlmostEqual(value(m.x), 8) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_boolean_vars_on_disjuncts(self): m = models.makeBooleanVarsOnDisjuncts() SolverFactory('gdpopt.ric').solve(m, mip_solver=mip_solver, @@ -1024,7 +1016,6 @@ def test_RIC_8PP_default_init(self): tee=False) ct.check_8PP_solution(self, eight_process, results) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_RIC_8PP_logical_default_init(self): """Test logic-based outer approximation with 8PP.""" exfile = import_file( @@ -1070,7 +1061,6 @@ def test_RIC_strip_pack_default_init(self): self.assertTrue( fabs(value(strip_pack.total_length.expr) - 11) <= 1E-2) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_RIC_strip_pack_default_init_logical_constraints(self): """Test logic-based outer approximation with strip packing with logical constraints.""" @@ -1130,7 +1120,6 @@ def test_RIC_strip_pack_maxBinary(self): self.assertTrue( fabs(value(strip_pack.total_length.expr) - 11) <= 1E-2) - @unittest.skipUnless(sympy_available, "Sympy not available") def test_RIC_strip_pack_maxBinary_logical_constraints(self): """Test RIC with strip packing using max_binary initialization and including logical constraints.""" @@ -1348,18 +1337,16 @@ def test_GDP_integer_vars_infeasible(self): self.assertEqual(res.solver.termination_condition, TerminationCondition.infeasible) - @unittest.skipUnless(license_available and sympy_available, - "Global NLP solver license not available or sympy " - "not available.") + @unittest.skipUnless(license_available, + "Global NLP solver license not available") def test_logical_constraints_on_disjuncts(self): m = models.makeLogicalConstraintsOnDisjuncts() SolverFactory('gdpopt.gloa').solve(m, mip_solver=mip_solver, nlp_solver=global_nlp_solver) self.assertAlmostEqual(value(m.x), 8) - @unittest.skipUnless(license_available and sympy_available, - "Global NLP solver license not available or sympy " - "not available.") + @unittest.skipUnless(license_available, + "Global NLP solver license not available") def test_boolean_vars_on_disjuncts(self): m = models.makeBooleanVarsOnDisjuncts() SolverFactory('gdpopt.gloa').solve(m, mip_solver=mip_solver, @@ -1381,9 +1368,8 @@ def test_GLOA_8PP(self): ) ct.check_8PP_solution(self, eight_process, results) - @unittest.skipUnless(license_available and sympy_available, - "Global NLP solver license not available or sympy " - "not available.") + @unittest.skipUnless(license_available, + "Global NLP solver license not available") def test_GLOA_8PP_logical(self): """Test the global logic-based outer approximation algorithm.""" exfile = import_file( @@ -1428,9 +1414,8 @@ def test_GLOA_strip_pack_default_init(self): self.assertTrue( fabs(value(strip_pack.total_length.expr) - 11) <= 1E-2) - @unittest.skipUnless(license_available and sympy_available, - "Global NLP solver license not available or sympy " - "not available.") + @unittest.skipUnless(license_available, + "Global NLP solver license not available") def test_GLOA_strip_pack_default_init_logical_constraints(self): """Test logic-based outer approximation with strip packing.""" exfile = import_file( diff --git a/pyomo/core/expr/compare.py b/pyomo/core/expr/compare.py index e619d69fea9..d1c9408a12f 100644 --- a/pyomo/core/expr/compare.py +++ b/pyomo/core/expr/compare.py @@ -270,7 +270,9 @@ def assertExpressionsStructurallyEqual(test, a, b, include_named_exprs=True): try: test.assertEqual(len(prefix_a), len(prefix_b)) for _a, _b in zip(prefix_a, prefix_b): - test.assertIs(_a.__class__, _b.__class__) + if _a.__class__ not in native_types and \ + _b.__class__ not in native_types: + test.assertIs(_a.__class__, _b.__class__) test.assertEqual(_a, _b) except (PyomoException, AssertionError): test.fail(f"Expressions not structurally equal:\n\t" diff --git a/pyomo/core/plugins/transform/logical_to_linear.py b/pyomo/core/plugins/transform/logical_to_linear.py index dd0927c7353..17e0b5403f6 100644 --- a/pyomo/core/plugins/transform/logical_to_linear.py +++ b/pyomo/core/plugins/transform/logical_to_linear.py @@ -26,8 +26,9 @@ from pyomo.core.plugins.transform.hierarchy import IsomorphicTransformation from pyomo.core.util import target_list -@TransformationFactory.register("core.logical_to_linear", - doc="Convert logic to linear constraints") +@TransformationFactory.register( + "core.logical_to_linear", + doc="Convert logic to linear constraints") class LogicalToLinear(IsomorphicTransformation): """ Re-encode logical constraints as linear constraints, diff --git a/pyomo/core/tests/unit/test_logical_constraint.py b/pyomo/core/tests/unit/test_logical_constraint.py index cb9d7dd3be8..2702706c91a 100644 --- a/pyomo/core/tests/unit/test_logical_constraint.py +++ b/pyomo/core/tests/unit/test_logical_constraint.py @@ -43,12 +43,21 @@ def check_lor_on_disjunct(self, model, disjunct, x1, x2): self.assertEqual(repn.linear_coefs[1], 1) @unittest.skipUnless(sympy_available, "Sympy not available") - def test_statement_in_Disjunct(self): + def test_statement_in_Disjunct_with_logical_to_linear(self): + # This is an old test that originally tested that GDP's + # BigM/Hull correctly handled Disjuncts with LogicalConstraints + # (implicitly calling logical_to_linear to leave the transformed + # algebraic constraints on the Disjuncts). That is no longer + # the default behavior. However, we will preserve this (with an + # explicit call to logical_to_linear) for posterity model = self.create_model() model.disj = Disjunction(expr=[ [model.x.lor(model.y)], [model.y.lor(model.z)] ]) + TransformationFactory('core.logical_to_linear').apply_to( + model, targets=model.disj.disjuncts) + bigmed = TransformationFactory('gdp.bigm').create_using(model) # check that the algebraic versions are living on the Disjuncts self.check_lor_on_disjunct(bigmed, bigmed.disj.disjuncts[0], bigmed.x, diff --git a/pyomo/environ/__init__.py b/pyomo/environ/__init__.py index 2b716323b55..a0cf8986391 100644 --- a/pyomo/environ/__init__.py +++ b/pyomo/environ/__init__.py @@ -36,6 +36,7 @@ def _do_import(pkg_name): 'pyomo.contrib.ampl_function_demo', 'pyomo.contrib.appsi', 'pyomo.contrib.community_detection', + 'pyomo.contrib.cp', 'pyomo.contrib.example', 'pyomo.contrib.fme', 'pyomo.contrib.gdp_bounds', diff --git a/pyomo/gdp/plugins/bigm.py b/pyomo/gdp/plugins/bigm.py index c7c5ef00bc7..d097c440ae1 100644 --- a/pyomo/gdp/plugins/bigm.py +++ b/pyomo/gdp/plugins/bigm.py @@ -18,13 +18,13 @@ from pyomo.common.log import is_debug_set from pyomo.common.modeling import unique_component_name from pyomo.common.deprecation import deprecated, deprecation_warning +from pyomo.contrib.cp.transform.logical_to_disjunctive_program import ( + LogicalToDisjunctive) from pyomo.contrib.fbbt.fbbt import compute_bounds_on_expr from pyomo.core import ( Block, BooleanVar, Connector, Constraint, Param, Set, SetOf, Suffix, Var, - Expression, TraversalStrategy, value, RangeSet, + Expression, SortComponents, TraversalStrategy, value, RangeSet, NonNegativeIntegers, Binary, Any) -from pyomo.core.base.boolean_var import ( - _DeprecatedImplicitAssociatedBinaryVariable) from pyomo.core.base.external import ExternalFunction from pyomo.core.base import Transformation, TransformationFactory, Reference import pyomo.core.expr.current as EXPR @@ -210,6 +210,7 @@ def _apply_to_impl(self, instance, **kwds): # up. config.set_value(kwds) bigM = config.bigM + self.assume_fixed_vars_permanent = config.assume_fixed_vars_permanent targets = config.targets # We need to check that all the targets are in fact on instance. As we @@ -237,6 +238,23 @@ def _filter_inactive(targets): yield t targets = list(_filter_inactive(targets)) + # transform any logical constraints that might be anywhere on the stuff + # we're about to transform. We do this before we preprocess targets + # because we will likely create more disjunctive components that will + # need transformation. + disj_targets = [] + for t in targets: + disj_datas = t.values() if t.is_indexed() else [t,] + if t.ctype is Disjunct: + disj_targets.extend(disj_datas) + if t.ctype is Disjunction: + disj_targets.extend([d for disjunction in disj_datas for d in + disjunction.disjuncts]) + TransformationFactory('contrib.logical_to_disjunctive').apply_to( + instance, + targets=[blk for blk in targets if blk.ctype is Block] + + disj_targets) + # we need to preprocess targets to make sure that if there are any # disjunctions in targets that their disjuncts appear before them in # the list. @@ -245,15 +263,6 @@ def _filter_inactive(targets): knownBlocks, gdp_tree=gdp_tree) - # transform any logical constraints that might be anywhere on the stuff - # we're about to transform. - TransformationFactory('core.logical_to_linear').apply_to( - instance, - targets=[blk for blk in targets if blk.ctype is Block] + - [disj for disj in preprocessed_targets if disj.ctype is Disjunct]) - - self.assume_fixed_vars_permanent = config.assume_fixed_vars_permanent - for t in preprocessed_targets: if t.ctype is Disjunction: self._transform_disjunctionData( @@ -394,24 +403,6 @@ def _transform_disjunct(self, obj, bigM, root_disjunct): def _transform_block_components(self, block, disjunct, bigM, arg_list, suffix_list): - # We don't know where all the BooleanVars are going to be used, so if - # there are any that the logical_to_linear transformation didn't - # transform, we need to do it now, so that the Reference gets moved - # up. This won't be necessary when the writers are willing to find Vars - # not in the active subtree. - for boolean in block.component_data_objects(BooleanVar, - descend_into=Block, - active=None): - if isinstance(boolean._associated_binary, - _DeprecatedImplicitAssociatedBinaryVariable): - parent_block = boolean.parent_block() - new_var = Var(domain=Binary) - parent_block.add_component( - unique_component_name(parent_block, - boolean.local_name + "_asbinary"), - new_var) - boolean.associate_binary_var(new_var) - # Find all the variables declared here (including the indicator_var) and # add a reference on the transformation block so these will be # accessible when the Disjunct is deactivated. diff --git a/pyomo/gdp/plugins/fix_disjuncts.py b/pyomo/gdp/plugins/fix_disjuncts.py index 69768a43022..d00627b5cc9 100644 --- a/pyomo/gdp/plugins/fix_disjuncts.py +++ b/pyomo/gdp/plugins/fix_disjuncts.py @@ -4,7 +4,9 @@ import logging from math import fabs -from pyomo.common.config import ConfigBlock, NonNegativeFloat +from pyomo.common.config import ConfigDict, ConfigValue, NonNegativeFloat +from pyomo.contrib.cp.transform .logical_to_disjunctive_program import ( + LogicalToDisjunctive) from pyomo.core.base import Transformation, TransformationFactory from pyomo.core.base.block import Block from pyomo.core.expr.numvalue import value @@ -13,6 +15,13 @@ logger = logging.getLogger('pyomo.gdp.fix_disjuncts') +def _is_transformation(transformation_name): + xform = TransformationFactory(transformation_name) + if xform is None: + raise ValueError( + "Expected valid name for a registered Pyomo transformation. " + "\n\tRecieved: %s" % transformation_name) + return transformation_name @TransformationFactory.register( 'gdp.fix_disjuncts', @@ -35,7 +44,20 @@ def __init__(self, **kwargs): # standardized. super(GDP_Disjunct_Fixer, self).__init__(**kwargs) - CONFIG = ConfigBlock("gdp.fix_disjuncts") + CONFIG = ConfigDict("gdp.fix_disjuncts") + CONFIG.declare("GDP_to_MIP_transformation", ConfigValue( + default='gdp.bigm', + domain=_is_transformation, + description="The name of the transformation to call after the " + "'logical_to_disjunctive' transformation in order to finish " + "transforming to a MI(N)LP.", + doc=""" + If there are no logical constraints on the model being transformed, + this option is not used. However, if there are logical constraints + that involve mixtures of Boolean and integer variables, this option + specifies the transformation to use to transform the model with fixed + Disjuncts to a MI(N)LP. Uses 'gdp.bigm' by default. + """)) def _apply_to(self, model, **kwds): """Fix all disjuncts in the given model and reclassify them to Blocks.""" @@ -52,7 +74,9 @@ def _apply_to(self, model, **kwds): disjunct_object, Block) # Transform any remaining logical stuff - TransformationFactory('core.logical_to_linear').apply_to(model) + TransformationFactory('contrib.logical_to_disjunctive').apply_to(model) + # Transform anything disjunctive that the above created: + TransformationFactory(config.GDP_to_MIP_transformation).apply_to(model) def _transformContainer(self, obj): """Find all disjuncts in the container and transform them.""" diff --git a/pyomo/gdp/plugins/hull.py b/pyomo/gdp/plugins/hull.py index 7470b426e8a..eb52fab2307 100644 --- a/pyomo/gdp/plugins/hull.py +++ b/pyomo/gdp/plugins/hull.py @@ -23,8 +23,6 @@ Block, BooleanVar, Connector, Constraint, Param, Set, SetOf, Suffix, Var, Expression, SortComponents, TraversalStrategy, Any, RangeSet, Reals, value, NonNegativeIntegers, Binary ) -from pyomo.core.base.boolean_var import ( - _DeprecatedImplicitAssociatedBinaryVariable) from pyomo.gdp import Disjunct, Disjunction, GDP_Error from pyomo.gdp.transformed_disjunct import _TransformedDisjunct from pyomo.gdp.util import ( @@ -273,6 +271,23 @@ def _filter_inactive(targets): yield t targets = list(_filter_inactive(targets)) + # transform any logical constraints that might be anywhere on the stuff + # we're about to transform. We do this before we preprocess targets + # because we will likely create more disjunctive components that will + # need transformation. + disj_targets = [] + for t in targets: + disj_datas = t.values() if t.is_indexed() else [t,] + if t.ctype is Disjunct: + disj_targets.extend(disj_datas) + if t.ctype is Disjunction: + disj_targets.extend([d for disjunction in disj_datas for d in + disjunction.disjuncts]) + TransformationFactory('contrib.logical_to_disjunctive').apply_to( + instance, + targets=[blk for blk in targets if blk.ctype is Block] + + disj_targets) + # we need to preprocess targets to make sure that if there are any # disjunctions in targets that they appear before disjunctions that are # contained in their disjuncts. That is, in hull, we will transform from @@ -284,12 +299,6 @@ def _filter_inactive(targets): preprocessed_targets = gdp_tree.topological_sort() self._targets = set(preprocessed_targets) - # transform any logical constraints that might be anywhere on the stuff - # we're about to transform. - TransformationFactory('core.logical_to_linear').apply_to( - instance, - targets=[blk for blk in targets if blk.ctype is Block] + - [disj for disj in preprocessed_targets if disj.ctype is Disjunct]) for t in preprocessed_targets: if t.ctype is Disjunction: @@ -665,23 +674,6 @@ def _transform_disjunct(self, obj, transBlock, varSet, localVars, def _transform_block_components( self, block, disjunct, var_substitute_map, zero_substitute_map): - # We don't know where all the BooleanVars are used, so if there are any - # that the above transformation didn't transform, we need to do it now, - # so that the Reference gets moved up. This won't be necessary when the - # writers are willing to find Vars not in the active subtree. - for boolean in block.component_data_objects(BooleanVar, - descend_into=Block, - active=None): - if isinstance(boolean._associated_binary, - _DeprecatedImplicitAssociatedBinaryVariable): - parent_block = boolean.parent_block() - new_var = Var(domain=Binary) - parent_block.add_component( - unique_component_name(parent_block, - boolean.local_name + "_asbinary"), - new_var) - boolean.associate_binary_var(new_var) - # add references to all local variables on block (including the # indicator_var). Note that we do this after we have moved up the # transformation blocks for nested disjunctions, so that we don't have diff --git a/pyomo/gdp/plugins/multiple_bigm.py b/pyomo/gdp/plugins/multiple_bigm.py index f306aad8279..91e3c6f5615 100644 --- a/pyomo/gdp/plugins/multiple_bigm.py +++ b/pyomo/gdp/plugins/multiple_bigm.py @@ -23,8 +23,6 @@ Param, RangeSet, Set, SetOf, SortComponents, Suffix, value, Var ) from pyomo.core.base import Reference, Transformation, TransformationFactory -from pyomo.core.base.boolean_var import ( - _DeprecatedImplicitAssociatedBinaryVariable) import pyomo.core.expr.current as EXPR from pyomo.core.util import target_list @@ -228,6 +226,23 @@ def _apply_to_impl(self, instance, **kwds): if targets is None: targets = (instance, ) + # transform any logical constraints that might be anywhere on the stuff + # we're about to transform. We do this before we preprocess targets + # because we will likely create more disjunctive components that will + # need transformation. + disj_targets = [] + for t in targets: + disj_datas = t.values() if t.is_indexed() else [t,] + if t.ctype is Disjunct: + disj_targets.extend(disj_datas) + if t.ctype is Disjunction: + disj_targets.extend([d for disjunction in disj_datas for d in + disjunction.disjuncts]) + TransformationFactory('contrib.logical_to_disjunctive').apply_to( + instance, + targets=[blk for blk in targets if blk.ctype is Block] + + disj_targets) + # We don't allow nested, so it doesn't much matter which way we sort # this. But transforming from leaf to root makes the error checking for # complaining about nested smoother, so we do that. We have to transform @@ -236,13 +251,6 @@ def _apply_to_impl(self, instance, **kwds): gdp_tree = get_gdp_tree(targets, instance, knownBlocks) preprocessed_targets = gdp_tree.reverse_topological_sort() - # transform any logical constraints that might be anywhere on the stuff - # we're about to transform. - TransformationFactory('core.logical_to_linear').apply_to( - instance, - targets=[blk for blk in targets if blk.ctype is Block] + - [disj for disj in preprocessed_targets if disj.ctype is Disjunct]) - for t in preprocessed_targets: if t.ctype is Disjunction: self._transform_disjunctionData( @@ -366,21 +374,6 @@ def _transform_disjunct(self, obj, transBlock, active_disjuncts, Ms): obj._deactivate_without_fixing_indicator() def _transform_block_components(self, disjunct, active_disjuncts, Ms): - # We don't know where all the BooleanVars are used, so if there are any - # that logical_to_linear didn't transform, we need to do it now - for boolean in disjunct.component_data_objects(BooleanVar, - descend_into=Block, - active=None): - if isinstance(boolean._associated_binary, - _DeprecatedImplicitAssociatedBinaryVariable): - parent_block = boolean.parent_block() - new_var = Var(domain=Binary) - parent_block.add_component( - unique_component_name(parent_block, - boolean.local_name + "_asbinary"), - new_var) - boolean.associate_binary_var(new_var) - # add references to all local variables on block (including the # indicator_var). We won't have to do this when the writers can find # Vars not in the active subtree. diff --git a/pyomo/gdp/tests/common_tests.py b/pyomo/gdp/tests/common_tests.py index 3e9b9463503..291b63047a6 100644 --- a/pyomo/gdp/tests/common_tests.py +++ b/pyomo/gdp/tests/common_tests.py @@ -1715,7 +1715,7 @@ def check_pprint_equal(self, m, unpickle): self.assertMultiLineEqual(m_output, unpickle_output) def check_transformed_model_pickles(self, transformation): - # Do a model where we'll have to call logical_to_linear too. + # Do a model where we'll have to call logical_to_disjunctive too. m = models.makeLogicalConstraintsOnDisjuncts_NonlinearConvex() trans = TransformationFactory('gdp.%s' % transformation) trans.apply_to(m) diff --git a/pyomo/gdp/tests/test_bigm.py b/pyomo/gdp/tests/test_bigm.py index 4418f48428b..d7d83b78474 100644 --- a/pyomo/gdp/tests/test_bigm.py +++ b/pyomo/gdp/tests/test_bigm.py @@ -18,7 +18,8 @@ Any, value) from pyomo.gdp import Disjunct, Disjunction, GDP_Error from pyomo.core.base import constraint, _ConstraintData -from pyomo.core.expr.sympy_tools import sympy_available +from pyomo.core.expr.compare import ( + assertExpressionsEqual, assertExpressionsStructurallyEqual) from pyomo.repn import generate_standard_repn from pyomo.common.log import LoggingIntercept import logging @@ -2563,7 +2564,6 @@ def test_solution_maximize(self): def test_solution_minimize(self): ct.check_network_disjuncts(self, minimize=True, transformation='bigm') -@unittest.skipUnless(sympy_available, "Sympy not available") class LogicalConstraintsOnDisjuncts(unittest.TestCase): def test_logical_constraints_transformed(self): m = models.makeLogicalConstraintsOnDisjuncts() @@ -2577,56 +2577,169 @@ def test_logical_constraints_transformed(self): # first d[1]: cons = bigm.get_transformed_constraints( - m.d[1].logic_to_linear.transformed_constraints[1]) + m.d[1]._logical_to_disjunctive.transformed_constraints[1]) + # big-M transformation of z = 1 - y1: + # z <= 1 - y1 + (1 - d[1].indicator_var) + # z >= 1 - y1 - (1 - d[1].indicator_var) + z = m.d[1]._logical_to_disjunctive.auxiliary_vars[1] self.assertEqual(len(cons), 2) leq = cons[0] - self.assertEqual(leq.lower, 1) + self.assertEqual(leq.lower, 0) self.assertIsNone(leq.upper) repn = generate_standard_repn(leq.body) self.assertTrue(repn.is_linear()) - self.assertEqual(repn.constant, 2) - self.assertEqual(len(repn.linear_vars), 2) - ct.check_linear_coef(self, repn, y1, -1) - ct.check_linear_coef(self, repn, m.d[1].binary_indicator_var, -1) - # this is a stupid constraint + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + z + y1 - m.d[1].binary_indicator_var) geq = cons[1] + self.assertEqual(geq.upper, 0) self.assertIsNone(geq.lower) - self.assertEqual(geq.upper, 1) repn = generate_standard_repn(geq.body) self.assertTrue(repn.is_linear()) - self.assertEqual(repn.constant, 1) - self.assertEqual(len(repn.linear_vars), 1) - ct.check_linear_coef(self, repn, y1, -1) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + z + y1 + m.d[1].binary_indicator_var - 2) # then d[4]: - # 1 <= 1 - Y[2] + Y[1] + z1 = m.d[4]._logical_to_disjunctive.auxiliary_vars[1] + z2 = m.d[4]._logical_to_disjunctive.auxiliary_vars[2] + z3 = m.d[4]._logical_to_disjunctive.auxiliary_vars[3] # fixed True cons = bigm.get_transformed_constraints( - m.d[4].logic_to_linear.transformed_constraints[1]) + m.d[4]._logical_to_disjunctive.transformed_constraints[1]) self.assertEqual(len(cons), 1) - leq = cons[0] - self.assertEqual(leq.lower, 1) - self.assertIsNone(leq.upper) - repn = generate_standard_repn(leq.body) + c = cons[0] + # (1 - z1) + (1 - y1) + y2 >= 1 - (1 - d4.ind_var) + self.assertIsNone(c.upper) + self.assertEqual(c.lower, 1) + repn = generate_standard_repn(c.body) self.assertTrue(repn.is_linear()) - self.assertEqual(repn.constant, 2) - self.assertEqual(len(repn.linear_vars), 3) - ct.check_linear_coef(self, repn, y1, 1) - ct.check_linear_coef(self, repn, y2, -1) - ct.check_linear_coef(self, repn, m.d[4].binary_indicator_var, -1) - # 1 <= 1 - Y[1] + Y[2] + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + -z1 - y1 + y2 - m.d[4].binary_indicator_var + 3) cons = bigm.get_transformed_constraints( - m.d[4].logic_to_linear.transformed_constraints[2]) + m.d[4]._logical_to_disjunctive.transformed_constraints[2]) self.assertEqual(len(cons), 1) - leq = cons[0] - self.assertEqual(leq.lower, 1) - self.assertIsNone(leq.upper) - repn = generate_standard_repn(leq.body) + c = cons[0] + # z1 + 1 - (1 - y1) >= 1 - (1 - d4.ind_var) + self.assertIsNone(c.upper) + self.assertEqual(c.lower, 1) + repn = generate_standard_repn(c.body) self.assertTrue(repn.is_linear()) - self.assertEqual(repn.constant, 2) - self.assertEqual(len(repn.linear_vars), 3) - ct.check_linear_coef(self, repn, y2, 1) - ct.check_linear_coef(self, repn, y1, -1) - ct.check_linear_coef(self, repn, m.d[4].binary_indicator_var, -1) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + y1 + z1 - m.d[4].binary_indicator_var + 1) + cons = bigm.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[3]) + self.assertEqual(len(cons), 1) + c = cons[0] + # z1 + (1 - y2) >= 1 - (1 - d4.ind_var) + self.assertIsNone(c.upper) + self.assertEqual(c.lower, 1) + repn = generate_standard_repn(c.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + - y2 + z1 - m.d[4].binary_indicator_var + 2) + cons = bigm.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[4]) + self.assertEqual(len(cons), 1) + c = cons[0] + # (1 - z2) + y1 + (1 - y2) >= 1 - (1 - d4.ind_var) + self.assertIsNone(c.upper) + self.assertEqual(c.lower, 1) + repn = generate_standard_repn(c.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + -z2 - y2 + y1 - m.d[4].binary_indicator_var + 3) + cons = bigm.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[5]) + self.assertEqual(len(cons), 1) + c = cons[0] + # z2 + (1 - y1) >= 1 - (1 - d4.ind_var) + self.assertIsNone(c.upper) + self.assertEqual(c.lower, 1) + repn = generate_standard_repn(c.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + - y1 + z2 - m.d[4].binary_indicator_var + 2) + cons = bigm.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[6]) + self.assertEqual(len(cons), 1) + c = cons[0] + # z2 + 1 - (1 - y2) >= 1 - (1 - d4.ind_var) + self.assertIsNone(c.upper) + self.assertEqual(c.lower, 1) + repn = generate_standard_repn(c.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + y2 + z2 - m.d[4].binary_indicator_var + 1) + cons = bigm.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[7]) + self.assertEqual(len(cons), 1) + c = cons[0] + # z3 <= z1 + (1 - d4.ind_var) + self.assertIsNone(c.lower) + self.assertEqual(c.upper, 0) + repn = generate_standard_repn(c.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + z3 - z1 + m.d[4].binary_indicator_var - 1) + cons = bigm.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[8]) + self.assertEqual(len(cons), 1) + c = cons[0] + # z3 <= z2 + (1 - d4.ind_var) + self.assertIsNone(c.lower) + self.assertEqual(c.upper, 0) + repn = generate_standard_repn(c.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + z3 - z2 + m.d[4].binary_indicator_var - 1) # check that the global logical constraints were also transformed. self.assertFalse(m.p.active) diff --git a/pyomo/gdp/tests/test_hull.py b/pyomo/gdp/tests/test_hull.py index 163e7e49428..425578306e8 100644 --- a/pyomo/gdp/tests/test_hull.py +++ b/pyomo/gdp/tests/test_hull.py @@ -18,8 +18,8 @@ RealSet, ComponentMap, value, log, ConcreteModel, Any, Suffix, SolverFactory, RangeSet, Param, Objective, TerminationCondition, Reference) +from pyomo.core.expr.compare import assertExpressionsStructurallyEqual from pyomo.core.base import constraint -from pyomo.core.expr.sympy_tools import sympy_available from pyomo.repn import generate_standard_repn from pyomo.gdp import Disjunct, Disjunction, GDP_Error @@ -2141,7 +2141,6 @@ def test_solution_maximize(self): def test_solution_minimize(self): ct.check_network_disjuncts(self, minimize=True, transformation='hull') -@unittest.skipUnless(sympy_available, "Sympy not available") class LogicalConstraintsOnDisjuncts(unittest.TestCase): def test_logical_constraints_transformed(self): m = models.makeLogicalConstraintsOnDisjuncts() @@ -2155,25 +2154,50 @@ def test_logical_constraints_transformed(self): # first d[1]: cons = hull.get_transformed_constraints( - m.d[1].logic_to_linear.transformed_constraints[1]) + m.d[1]._logical_to_disjunctive.transformed_constraints[1]) + dis_z1 = hull.get_disaggregated_var( + m.d[1]._logical_to_disjunctive.auxiliary_vars[1], m.d[1]) + dis_y1 = hull.get_disaggregated_var(y1, m.d[1]) + self.assertEqual(len(cons), 1) # this simplifies because the dissaggregated variable is *always* 0 c = cons[0] + # hull transformation of z1 = 1 - y1: + # dis_z1 + dis_y1 = d[1].ind_var self.assertEqual(c.lower, 0) self.assertEqual(c.upper, 0) repn = generate_standard_repn(c.body) self.assertTrue(repn.is_linear()) - self.assertEqual(repn.constant, 0) - self.assertEqual(len(repn.linear_vars), 1) - ct.check_linear_coef(self, repn, hull.get_disaggregated_var(y1, m.d[1]), - -1) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + dis_z1 + dis_y1 - m.d[1].binary_indicator_var) + + cons = hull.get_transformed_constraints( + m.d[1]._logical_to_disjunctive.transformed_constraints[2]) + self.assertEqual(len(cons), 1) + c = cons[0] + # hull transformation of z1 >= 1 + assertExpressionsStructurallyEqual(self, c.expr, dis_z1 >= + m.d[1].binary_indicator_var) # then d[4]: y1d = hull.get_disaggregated_var(y1, m.d[4]) y2d = hull.get_disaggregated_var(y2, m.d[4]) - # 1 <= 1 - Y[2] + Y[1] + z1d = hull.get_disaggregated_var( + m.d[4]._logical_to_disjunctive.auxiliary_vars[1], m.d[4]) + z2d = hull.get_disaggregated_var( + m.d[4]._logical_to_disjunctive.auxiliary_vars[2], m.d[4]) + z3d = hull.get_disaggregated_var( + m.d[4]._logical_to_disjunctive.auxiliary_vars[3], m.d[4]) + + # hull transformation of (1 - z1) + (1 - y1) + y2 >= 1: + # dz1 + dy1 - dy2 <= m.d[4].ind_var cons = hull.get_transformed_constraints( - m.d[4].logic_to_linear.transformed_constraints[1]) + m.d[4]._logical_to_disjunctive.transformed_constraints[1]) # these also are simple because it's really an equality, and since both # disaggregated variables will be 0 when the disjunct isn't selected, it # doesn't even need big-Ming. @@ -2183,24 +2207,147 @@ def test_logical_constraints_transformed(self): self.assertEqual(cons.upper, 0) repn = generate_standard_repn(cons.body) self.assertTrue(repn.is_linear()) - self.assertEqual(repn.constant, 0) - self.assertEqual(len(repn.linear_vars), 2) - ct.check_linear_coef(self, repn, y2d, 1) - ct.check_linear_coef(self, repn, y1d, -1) - - # 1 <= 1 - Y[1] + Y[2] + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + - m.d[4].binary_indicator_var + z1d + y1d - y2d) + + # hull transformation of z1 + 1 - (1 - y1) >= 1 + # -y1d - z1d <= -d[4].ind_var cons = hull.get_transformed_constraints( - m.d[4].logic_to_linear.transformed_constraints[2]) + m.d[4]._logical_to_disjunctive.transformed_constraints[2]) self.assertEqual(len(cons), 1) cons = cons[0] self.assertIsNone(cons.lower) self.assertEqual(cons.upper, 0) repn = generate_standard_repn(cons.body) self.assertTrue(repn.is_linear()) - self.assertEqual(repn.constant, 0) - self.assertEqual(len(repn.linear_vars), 2) - ct.check_linear_coef(self, repn, y2d, -1) - ct.check_linear_coef(self, repn, y1d, 1) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + m.d[4].binary_indicator_var - y1d - z1d) + + # hull transformation of z1 + (1 - y2) >= 1 + # y2d - z1d <= 0 + cons = hull.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[3]) + self.assertEqual(len(cons), 1) + cons = cons[0] + self.assertIsNone(cons.lower) + self.assertEqual(cons.upper, 0) + repn = generate_standard_repn(cons.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + y2d - z1d) + + # hull transformation of (1 - z2) + y1 + (1 - y2) >= 1 + # z2d - y1d + y2d <= m.d[4].ind_var + cons = hull.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[4]) + self.assertEqual(len(cons), 1) + cons = cons[0] + self.assertIsNone(cons.lower) + self.assertEqual(cons.upper, 0) + repn = generate_standard_repn(cons.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + - m.d[4].binary_indicator_var + z2d + y2d - y1d) + + # hull transformation of z2 + (1 - y1) >= 1 + # y1d - z2d <= 0 + cons = hull.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[5]) + self.assertEqual(len(cons), 1) + cons = cons[0] + self.assertIsNone(cons.lower) + self.assertEqual(cons.upper, 0) + repn = generate_standard_repn(cons.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + y1d - z2d) + + # hull transformation of z2 + 1 - (1 - y2) >= 1 + # -y2d - z2d <= -d[4].ind_var + cons = hull.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[6]) + self.assertEqual(len(cons), 1) + cons = cons[0] + self.assertIsNone(cons.lower) + self.assertEqual(cons.upper, 0) + repn = generate_standard_repn(cons.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + m.d[4].binary_indicator_var - y2d - z2d) + + # hull transformation of z3 <= z1 + # z3d - z1d <= 0 + cons = hull.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[7]) + self.assertEqual(len(cons), 1) + cons = cons[0] + self.assertIsNone(cons.lower) + self.assertEqual(cons.upper, 0) + repn = generate_standard_repn(cons.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + z3d - z1d) + + # hull transformation of z3 <= z2 + # z3d - z2d <= 0 + cons = hull.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[8]) + self.assertEqual(len(cons), 1) + cons = cons[0] + self.assertIsNone(cons.lower) + self.assertEqual(cons.upper, 0) + repn = generate_standard_repn(cons.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + z3d - z2d) + + # hull transformation of z3 >= 1 + cons = hull.get_transformed_constraints( + m.d[4]._logical_to_disjunctive.transformed_constraints[9]) + self.assertEqual(len(cons), 1) + cons = cons[0] + assertExpressionsStructurallyEqual( + self, cons.expr, z3d >= m.d[4].binary_indicator_var) self.assertFalse(m.bwahaha.active) self.assertFalse(m.p.active) diff --git a/pyomo/gdp/tests/test_mbigm.py b/pyomo/gdp/tests/test_mbigm.py index 6d2208a5035..367f4c67037 100644 --- a/pyomo/gdp/tests/test_mbigm.py +++ b/pyomo/gdp/tests/test_mbigm.py @@ -16,7 +16,9 @@ from pyomo.common.fileutils import import_file, PYOMO_ROOT_DIR from pyomo.common.log import LoggingIntercept import pyomo.common.unittest as unittest -from pyomo.core.expr.compare import assertExpressionsEqual +from pyomo.core.expr.compare import ( + assertExpressionsEqual, assertExpressionsStructurallyEqual +) from pyomo.environ import ( BooleanVar, ConcreteModel, Constraint, LogicalConstraint, @@ -619,31 +621,73 @@ def test_logical_constraints_on_disjuncts(self): m = self.make_model() m.d1.Y = BooleanVar() m.d1.Z = BooleanVar() - # include an unused boolean var because it's our job to move this out of - # what will be deactivated. - m.d1.W = BooleanVar() m.d1.logical = LogicalConstraint(expr=m.d1.Y.implies(m.d1.Z)) mbm = TransformationFactory('gdp.mbigm') mbm.apply_to(m, bigM=self.get_Ms(m), reduce_bound_constraints=False) - transformed = mbm.get_transformed_constraints( - m.d1.logic_to_linear.transformed_constraints[1]) - self.assertEqual(len(transformed), 1) y = m.d1.Y.get_associated_binary() z = m.d1.Z.get_associated_binary() + z1 = m.d1._logical_to_disjunctive.auxiliary_vars[3] + + # MbigM transformation of: (1 - z1) + (1 - y) + z >= 1 + # (1 - z1) + (1 - y) + z >= 1 - d2.ind_var - d3.ind_var + transformed = mbm.get_transformed_constraints( + m.d1._logical_to_disjunctive.transformed_constraints[1]) + self.assertEqual(len(transformed), 1) c = transformed[0] check_obj_in_active_tree(self, c) self.assertIsNone(c.lower) self.assertEqual(value(c.upper), 0) repn = generate_standard_repn(c.body) self.assertTrue(repn.is_linear()) - self.assertEqual(repn.constant, 0) - self.assertEqual(len(repn.linear_vars), 4) - check_linear_coef(self, repn, y, 1) - check_linear_coef(self, repn, z, -1) - check_linear_coef(self, repn, m.d2.binary_indicator_var, -1) - check_linear_coef(self, repn, m.d3.binary_indicator_var, -1) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + - m.d2.binary_indicator_var - m.d3.binary_indicator_var + z1 + + y - z - 1) + + # MbigM transformation of: z1 + 1 - (1 - y) >= 1 + # z1 + y >= 1 - d2.ind_var - d3.ind_var + transformed = mbm.get_transformed_constraints( + m.d1._logical_to_disjunctive.transformed_constraints[2]) + self.assertEqual(len(transformed), 1) + c = transformed[0] + check_obj_in_active_tree(self, c) + self.assertIsNone(c.lower) + self.assertEqual(value(c.upper), 0) + repn = generate_standard_repn(c.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + - m.d2.binary_indicator_var - m.d3.binary_indicator_var - y - z1 + + 1) + + # MbigM transformation of: z1 + 1 - z >= 1 + # z1 + 1 - z >= 1 - d2.ind_var - d3.ind_var + transformed = mbm.get_transformed_constraints( + m.d1._logical_to_disjunctive.transformed_constraints[3]) + self.assertEqual(len(transformed), 1) + c = transformed[0] + check_obj_in_active_tree(self, c) + self.assertIsNone(c.lower) + self.assertEqual(value(c.upper), 0) + repn = generate_standard_repn(c.body) + self.assertTrue(repn.is_linear()) + simplified = repn.constant + sum( + repn.linear_coefs[i]*repn.linear_vars[i] + for i in range(len(repn.linear_vars))) + assertExpressionsStructurallyEqual( + self, + simplified, + - m.d2.binary_indicator_var - m.d3.binary_indicator_var + z - z1) def check_traditionally_bigmed_constraints(self, m, mbm, Ms): cons = mbm.get_transformed_constraints(m.d1.func) diff --git a/pyomo/gdp/tests/test_partition_disjuncts.py b/pyomo/gdp/tests/test_partition_disjuncts.py index ab840ac3302..8914f559afa 100644 --- a/pyomo/gdp/tests/test_partition_disjuncts.py +++ b/pyomo/gdp/tests/test_partition_disjuncts.py @@ -16,7 +16,6 @@ from pyomo.core.expr.logical_expr import ( EquivalenceExpression, NotExpression, AndExpression, ExactlyExpression ) -from pyomo.core.expr.sympy_tools import sympy_available from pyomo.gdp import Disjunct, Disjunction from pyomo.gdp.util import GDP_Error, check_model_algebraic from pyomo.gdp.plugins.partition_disjuncts import ( @@ -1872,7 +1871,6 @@ def test_untransformed_network_raises_GDPError(self): num_partitions=2) @unittest.skipUnless(ct.linear_solvers, "Could not find a linear solver") - @unittest.skipUnless(sympy_available, "Sympy not available") def test_network_disjuncts(self): ct.check_network_disjuncts(self, True, 'between_steps', num_partitions=2)