Using discopt as a Pyomo solver

discopt ships an optional Pyomo [Hart et al., 2017] plugin so that models written with Pyomo’s algebraic modeling layer can be solved with discopt through the familiar SolverFactory interface — no rewrite required.

Install the extra:

pip install discopt[pyomo]

Under the hood the bridge round-trips the model through a temporary AMPL .nl file [Fourer et al., 2003]: Pyomo’s .nl writer and discopt’s from_nl agree on the variable column order, so the solution maps back by index. This notebook shows a basic solve, solver options, and dual values.

Activate the plugin

Importing discopt.pyomo — or calling discopt.pyomo.register() explicitly — registers the 'discopt' solver with Pyomo’s SolverFactory. A pyomo.solvers entry point is also declared, so discovery can happen automatically after installation.

import pyomo.environ as pyo
import discopt.pyomo  # registers the 'discopt' solver

print("discopt registered:", "discopt" in pyo.SolverFactory)
opt = pyo.SolverFactory("discopt")
print("available:", opt.available())

discopt registered: True
available: True

A small MINLP

A mixed-integer nonlinear program with one continuous variable, one binary variable, a quadratic objective, and a linear constraint:

\[\min_{x,\,y}\;(x-3)^2 + 2y \quad\text{s.t.}\quad x + 5y \ge 4,\; x \in [0,10],\; y \in \{0,1\}.\]

The optimum is \(x=4,\,y=0\) with objective \(1\).

m = pyo.ConcreteModel()
m.x = pyo.Var(bounds=(0, 10))
m.y = pyo.Var(domain=pyo.Binary)
m.obj = pyo.Objective(expr=(m.x - 3) ** 2 + 2 * m.y, sense=pyo.minimize)
m.c = pyo.Constraint(expr=m.x + 5 * m.y >= 4)

results = opt.solve(m)
print("termination:", results.solver.termination_condition)
print("objective:  ", pyo.value(m.obj))
print("x =", m.x.value, "  y =", m.y.value)

OA: generating OA cuts for 0 of 1 rows (1 convex; the rest are either nonconvex or already exact in the master)

termination: optimal
objective:   1.000000002507477
x = 4.0000000012537384   y = 0.0

discopt solved the model and loaded the solution back onto the Pyomo variables. The binary variable comes back as an exact integer (integrality drift is rounded away on load).

Solver options

Pyomo solver options map onto Model.solve keyword arguments:

Pyomo option	discopt solve() kwarg
timelimit (or time_limit)	time_limit
mipgap / gap	gap_tolerance
threads	threads

Unrecognised options are forwarded verbatim to Model.solve, so any discopt solve keyword works without a plugin change. Options can be given per solve (options={...}) or set persistently on the solver object (opt.options[...]).

tee=True is not one of the forwarded options. For the duration of the solve the plugin attaches a stdout handler to the discopt logger (lifting it to INFO if it is quieter), so the solver log streams to the terminal the way Pyomo users expect, and then restores the logger. It is deliberately not mapped onto Model.solve(stream=True): that keyword asks for an iterator of SolveUpdate objects instead of a SolveResult and is not implemented, which used to make every tee=True solve return termination error with no solution loaded (issue #1178).

results = opt.solve(m, options={"timelimit": 30, "mipgap": 1e-4})
print("termination:", results.solver.termination_condition)
print("wall time (s):", results.solver.wallclock_time)

# The same options set persistently, applying to every later solve on this object:
opt.options["timelimit"] = 60
print("persistent options:", dict(opt.options))

OA: generating OA cuts for 0 of 1 rows (1 convex; the rest are either nonconvex or already exact in the master)

termination: optimal
wall time (s): 0.00956199993379414

Dual values and reduced costs

Declare an import Suffix and discopt’s KKT multipliers are loaded best-effort. Consider a convex NLP whose optimum sits on an active constraint:

\[\min_x\;(x-3)^2 \quad\text{s.t.}\quad x \ge 4.\]

The optimum is \(x=4\), and the constraint multiplier is \(2\) (the AMPL/Pyomo sign convention, matching e.g. Ipopt).

n = pyo.ConcreteModel()
n.x = pyo.Var(bounds=(0, 10))
n.obj = pyo.Objective(expr=(n.x - 3) ** 2, sense=pyo.minimize)
n.c = pyo.Constraint(expr=n.x >= 4)
n.dual = pyo.Suffix(direction=pyo.Suffix.IMPORT)
n.rc = pyo.Suffix(direction=pyo.Suffix.IMPORT)

opt.solve(n)
print("x =", round(n.x.value, 6))
print("dual[c] =", n.dual.get(n.c))
print("rc[x]   =", n.rc.get(n.x))

x = 4.0
dual[c] = 2.0000000022985165
rc[x]   = 2.0898051791981e-10

Duals are populated only when discopt exposes multipliers (the NLP/KKT path); for problems solved by another route the Suffix is simply left empty — values are never fabricated.

How it works and limitations

Pyomo’s .nl writer emits the model in its original variable space (presolve and scaling disabled), discopt’s from_nl reads it with the same column/row order, and the solution maps back by index. As a result:

• Variable and constraint values map by column/row order, not by name (Pyomo and discopt name objects differently).

• Constraint names are not preserved through .nl (cosmetic; solving and duals are unaffected).

• A model whose operators discopt’s reader does not support returns a structured error result rather than crashing.

Importing a Pyomo model as a discopt Model

To get a native discopt Model out of a Pyomo ConcreteModel — rather than just solving it in place — use discopt.modeling.from_pyomo. That is the entry point when you want discopt’s own API afterwards (relaxations, transformations, serialization) instead of Pyomo’s result object.

import discopt.modeling as dm

p = pyo.ConcreteModel()
p.x = pyo.Var(bounds=(0, 10))
p.obj = pyo.Objective(expr=(p.x - 3) ** 2)

dmodel = dm.from_pyomo(p)  # a native discopt Model
print(type(dmodel).__name__, "with", len(dmodel._variables), "variable(s)")

res = dmodel.solve()
print("status:", res.status, " objective:", round(res.objective, 6))

from_pyomo reuses the same temporary .nl round-trip as the solver plugin, so variables and constraints come back in .nl column/row order and their names differ from the Pyomo model’s. A future direct in-memory translator may replace the round-trip without changing this interface.