Outer Approximation Solver

Outer Approximation Solver#

This tutorial introduces discopt’s general-purpose Outer Approximation (OA) solver for Mixed-Integer Nonlinear Programming (MINLP). OA decomposes a MINLP into alternating NLP subproblems and MILP master problems, and is especially effective for convex MINLP problems with many binary/integer variables [Duran and Grossmann, 1986].

When to use OA vs. Branch & Bound#

Criterion	OA	Spatial B&B
Convex MINLP	Preferred — finite convergence to global optimum	Valid but often slower
Non-convex MINLP	Local optimum only — cuts may be invalid	Global optimum with relaxation guarantees
Many binary variables	Strong — MILP solver handles combinatorics	Can struggle with large B&B tree
Few integers, hard NLP	Overkill — NLP dominates	Preferred

Algorithm Overview#

The OA algorithm [Duran and Grossmann, 1986, Fletcher and Leyffer, 1994] works as follows:

Initialize: Solve continuous NLP relaxation → first linearization point
Loop:
- Solve MILP master (linear constraints + accumulated OA cuts) → lower bound
- Fix integers at master solution, solve NLP subproblem → upper bound
- Generate OA cuts (tangent hyperplanes) at NLP solution
- If NLP infeasible: add feasibility cuts [Fletcher and Leyffer, 1994] or no-good cuts
Terminate when gap converges or master becomes infeasible

import os

os.environ["JAX_PLATFORMS"] = "cpu"
os.environ["JAX_ENABLE_X64"] = "1"

import discopt.modeling as dm

Example 1: Simple Convex MINLP#

\[\min_{x,y} \quad x_1^2 + x_2^2 + x_3\]

\[\text{s.t.} \quad x_1 + x_2 \geq 1, \quad x_1^2 + x_2 \leq 3, \quad x_3 \in \{0,1\}\]

The optimal solution is \(x_1 = x_2 = 0.5\), \(x_3 = 0\), with objective 0.5.

from discopt.modeling.examples import example_simple_minlp

m = example_simple_minlp()
result = m.solve(gdp_method="oa", time_limit=60)

print(f"Status: {result.status}")
print(f"Objective: {result.objective:.4f}")
print(f"Solution: {result.x}")
print(f"Wall time: {result.wall_time:.2f}s")

Model: textbook
  Variables: 3 (2 continuous, 1 integer/binary)
  Constraints: 2
  Objective: minimize (((x1 ** 2) + (x2 ** 2)) + x3)
  Parameters: 0

Status: optimal
Objective: 0.5000
Solution: {'x1': array(0.5), 'x2': array(0.5), 'x3': array(-1.e-12)}
Wall time: 1.16s

Example 2: Quadratic Objective with Binary Activation#

A common pattern in process design: a quadratic cost with binary decisions that activate or deactivate options.

\[\min \quad (x - 3)^2 + 2y_1 + 3y_2\]

\[\text{s.t.} \quad y_1 + y_2 \leq 1, \quad x \leq 2y_1 + 4y_2, \quad x \in [0, 5],\; y_i \in \{0,1\}\]

m = dm.Model("activation")
x = m.continuous("x", lb=0, ub=5)
y1 = m.binary("y1")
y2 = m.binary("y2")

m.minimize((x - 3) ** 2 + 2 * y1 + 3 * y2)
m.subject_to(y1 + y2 <= 1)
m.subject_to(x <= 2 * y1 + 4 * y2)

result = m.solve(gdp_method="oa", time_limit=60)
print(f"Status: {result.status}, Objective: {result.objective:.4f}")
print(f"x={result.x['x']:.3f}, y1={result.x['y1']:.0f}, y2={result.x['y2']:.0f}")

Status: feasible, Objective: 3.0000
x=2.000, y1=1, y2=-0

Comparing OA with Branch & Bound#

Both solvers find the same optimum for convex MINLP, but OA is typically faster when the integer structure dominates.

m_oa = example_simple_minlp()
m_bb = example_simple_minlp()

