This basic example of the use of GDP (General Disjunctive Modeling) is taken from the pyomodoc. The GDPopt solver in Pyomo allows users to solve nonlinear Generalized Disjunctive Programming (GDP) models using logic-based decomposition approaches, as opposed to the conventional approach via reformulation to a Mixed Integer Nonlinear Programming (MINLP) model.
\[
\begin{aligned}
\min_{x,\, y} \quad & x + 0.1\, y \\[6pt]
\text{s.t.} \quad & x + y = 1 \\[4pt]
& \begin{bmatrix} x = 0 \end{bmatrix} \;\vee\; \begin{bmatrix} y = 0 \end{bmatrix} \\[4pt]
& -1.2 \leq x \leq 2 \\
& -10 \leq y \leq 10
\end{aligned}
\]
The disjunction \([\,x=0\,] \vee [\,y=0\,]\) is a GDP (Generalized Disjunctive Program) constraint: at least one of the two bracketed sub-problems must hold. Combined with the equality \(x + y = 1\), the two feasible cases are:
Active disjunct
\(x\)
\(y\)
Objective
\(x = 0\)
\(0\)
\(1\)
\(0 + 0.1(1) = \mathbf{0.1}\)
\(y = 0\)
\(1\)
\(0\)
\(1 + 0.1(0) = 1.0\)
So the optimal solution is \(x^* = 0,\ y^* = 1\) with objective value \({0.1}\), achieved by activating the first disjunct.
The LOA (Logic-based Outer Approximation) solver finds this by introducing binary indicator variables \(z_1, z_2 \in \{0,1\}\) with \(z_1 + z_2 \geq 1\) and enforcing the disjunct constraints only when the corresponding indicator is active.
Starting GDPopt version 22.5.13 using LOA algorithm
iterlim: None
time_limit: None
tee: true
logger: <Logger pyomo.contrib.gdpopt (INFO)>
integer_tolerance: 1.0e-05
constraint_tolerance: 1.0e-06
variable_tolerance: 1.0e-08
subproblem_initialization_method: <function restore_vars_to_original_values at 0x108c374c0>
call_before_subproblem_solve: <class 'type'>
call_after_subproblem_solve: <class 'type'>
call_after_subproblem_feasible: <class 'type'>
force_subproblem_nlp: false
subproblem_presolve: true
tighten_nlp_var_bounds: false
round_discrete_vars: true
max_fbbt_iterations: 3
init_strategy: None
init_algorithm: set_covering
custom_init_disjuncts: []
max_slack: 1000.0
OA_penalty_factor: 1000.0
set_cover_iterlim: 8
discrete_problem_transformation: gdp.bigm
call_before_discrete_problem_solve: <class 'type'>
call_after_discrete_problem_solve: <class 'type'>
call_before_master_solve: <class 'type'>
call_after_master_solve: <class 'type'>
mip_presolve: true
calc_disjunctive_bounds: false
obbt_disjunctive_bounds: false
mip_solver: glpk
mip_solver_args:
nlp_solver: ipopt
nlp_solver_args:
minlp_solver: baron
minlp_solver_args:
local_minlp_solver: bonmin
local_minlp_solver_args:
small_dual_tolerance: 1.0e-08
bound_tolerance: 1.0e-06
If you use this software, you may cite the following:
- Implementation:
Chen, Q; Johnson, ES; Bernal, DE; Valentin, R; Kale, S;
Bates, J; Siirola, JD; Grossmann, IE.
Pyomo.GDP: an ecosystem for logic based modeling and optimization
development.
Optimization and Engineering, 2021.
- LOA algorithm:
Türkay, M; Grossmann, IE.
Logic-based MINLP algorithms for the optimal synthesis of process
networks. Comp. and Chem. Eng. 1996, 20(8), 959–978.
DOI: 10.1016/0098-1354(95)00219-7.
Original model has 3 constraints (0 nonlinear) and 1 disjunctions, with 4 variables, of which 2 are binary, 0 are integer, and 2 are continuous.
---Starting GDPopt initialization---
Starting set covering initialization.
=============================================================================================
Iteration | Subproblem Type | Lower Bound | Upper Bound | Gap | Time(s)
Initialization complete.
Finished discrete problem initialization in 0.00s and 0 iterations
=============================================================================================
Iteration | Subproblem Type | Lower Bound | Upper Bound | Gap | Time(s)
1 discrete 0.10000 inf nan% 0.06
1 subproblem 0.10000 0.10000 0.00% 0.07 *
1 0.10000 0.10000 0.00% 0.07
GDPopt exiting--bounds have converged or crossed.
Solved in 1 iterations and 0.07433 seconds
Optimal objective value 0.1000000000
Relative optimality gap 0.00000%
Model LOA example
Variables:
x : Size=1, Index=None
Key : Lower : Value : Upper : Fixed : Stale : Domain
None : -1.2 : 0 : 2 : False : False : Reals
y : Size=1, Index=None
Key : Lower : Value : Upper : Fixed : Stale : Domain
None : -10 : 1 : 10 : False : False : Reals
Objectives:
objective : Size=1, Index=None, Active=True
Key : Active : Value
None : True : 0.1
Constraints:
c : Size=1
Key : Lower : Body : Upper
None : 1.0 : 1 : 1.0
model.pprint()
2 Var Declarations
x : Size=1, Index=None
Key : Lower : Value : Upper : Fixed : Stale : Domain
None : -1.2 : 0 : 2 : False : False : Reals
y : Size=1, Index=None
Key : Lower : Value : Upper : Fixed : Stale : Domain
None : -10 : 1 : 10 : False : False : Reals
1 Objective Declarations
objective : Size=1, Index=None, Active=True
Key : Active : Sense : Expression
None : True : minimize : x + 0.1*y
1 Constraint Declarations
c : Size=1, Index=None, Active=True
Key : Lower : Body : Upper : Active
None : 1.0 : x + y : 1.0 : True
2 Disjunct Declarations
fix_x : Size=1, Index=None, Active=True
1 Var Declarations
binary_indicator_var : Size=1, Index=None
Key : Lower : Value : Upper : Fixed : Stale : Domain
None : 0 : 1 : 1 : False : False : Binary
1 Constraint Declarations
c : Size=1, Index=None, Active=True
Key : Lower : Body : Upper : Active
None : 0.0 : x : 0.0 : True
1 BooleanVar Declarations
indicator_var : Size=1, Index=None
Key : Value : Fixed : Stale
None : True : False : False
3 Declarations: indicator_var binary_indicator_var c
fix_y : Size=1, Index=None, Active=True
1 Var Declarations
binary_indicator_var : Size=1, Index=None
Key : Lower : Value : Upper : Fixed : Stale : Domain
None : 0 : 0 : 1 : False : False : Binary
1 Constraint Declarations
c : Size=1, Index=None, Active=True
Key : Lower : Body : Upper : Active
None : 0.0 : y : 0.0 : True
1 BooleanVar Declarations
indicator_var : Size=1, Index=None
Key : Value : Fixed : Stale
None : False : False : False
3 Declarations: indicator_var binary_indicator_var c
1 Disjunction Declarations
d : Size=1, Index=None, Active=True
Key : Disjuncts : Active : XOR
None : ['fix_x', 'fix_y'] : True : True
7 Declarations: x y c fix_x fix_y d objective
results
{'Problem': [{'Name': 'LOA example', 'Lower bound': 0.1, 'Upper bound': 0.1, 'Number of objectives': 1, 'Number of constraints': 3, 'Number of variables': 4, 'Number of binary variables': 2, 'Number of integer variables': 0, 'Number of continuous variables': 2, 'Number of nonzeros': None, 'Sense': 'minimize', 'Number of disjunctions': 1}], 'Solver': [{'Name': 'GDPopt (22, 5, 13) - LOA', 'Status': 'ok', 'User time': 0.0742967079859227, 'Wallclock time': 0.0742967079859227, 'Termination condition': 'optimal', 'Iterations': 1, 'Timing': Bunch(integer cut generation = 0.00014595797983929515, main_timer_start_time = 338634.26841425, mip = 0.06020570796681568, nlp = 0.008367584028746933, total = 0.0742967079859227)}]}
# ── Pretty-print results ──────────────────────────────────────────────────────sep ="─"*50print(f"\n{'═'*50}")print(f" Solver Summary: {model.name}")print(f"{'═'*50}")# Solver metadatasolver_info = results.solverprint(f" Solver status : {solver_info.status}")print(f" Termination cond: {solver_info.termination_condition}")print(sep)# Objective & decision variablesprint(f" Objective value : {pyo.value(model.objective):.6f}")print(sep)print(f" Decision Variables")print(f" x = {pyo.value(model.x):.6f}")print(f" y = {pyo.value(model.y):.6f}")print(sep)# Active disjunctactive ="fix_x (x = 0)"if pyo.value(model.fix_x.binary_indicator_var) >0.5\else"fix_y (y = 0)"print(f" Active disjunct : {active}")print(f"{'═'*50}\n")
══════════════════════════════════════════════════
Solver Summary: LOA example
══════════════════════════════════════════════════
Solver status : ok
Termination cond: optimal
──────────────────────────────────────────────────
Objective value : 0.100000
──────────────────────────────────────────────────
Decision Variables
x = 0.000000
y = 1.000000
──────────────────────────────────────────────────
Active disjunct : fix_x (x = 0)
══════════════════════════════════════════════════
problem_info = results.problem# Timingprint(sep)solver_time =getattr(solver_info, 'user_time', None)wall_time =getattr(solver_info, 'wallclock_time', None)print(f" Timing")print(f" CPU time : {solver_time:.4f} s"if solver_time isnotNoneelse" CPU time : n/a")print(f" Wall time : {wall_time:.4f} s"if wall_time isnotNoneelse" Wall time : n/a")print(sep)# Problem sizen_vars =getattr(problem_info, 'number_of_variables', None)n_bin =getattr(problem_info, 'number_of_binary_variables', None)n_cont =getattr(problem_info, 'number_of_continuous_variables',None)n_cons =getattr(problem_info, 'number_of_constraints', None)print(f" Problem Size")print(f" Variables : {n_vars}"if n_vars isnotNoneelse" Variables : n/a")print(f" Binary : {n_bin}"if n_bin isnotNoneelse" -Binary : n/a")print(f" Continuous: {n_cont}"if n_cont isnotNoneelse" -Continuous: n/a")print(f" Constraints : {n_cons}"if n_cons isnotNoneelse" Constraints : n/a")print(sep)
──────────────────────────────────────────────────
Timing
CPU time : 0.0743 s
Wall time : 0.0743 s
──────────────────────────────────────────────────
Problem Size
Variables : 4
Binary : 2
Continuous: 2
Constraints : 3
──────────────────────────────────────────────────