Modeling in Pyomo.GDP

Disjunctions

To demonstrate modeling with disjunctions in Pyomo.GDP, we revisit the small example from the previous page.

\[\begin{split}\left[\begin{gathered} Y_1 \\ \exp(x_2) - 1 = x_1 \\ x_3 = x_4 = 0 \end{gathered} \right] \bigvee \left[\begin{gathered} Y_2 \\ \exp\left(\frac{x_4}{1.2}\right) - 1 = x_3 \\ x_1 = x_2 = 0 \end{gathered} \right]\end{split}\]

Explicit syntax: more descriptive

Pyomo.GDP explicit syntax (see below) provides more clarity in the declaration of each modeling object, and gives the user explicit control over the Disjunct names. Assuming the ConcreteModel object m and variables have been defined, lines 1 and 5 declare the Disjunct objects corresponding to selection of unit 1 and 2, respectively. Lines 2 and 6 define the input-output relations for each unit, and lines 3-4 and 7-8 enforce zero flow through the unit that is not selected. Finally, line 9 declares the logical disjunction between the two disjunctive terms.

m.unit1 = Disjunct()
m.unit1.inout = Constraint(expr=exp(m.x[2]) - 1 == m.x[1])
m.unit1.no_unit2_flow1 = Constraint(expr=m.x[3] == 0)
m.unit1.no_unit2_flow2 = Constraint(expr=m.x[4] == 0)
m.unit2 = Disjunct()
m.unit2.inout = Constraint(expr=exp(m.x[4] / 1.2) - 1 == m.x[3])
m.unit2.no_unit1_flow1 = Constraint(expr=m.x[1] == 0)
m.unit2.no_unit1_flow2 = Constraint(expr=m.x[2] == 0)
m.use_unit1or2 = Disjunction(expr=[m.unit1, m.unit2])

The indicator variables for each disjunct \(Y_1\) and \(Y_2\) are automatically generated by Pyomo.GDP, accessible via m.unit1.indicator_var and m.unit2.indicator_var.

Compact syntax: more concise

For more advanced users, a compact syntax is also available below, taking advantage of the ability to declare disjuncts and constraints implicitly. When the Disjunction object constructor is passed a list of lists, the outer list defines the disjuncts and the inner list defines the constraint expressions associated with the respective disjunct.

m.use1or2 = Disjunction(expr=[
    # First disjunct
    [exp(m.x[2])-1 == m.x[1],
     m.x[3] == 0, m.x[4] == 0],
    # Second disjunct
    [exp(m.x[4]/1.2)-1 == m.x[3],
     m.x[1] == 0, m.x[2] == 0]])

Note

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, \(Y_1 \veebar Y_2\).

Logical Propositions

Pyomo.GDP also supports the use of logical propositions through the use of the BooleanVar and LogicalConstraint objects. The BooleanVar object in Pyomo represents Boolean variables, analogous to Var for numeric variables. BooleanVar can be indexed over a Pyomo Set, as below:

>>> m = ConcreteModel()
>>> m.my_set = RangeSet(4)
>>> m.Y = BooleanVar(m.my_set)
>>> m.Y.display()
Y : Size=4, Index=my_set
    Key : Value : Fixed : Stale
      1 :  None : False :  True
      2 :  None : False :  True
      3 :  None : False :  True
      4 :  None : False :  True

Using these Boolean variables, we can define LogicalConstraint objects, analogous to algebraic Constraint objects.

>>> m.p = LogicalConstraint(expr=m.Y[1].implies(m.Y[2] & m.Y[3]) | m.Y[4])
>>> m.p.pprint()
p : Size=1, Index=None, Active=True
    Key  : Body                          : Active
    None : (Y[1] --> Y[2] ∧ Y[3]) ∨ Y[4] :   True

Supported Logical Operators

Pyomo.GDP logical expression system supported operators and their usage are listed in the table below.

Operator	Operator	Method	Function
Negation	~Y[1]		lnot(Y[1])
Conjunction	Y[1] & Y[2]	Y[1].land(Y[2])	land(Y[1],Y[2])
Disjunction	Y[1] | Y[2]	Y[1].lor(Y[2])	lor(Y[1],Y[2])
Exclusive OR	Y[1] ^ Y[2]	Y[1].xor(Y[2])	xor(Y[1], Y[2])
Implication		Y[1].implies(Y[2])	implies(Y[1], Y[2])
Equivalence		Y[1].equivalent_to(Y[2])	equivalent(Y[1], Y[2])

Note