result_oa = m_oa.solve(gdp_method="oa", time_limit=60)
result_bb = m_bb.solve(time_limit=60)

oa_obj = result_oa.objective
oa_t = result_oa.wall_time
bb_obj = result_bb.objective
bb_t = result_bb.wall_time
print(f"OA:  obj={oa_obj:.4f}, time={oa_t:.2f}s, status={result_oa.status}")
print(f"B&B: obj={bb_obj:.4f}, time={bb_t:.2f}s, status={result_bb.status}")

Model: textbook
  Variables: 3 (2 continuous, 1 integer/binary)
  Constraints: 2
  Objective: minimize (((x1 ** 2) + (x2 ** 2)) + x3)
  Parameters: 0
Model: textbook
  Variables: 3 (2 continuous, 1 integer/binary)
  Constraints: 2
  Objective: minimize (((x1 ** 2) + (x2 ** 2)) + x3)
  Parameters: 0

******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
 Ipopt is released as open source code under the Eclipse Public License (EPL).
         For more information visit https://github.com/coin-or/Ipopt
******************************************************************************

OA:  obj=0.5000, time=0.83s, status=optimal
B&B: obj=0.5000, time=1.14s, status=optimal

Extended Cutting Plane (ECP) Mode#

The Extended Cutting Plane method [Westerlund and Pettersson, 1995] avoids NLP subproblem solves entirely. Instead, it generates OA cuts at the MILP master solution for violated nonlinear constraints. This is simpler but converges more slowly.

m = example_simple_minlp()
result_ecp = m.solve(gdp_method="oa", ecp_mode=True, time_limit=60)

print(f"ECP Status: {result_ecp.status}")
print(f"ECP Objective: {result_ecp.objective}")

Model: textbook
  Variables: 3 (2 continuous, 1 integer/binary)
  Constraints: 2
  Objective: minimize (((x1 ** 2) + (x2 ** 2)) + x3)
  Parameters: 0

ECP Status: optimal
ECP Objective: 0.5000000863278043

Equality Relaxation#

For problems with nonlinear equality constraints, the OA linearizations can make the MILP master infeasible. The Equality Relaxation strategy [Viswanathan and Grossmann, 1990] relaxes these to inequalities in the OA cuts, restoring master feasibility.

m = dm.Model("eq_relax")
x = m.continuous("x", lb=0, ub=2)
y = m.binary("y")
m.minimize(x**2 + y)
m.subject_to(x**2 - y == 0)  # nonlinear equality

result_er = m.solve(gdp_method="oa", equality_relaxation=True, time_limit=60)
print(f"Status: {result_er.status}, Objective: {result_er.objective}")

OA: generating OA cuts only for 0 of 1 constraints classified convex

Status: feasible, Objective: 3.333740079851804e-12

Configuration Options#

Parameter	Default	Description
`gdp_method="oa"`	—	Select the OA solver
`time_limit`	3600	Wall-clock time limit (seconds)
`gap_tolerance`	1e-4	Relative optimality gap
`max_nodes`	100000	Maximum OA iterations
`ecp_mode`	False	Extended Cutting Plane mode
`equality_relaxation`	False	Relax nonlinear equalities
`feasibility_cuts`	True	Gradient-based feasibility cuts
`nlp_solver`	“ipm”	NLP backend: “ipm”, “ipopt”, “ripopt”

Limitations#

Non-convex problems: OA only guarantees local optimality. For global optimality on non-convex MINLP, use the default spatial B&B solver.
Single-tree: The current implementation uses a multi-tree approach (re-solves the MILP from scratch each iteration). Single-tree LP/NLP Branch-and-Bound [Quesada and Grossmann, 1992] is a future enhancement.
Regularization: Level-method regularization [Kronqvist et al., 2020] to stabilize integer selections is not yet implemented.

For a comprehensive review and comparison of convex MINLP solvers, see [Kronqvist et al., 2019]. The Bonmin solver [Bonami et al., 2008] provides multiple algorithmic variants including OA, NLP-BB, and hybrid approaches.