NLP Backend Comparison: POUNCE vs Ipopt#
discopt solves continuous relaxations through two distinct NLP engines:
Backend |
Flag |
Description |
|---|---|---|
POUNCE |
|
Pure-Rust primal-dual interior-point method [Nocedal and Wright, 2006] via PyO3 bindings; the universal default |
Ipopt |
|
Industry-standard NLP solver [Wächter and Biegler, 2006] via cyipopt |
The comparison below runs a third selector, nlp_solver="ipm". This is a
back-compat alias, not a third engine: the pure-JAX interior-point method it
used to name has been retired, and on the NLP/MINLP path "ipm" now resolves to
POUNCE. It is kept in the sweep deliberately — the pounce and IPM columns
below agree to every printed digit, and to the node on the MINLPs, which is what
that alias resolving to the same engine looks like from the outside.
This notebook measures:
Accuracy: objective value agreement with known optima
Speed: wall-clock solve time
Robustness: convergence status
import time
import discopt.modeling as dm
import numpy as np
print("discopt loaded successfully")
discopt loaded successfully
Test Problems#
We define a mix of problem types built with dm.Model() and solve with each backend selector.
Both backends support dm.Model()-built models and .nl file models alike. Derivatives
come from the same source in either case: a Rust automatic-differentiation tape built by
discopt._nl_expr_compiler and consumed through TapeNLPEvaluator, so the objective,
gradient, Jacobian, and Lagrangian Hessian are exact rather than finite-differenced.
def make_unconstrained_quadratic():
"""Simple unconstrained quadratic. Optimal: f*=0 at (3, -1)."""
m = dm.Model("quadratic")
x = m.continuous("x", lb=-10, ub=10)
y = m.continuous("y", lb=-10, ub=10)
m.minimize((x - 3) ** 2 + (y + 1) ** 2)
return m, 0.0
def make_rosenbrock():
"""Rosenbrock function. Optimal: f*=0 at (1, 1)."""
m = dm.Model("rosenbrock")
x = m.continuous("x", lb=-5, ub=5)
y = m.continuous("y", lb=-5, ub=5)
m.minimize((1 - x) ** 2 + 100 * (y - x**2) ** 2)
return m, 0.0
def make_constrained_quadratic():
"""Constrained quadratic. Optimal: f*=0.5 at (0.5, 0.5)."""
m = dm.Model("constrained_quad")
x = m.continuous("x", lb=-5, ub=5)
y = m.continuous("y", lb=-5, ub=5)
m.minimize(x**2 + y**2)
m.subject_to(x + y >= 1)
return m, 0.5
def make_nonlinear_equality():
"""Nonlinear equality constraint. Optimal: f*=2 at (1, 1)."""
m = dm.Model("nl_equality")
x = m.continuous("x", lb=0.1, ub=5)
y = m.continuous("y", lb=0.1, ub=5)
m.minimize(x**2 + y**2)
m.subject_to(x * y == 1)
return m, 2.0
def make_exp_nlp():
"""Exponential objective with constraint. Optimal: f*~1.8395."""
from scipy.optimize import fsolve
# Solve KKT: exp(x)=2y, x+y=1 => exp(x)=2(1-x)
x_star = fsolve(lambda x: np.exp(x) - 2 * (1 - x), 0.0)[0]
y_star = 1 - x_star
opt = np.exp(x_star) + y_star**2
m = dm.Model("exp_nlp")
x = m.continuous("x", lb=-2, ub=2)
y = m.continuous("y", lb=-2, ub=2)
m.minimize(dm.exp(x) + y**2)
m.subject_to(x + y >= 1)
return m, float(opt)
def make_log_nlp():
"""Logarithmic objective. Optimal: f*=2*ln(2)~1.3863 at (0.5, 0.5)."""
import math
m = dm.Model("log_nlp")
x = m.continuous("x", lb=0.01, ub=0.99)
y = m.continuous("y", lb=0.01, ub=0.99)
m.minimize(-dm.log(x) - dm.log(y))
m.subject_to(x + y <= 1)
return m, 2 * math.log(2)
def make_multi_constraint():
"""Multiple constraints. Optimal: f*=0.5 at (1.5, 1.5)."""
m = dm.Model("multi_constraint")
x = m.continuous("x", lb=0, ub=5)
y = m.continuous("y", lb=0, ub=5)
m.minimize((x - 2) ** 2 + (y - 2) ** 2)
m.subject_to(x + y <= 3)
return m, 0.5
def make_bilinear():
"""Bilinear constraints. Optimal: f*=2 at (1, 1)."""
m = dm.Model("bilinear")
x = m.continuous("x", lb=0.1, ub=5)
y = m.continuous("y", lb=0.1, ub=5)
m.minimize(x**2 + y**2)
m.subject_to(x * y >= 1)
m.subject_to(x + y <= 4)
return m, 2.0
def make_trig_nlp():
"""Trigonometric objective. Optimal: f*~-1.0 at x~pi/2."""
m = dm.Model("trig_nlp")
x = m.continuous("x", lb=0.1, ub=3.0)
y = m.continuous("y", lb=0.1, ub=3.0)
m.minimize(-dm.sin(x) * dm.cos(y))
m.subject_to(x + y <= 4)
return m, -1.0 # approximate
problems = [
("Unconstrained quadratic", make_unconstrained_quadratic),
("Rosenbrock", make_rosenbrock),
("Constrained quadratic", make_constrained_quadratic),
("Nonlinear equality", make_nonlinear_equality),
("Exponential NLP", make_exp_nlp),
("Logarithmic NLP", make_log_nlp),
("Multiple constraints", make_multi_constraint),
("Bilinear constraints", make_bilinear),
("Trigonometric NLP", make_trig_nlp),
]
print(f"Defined {len(problems)} test problems")
Defined 9 test problems
Run Comparison#
For each problem, solve with all three backends and record objective, time, and status.
backends = ["pounce", "ipm", "ipopt"]
results = []
for name, make_fn in problems:
model, expected = make_fn()
row = {"problem": name, "expected": expected}
for backend in backends:
m_copy, _ = make_fn() # fresh model each time
try:
t0 = time.perf_counter()
result = m_copy.solve(nlp_solver=backend, time_limit=30)
elapsed = time.perf_counter() - t0
row[f"{backend}_obj"] = result.objective
row[f"{backend}_time"] = elapsed
row[f"{backend}_status"] = result.status
if result.objective is not None:
row[f"{backend}_error"] = abs(result.objective - expected)
else:
row[f"{backend}_error"] = None
except Exception as e:
row[f"{backend}_obj"] = None
row[f"{backend}_time"] = None
row[f"{backend}_status"] = f"ERROR: {e}"
row[f"{backend}_error"] = None
results.append(row)
print(f"Solved {len(results)} problems with {len(backends)} backends")
******************************************************************************
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
******************************************************************************
Solved 9 problems with 3 backends
Results: Accuracy Comparison#
header = (
f"{'Problem':<26s} {'Expected':>10s} "
f"{'pounce obj':>12s} {'pounce err':>10s} "
f"{'IPM obj':>12s} {'IPM err':>10s} "
f"{'Ipopt obj':>12s} {'Ipopt err':>10s}"
)
print(header)
print("-" * len(header))
for r in results:
cols = []
cols.append(f"{r['problem']:<26s} {r['expected']:>10.6f}")
for b in backends:
obj = r.get(f"{b}_obj")
err = r.get(f"{b}_error")
obj_s = f"{obj:.6f}" if obj is not None else "FAIL"
err_s = f"{err:.2e}" if err is not None else "N/A"
cols.append(f"{obj_s:>12s} {err_s:>10s}")
print(" ".join(cols))
Problem Expected pounce obj pounce err IPM obj IPM err Ipopt obj Ipopt err
-------------------------------------------------------------------------------------------------------------
Unconstrained quadratic 0.000000 0.000000 3.55e-15 0.000000 3.55e-15 0.000000 3.55e-15
Rosenbrock 0.000000 -0.000000 1.83e-11 -0.000000 1.83e-11 -0.000000 1.83e-11
Constrained quadratic 0.500000 0.500000 2.51e-09 0.500000 2.51e-09 0.500000 2.51e-09
Nonlinear equality 2.000000 2.000000 0.00e+00 2.000000 0.00e+00 2.000000 0.00e+00
Exponential NLP 1.839484 1.839484 2.51e-09 1.839484 2.51e-09 1.839484 2.51e-09
Logarithmic NLP 1.386294 1.386294 2.51e-09 1.386294 2.51e-09 1.386294 2.51e-09
Multiple constraints 0.500000 0.500000 2.51e-09 0.500000 2.51e-09 0.500000 2.51e-09
Bilinear constraints 2.000000 2.000000 2.51e-09 2.000000 2.51e-09 2.000000 2.51e-09
Trigonometric NLP -1.000000 -0.995004 5.00e-03 -0.995004 5.00e-03 -0.995004 5.00e-03
Results: Speed Comparison#
header = f"{'Problem':<26s} {'pounce':>10s} {'IPM':>10s} {'Ipopt':>10s} {'Fastest':>8s}"
print(header)
print("-" * len(header))
totals = {b: 0.0 for b in backends}
for r in results:
times = {}
for b in backends:
t = r.get(f"{b}_time")
times[b] = t
if t is not None:
totals[b] += t
time_strs = []
for b in backends:
t = times[b]
time_strs.append(f"{t:.4f}s" if t is not None else "FAIL")
valid = {b: t for b, t in times.items() if t is not None}
fastest = min(valid, key=valid.get) if valid else "N/A"
print(f"{r['problem']:<26s} " + " ".join(f"{s:>10s}" for s in time_strs) + f" {fastest:>8s}")
print()
for b in backends:
print(f" {b} total: {totals[b]:.3f}s")
Problem pounce IPM Ipopt Fastest
--------------------------------------------------------------------
Unconstrained quadratic 0.1635s 0.0029s 0.0027s ipopt
Rosenbrock 15.3894s 14.6687s 15.7998s ipm
Constrained quadratic 0.0039s 0.0038s 0.0029s ipopt
Nonlinear equality 0.1750s 0.0056s 0.0044s ipopt
Exponential NLP 0.0035s 0.0038s 0.0035s ipopt
Logarithmic NLP 0.0026s 0.0027s 0.0027s pounce
Multiple constraints 0.0034s 0.0033s 0.0029s ipopt
Bilinear constraints 0.0046s 0.0051s 0.0044s ipopt
Trigonometric NLP 0.0168s 0.0186s 0.0184s pounce
pounce total: 15.763s
ipm total: 14.715s
ipopt total: 15.842s
Results: Status Summary#
header = f"{'Problem':<26s} {'pounce':<18s} {'IPM':<18s} {'Ipopt':<18s}"
print(header)
print("-" * len(header))
for r in results:
cols = [f"{r['problem']:<26s}"]
for b in backends:
s = r.get(f"{b}_status", "N/A")
cols.append(f"{s:<18s}")
print(" ".join(cols))
Problem pounce IPM Ipopt
-----------------------------------------------------------------------------------
Unconstrained quadratic optimal optimal optimal
Rosenbrock optimal optimal optimal
Constrained quadratic optimal optimal optimal
Nonlinear equality optimal optimal optimal
Exponential NLP optimal optimal optimal
Logarithmic NLP optimal optimal optimal
Multiple constraints optimal optimal optimal
Bilinear constraints optimal optimal optimal
Trigonometric NLP optimal optimal optimal
MINLP Comparison (Branch & Bound)#
For MINLP problems, the NLP backend is used at each Branch & Bound [Land and Doig, 1960] node. Let’s compare on a few small MINLPs.
def make_simple_minlp():
"""Simple MINLP. Optimal: f*=0.5."""
m = dm.Model("simple_minlp")
x1 = m.continuous("x1", lb=0, ub=5)
x2 = m.continuous("x2", lb=0, ub=5)
x3 = m.binary("x3")
m.minimize(x1**2 + x2**2 + x3)
m.subject_to(x1 + x2 >= 1)
m.subject_to(x1**2 + x2 <= 3)
return m, 0.5
def make_binary_knapsack():
"""Binary knapsack-like MINLP. Optimal: f*=5.0."""
m = dm.Model("binary_knapsack")
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(x <= 5 * y1)
m.subject_to(x <= 4 * y2)
return m, 5.0
def make_integer_quadratic():
"""Integer quadratic. Optimal: f*=0.0."""
m = dm.Model("integer_quad")
x = m.continuous("x", lb=0, ub=5)
n = m.integer("n", lb=0, ub=5)
m.minimize((x - 2.7) ** 2 + (n - 3) ** 2)
m.subject_to(x + n <= 6)
return m, 0.0
def make_exp_minlp():
"""Exponential MINLP. Binary switches continuous region."""
m = dm.Model("exp_minlp")
x = m.continuous("x", lb=0.1, ub=5)
y = m.continuous("y", lb=0.1, ub=5)
z = m.binary("z")
m.minimize(dm.exp(x) * y + z)
m.subject_to(x + y >= 2)
m.subject_to(x * y <= 5 * z + 1)
return m, None # optimal unknown a priori
minlp_problems = [
("Simple MINLP", make_simple_minlp),
("Binary knapsack", make_binary_knapsack),
("Integer quadratic", make_integer_quadratic),
("Exp MINLP", make_exp_minlp),
]
minlp_results = []
for name, make_fn in minlp_problems:
row = {"problem": name}
_, expected = make_fn()
row["expected"] = expected
for backend in backends:
m_copy, _ = make_fn()
try:
t0 = time.perf_counter()
result = m_copy.solve(nlp_solver=backend, time_limit=30, max_nodes=500)
elapsed = time.perf_counter() - t0
row[f"{backend}_obj"] = result.objective
row[f"{backend}_time"] = elapsed
row[f"{backend}_status"] = result.status
row[f"{backend}_nodes"] = result.node_count
except Exception as e:
row[f"{backend}_obj"] = None
row[f"{backend}_time"] = None
row[f"{backend}_status"] = f"ERROR: {e}"
row[f"{backend}_nodes"] = None
minlp_results.append(row)
# Print MINLP results
header = f"{'Problem':<20s} " + " ".join(
f"{b + ' obj':>12s} {b + ' time':>10s} {b + ' nodes':>8s}" for b in backends
)
print(header)
print("-" * len(header))
for r in minlp_results:
cols = [f"{r['problem']:<20s}"]
for b in backends:
obj = r.get(f"{b}_obj")
t = r.get(f"{b}_time")
n = r.get(f"{b}_nodes")
obj_s = f"{obj:.4f}" if obj is not None else "FAIL"
t_s = f"{t:.3f}s" if t is not None else "N/A"
n_s = str(n) if n is not None else "N/A"
cols.append(f"{obj_s:>12s} {t_s:>10s} {n_s:>8s}")
print(" ".join(cols))
OA: generating OA cuts for 1 of 2 rows (2 convex; the rest are either nonconvex or already exact in the master)
OA: generating OA cuts for 1 of 2 rows (2 convex; the rest are either nonconvex or already exact in the master)
OA: generating OA cuts for 1 of 2 rows (2 convex; the rest are either nonconvex or already exact in the master)
OA: generating OA cuts for 0 of 2 rows (2 convex; the rest are either nonconvex or already exact in the master)
OA: generating OA cuts for 0 of 2 rows (2 convex; the rest are either nonconvex or already exact in the master)
OA: generating OA cuts for 0 of 2 rows (2 convex; the rest are either nonconvex or already exact in the master)
OA: generating OA cuts for 0 of 1 rows (1 convex; the rest are either nonconvex or already exact in the master)
OA: generating OA cuts for 0 of 1 rows (1 convex; the rest are either nonconvex or already exact in the master)
OA: generating OA cuts for 0 of 1 rows (1 convex; the rest are either nonconvex or already exact in the master)
Problem pounce obj pounce time pounce nodes ipm obj ipm time ipm nodes ipopt obj ipopt time ipopt nodes
--------------------------------------------------------------------------------------------------------------------------------
Simple MINLP 0.5000 0.016s 0 0.5000 0.008s 0 0.5000 0.009s 0
Binary knapsack 5.0000 0.014s 0 5.0000 0.015s 0 5.0000 0.015s 0
Integer quadratic 0.0000 0.007s 0 0.0000 0.006s 0 0.0000 0.007s 0
Exp MINLP 0.6686 0.062s 5 0.6686 0.061s 5 0.6686 0.089s 5
Summary#
Key takeaways:
POUNCE (
nlp_solver="pounce") is the project’s own Rust interior-point solver, integrated via PyO3 bindings. It provides the Hessian of the full Lagrangian for better convergence with constraints, and it is the default for every problem class.Ipopt (
nlp_solver="ipopt") is the industry-standard NLP solver, and the more mature and battle-tested of the two. POUNCE is a Rust port of the same primal-dual interior-point algorithm, so the two agree closely on well-conditioned problems — as they do on every row above.nlp_solver="ipm"is an alias, not a backend. The pure-JAX IPM it named has been retired; the selector now resolves to POUNCE on the NLP/MINLP path (and to simplex-first matrix routing for LP/MILP). The identicalpounceandIPMcolumns above are the evidence, not a coincidence. Existing scripts passing"ipm"keep working; new code should say"pounce".Read the timing table with care. Each problem is a single untimed-warmup run, and the
pouncecolumn is measured first, so it absorbs one-off setup cost — that is the whole of the 0.1402s vs 0.0023s gap on the unconstrained quadratic, and of POUNCE’s larger total. It is not a POUNCE-vs-Ipopt result. For scaling measurements with median-of-3 timing, seebenchmarks_by_class.ipynb.POUNCE avoids the cyipopt/Ipopt system dependency and provides a fully self-contained solver stack; it ships as a core dependency of discopt.