We omit support for some infix operators, e.g. Y[1] >> Y[2], due to concerns about non-intuitive Python operator precedence. That is Y[1] | Y[2] >> Y[3] would translate to \(Y_1 \lor (Y_2 \Rightarrow Y_3)\) rather than \((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.

Usage:

• atleast(3, Y[1], Y[2], Y[3])

• atmost(3, Y)

• exactly(3, Y)

>>> m = ConcreteModel()
>>> m.my_set = RangeSet(4)
>>> m.Y = BooleanVar(m.my_set)
>>> m.p = LogicalConstraint(expr=atleast(3, m.Y))
>>> m.p.pprint()
p : Size=1, Index=None, Active=True
    Key  : Body                                 : Active
    None : atleast(3: [Y[1], Y[2], Y[3], Y[4]]) :   True
>>> 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={1}, 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 on the next page.

Indexed logical constraints

Like Constraint objects for algebraic expressions, LogicalConstraint objects can be indexed. An example of this usage may be found below for the expression:

\[Y_{i+1} \Rightarrow Y_{i}, \quad i \in \{1, 2, \dots, n-1\}\]

>>> m = ConcreteModel()
>>> n = 5
>>> m.I = RangeSet(n)
>>> m.Y = BooleanVar(m.I)

>>> @m.LogicalConstraint(m.I)
... def p(m, i):
...     return m.Y[i+1].implies(m.Y[i]) if i < n else Constraint.Skip

>>> m.p.pprint()
p : Size=4, Index=I, Active=True
    Key : Body          : Active
      1 : Y[2] --> Y[1] :   True
      2 : Y[3] --> Y[2] :   True
      3 : Y[4] --> Y[3] :   True
      4 : Y[5] --> Y[4] :   True

Integration with Disjunctions

Note

Historically, the indicator_var on Disjunct objects was implemented as a binary Var. Beginning in Pyomo 6.0, that has been changed to the more mathematically correct BooleanVar, with the associated binary variable available as binary_indicator_var.

The logical expression system is designed to augment the previously introduced Disjunct and Disjunction components. Mathematically, the disjunct indicator variable is Boolean, and can be used directly in logical propositions.

Here, we demonstrate this capability with a toy example:

\[\begin{array}{ll} \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 \veebar Y_2 \\ & Y_3 \veebar Y_4 \\ & Y_1 \Rightarrow Y_4 \end{array}\]

>>> m = ConcreteModel()
>>> m.s = RangeSet(4)
>>> m.ds = RangeSet(2)
>>> m.d = Disjunct(m.s)
>>> m.djn = Disjunction(m.ds)
>>> m.djn[1] = [m.d[1], m.d[2]]
>>> m.djn[2] = [m.d[3], m.d[4]]
>>> m.x = Var(bounds=(-2, 10))
>>> m.d[1].c = Constraint(expr=m.x >= 2)
>>> m.d[2].c = Constraint(expr=m.x >= 3)
>>> m.d[3].c = Constraint(expr=m.x <= 8)
>>> m.d[4].c = Constraint(expr=m.x == 2.5)
>>> m.o = Objective(expr=m.x)

>>> # Add the logical proposition
>>> m.p = LogicalConstraint(
...    expr=m.d[1].indicator_var.implies(m.d[4].indicator_var))
>>> # Note: the implicit XOR enforced by m.djn[1] and m.djn[2] still apply

>>> # Apply the Big-M reformulation: It will convert the logical
>>> # propositions to algebraic expressions.
>>> TransformationFactory('gdp.bigm').apply_to(m)

>>> # Before solve, Boolean vars have no value
>>> Reference(m.d[:].indicator_var).display()
IndexedBooleanVar : Size=4, Index=s, ReferenceTo=d[:].indicator_var
    Key : Value : Fixed : Stale
      1 :  None : False :  True
      2 :  None : False :  True
      3 :  None : False :  True
      4 :  None : False :  True

>>> # Solve the reformulated model
>>> run_data = SolverFactory('glpk').solve(m)
>>> Reference(m.d[:].indicator_var).display()
IndexedBooleanVar : Size=4, Index=s, ReferenceTo=d[:].indicator_var
    Key : Value : Fixed : Stale
      1 :  True : False : False
      2 : False : False : False
      3 : False : False : False
      4 :  True : False : False

Advanced LogicalConstraint Examples

Support for complex nested expressions is a key benefit of the logical expression system. Below are examples of expressions that we support, and with some, an explanation of their implementation.

Composition of standard operators

\[Y_1 \vee Y_2 \implies Y_3 \wedge \neg Y_4 \wedge (Y_5 \vee Y_6)\]

m.p = LogicalConstraint(expr=(m.Y[1] | m.Y[2]).implies(
    m.Y[3] & ~m.Y[4] & (m.Y[5] | m.Y[6]))
)

Expressions within CP-type operators

\[\text{atleast}(3, Y_1, Y_2 \vee Y_3, Y_4 \Rightarrow Y_5, Y_6)\]

Here, augmented variables may be automatically added to the model as follows:

\[\begin{split}\text{atleast}(3, Y_1, Y_A, Y_B, Y_6)\\ Y_A \Leftrightarrow Y_2 \vee Y_3\\ Y_B \Leftrightarrow (Y_4 \Rightarrow Y_5)\end{split}\]

m.p = LogicalConstraint(
    expr=atleast(3, m.Y[1], Or(m.Y[2], m.Y[3]), m.Y[4].implies(m.Y[5]), m.Y[6]))

Nested CP-style operators

\[\text{atleast}(2, Y_1, \text{exactly}(2, Y_2, Y_3, Y_4), Y_5, Y_6)\]

Here, we again need to add augmented variables:

\[\begin{split}\text{atleast}(2, Y_1, Y_A, Y_5, Y_6)\\ Y_A \Leftrightarrow \text{exactly}(2, Y_2, Y_3, Y_4)\end{split}\]

However, we also need to further interpret the second statement as a disjunction:

\begin{gather*} \text{atleast}(2, Y_1, Y_A, Y_5, Y_6)\\ \left[\begin{gathered}Y_A\\\text{exactly}(2, Y_2, Y_3, Y_4)\end{gathered}\right] \vee \left[\begin{gathered}\neg Y_A\\ \left[\begin{gathered}Y_B\\\text{atleast}(3, Y_2, Y_3, Y_4)\end{gathered}\right] \vee \left[\begin{gathered}Y_C\\\text{atmost}(1, Y_2, Y_3, Y_4)\end{gathered}\right] \end{gathered}\right] \end{gather*}

or equivalently,

\begin{gather*} \text{atleast}(2, Y_1, Y_A, Y_5, Y_6)\\ \text{exactly}(1, Y_A, Y_B, Y_C)\\ \left[\begin{gathered}Y_A\\\text{exactly}(2, Y_2, Y_3, Y_4)\end{gathered}\right] \vee \left[\begin{gathered}Y_B\\\text{atleast}(3, Y_2, Y_3, Y_4)\end{gathered}\right] \vee \left[\begin{gathered}Y_C\\\text{atmost}(1, Y_2, Y_3, Y_4)\end{gathered}\right] \end{gather*}

m.p = LogicalConstraint(
    expr=atleast(2, m.Y[1], exactly(2, m.Y[2], m.Y[3], m.Y[4]), m.Y[5], m.Y[6]))

In the logical_to_linear transformation, we automatically convert these special disjunctions to linear form using a Big M reformulation.

Additional Examples

The following models all work and are equivalent for \(\left[x = 0\right] \veebar \left[y = 0\right]\):

Option 1: Rule-based construction

>>> import pyomo.environ as pyo
>>> from pyomo.gdp import Disjunct, Disjunction
>>> model = pyo.ConcreteModel()

>>> model.x = pyo.Var()
>>> model.y = pyo.Var()

>>> # Two conditions
>>> def _d(disjunct, flag):
...    model = disjunct.model()
...    if flag:
...       # x == 0
...       disjunct.c = pyo.Constraint(expr=model.x == 0)
...    else:
...       # y == 0
...       disjunct.c = pyo.Constraint(expr=model.y == 0)
>>> model.d = Disjunct([0,1], rule=_d)

>>> # Define the disjunction
>>> def _c(model):
...    return [model.d[0], model.d[1]]
>>> model.c = Disjunction(rule=_c)

Option 2: Explicit disjuncts

>>> import pyomo.environ as pyo
>>> from pyomo.gdp import Disjunct, Disjunction
>>> model = pyo.ConcreteModel()

>>> model.x = pyo.Var()
>>> model.y = pyo.Var()

>>> model.fix_x = Disjunct()
>>> model.fix_x.c = pyo.Constraint(expr=model.x == 0)

>>> model.fix_y = Disjunct()
>>> model.fix_y.c = pyo.Constraint(expr=model.y == 0)

>>> model.c = Disjunction(expr=[model.fix_x, model.fix_y])

Option 3: Implicit disjuncts (disjunction rule returns a list of
expressions or a list of lists of expressions)

>>> import pyomo.environ as pyo
>>> from pyomo.gdp import Disjunction
>>> model = pyo.ConcreteModel()

>>> model.x = pyo.Var()
>>> model.y = pyo.Var()

>>> model.c = Disjunction(expr=[model.x == 0, model.y == 0])