discopt.solvers.convex_ipm_pounce

discopt.solvers.convex_ipm_pounce#

POUNCE’s dedicated convex LP/QP interior-point method (lp-ipm / qp-ipm).

The other POUNCE matrix backends (discopt.solvers.lp_pounce, discopt.solvers.qp_pounce) hand an LP or QP to POUNCE’s general NLP engine — the filter line-search interior-point method — through LP/QP-shaped callbacks. This module instead calls pounce.qp.solve_qp(), the specialised convex solver (Mehrotra predictor-corrector with presolve) that POUNCE’s own solver_selection=lp-ipm / qp-ipm reach. It backs Model.solve(solver="pounce") (#1533), whose purpose is that a pure LP or a convex QP is answered by that algorithm and by nothing else.

Both entry points follow the shared matrix contracts (discopt.solvers.lp_simplex.solve_lp() / qp_pounce.solve_qp): same signature, same LPResult / QPResult with HiGHS-convention duals, so solver._solve_lp_matrix and solver._solve_qp_matrix apply their primal-feasibility and KKT guards to the returned point unchanged.

Certificates#

POUNCE’s convex engine reports primal_infeasible and dual_infeasible (unbounded) from its own detection. Neither is turned into a certificate on that say-so alone, by the same rule the NLP-engine routes follow (#1309, #940):

  • primal_infeasible becomes INFEASIBLE only when the exact Rust simplex (lp_pounce._simplex_feasibility_verdict()) proves the linear system empty; otherwise ERROR.

  • dual_infeasible becomes UNBOUNDED only when an improving recession ray is verified exactly (lp_pounce._certify_unbounded_ray()) and the simplex exhibits a feasible point; otherwise ERROR.

Both checks only withhold a verdict; neither can produce an answer the IPM did not, so the route never silently becomes a simplex solve.

Bounds#

A declared bound at or beyond lp_pounce.finite_bound_threshold() is handed to the engine as infinite, exactly as on the other POUNCE routes. That includes the ±9.999e19 box a column declared without bounds receives, which, kept finite, costs the IPM several times the iterations. The caller passes relaxes_huge_bounds=True to the matrix route so an UNBOUNDED verdict over such a relaxed box is not certified (#850).

Objective scale#

The engine’s stopping test is not invariant to the scale of the objective: the same 2080-column production LP that stops at a primal residual of 5.5e-9 with its costs in dollars stops at 4.9e-6 with the costs in cents (x100) – 3.5x over the per-row test _matrix_solution_feasible applies – and the solve came back error (#1537). A QP 1e8*x + x**2/2 on [0, 10] ended in numerical_failure. _solve() therefore hands the engine sigma*P, sigma*c with sigma the power of two that brings the largest objective coefficient to at most 1 (objective_scale()), and maps the answer back: x is unchanged, objective and every dual divide by sigma. A power of two makes both directions exact, and sigma = 1 (no change at all) for every objective whose coefficients are already at most 1.

The engine’s tol is an absolute bound on its KKT residual (POUNCE QpOptions::tol, default 1e-8; a scale-relative arm opens only below the finite-precision floor). Dual infeasibility and complementarity are linear in the objective, so a solve of sigma*objective stopped at tol is stopped at tol/sigma in the caller’s units – 8192x looser on 3200*j**2, where the mapped-back complementarity of 6.6e-6 failed the #1384 stationarity guard and the solve came back error. tol*sigma (engine_tol()) is the engine tolerance that holds the caller’s tol in the caller’s units – down to a floor of 16*eps (_ENGINE_TOL_FLOOR), below which the scaled problem (coefficients at most 1) cannot resolve its residual and the engine runs to iteration_limit.

Handing it tol*sigma unconditionally cost iterations for nothing on LPs whose complementarity bottoms out at roundoff (#1658). The engine measures complementarity on the recomputed slack, so it cannot fall below about eps * |row activity| * |dual|; on a 10-flow blending LP with a contract row tight at every feasible point (no strict interior, sigma = 2**-7) that floor is ~5e-9 in engine units, tol*sigma asked for 7.8e-11, and the engine ran 59 iterations to its stall exit where 16 reached the same point (caller-unit complementarity 6.6e-7 vs 6.0e-7). The primal residual, which the issue suspected, was 1-2e-11 and passed. So _solve() runs the engine at the caller’s tol first and re-runs it at tol*sigma, warm-started from that point, only when the mapped-back residual fails caller_unit_converged(): primal absolute, dual infeasibility and complementarity against tol times the stationarity scale max(1, |c|, |Px|, |y|, |z|) (Ipopt’s s_d, the yardstick the #1384 guard applies). The 3200*j**2 class above fails that test and takes the re-run; the blending LP passes it and does not. The floor only ever loosens the request; the #1384 stationarity guard downstream still judges every point the engine returns.

Column scale#

qp-ipm’s iteration count depended on the units of the variables (#1683). The linear-MPC QP x[k+1] = 0.95 x[k] + (0.1/UA) u[k], u in [-2UA, 2UA], objective sum x**2 + 0.01 sum u**2 / UA**2 is the same problem for every UA, yet at N = 200 it took 13 iterations at UA = 1 and 115 at UA = 6e4 (94 vs 14 at N = 1200); pounce.qp.solve_qp has no scaling argument. For a QP, _solve() therefore hands the engine the substitution x = d * xt: D P D, d * c, A D, G D, lb / d, ub / d, with d from column_scale() (power-of-two Ruiz on the columns of [P; A; G]), and maps back x = d * xt, z_lb = zt_lb / d, z_ub = zt_ub / d; y, z and the objective are unchanged. Every d_j is a power of two, so both directions are exact. A warm start is mapped forward the same way. Measured after: 13 iterations at every UA in {1, 1e3, 6e4, 1e6} and at both N.

The scale is one-sided (d >= 1): only columns whose entries are small (a variable in large units) are scaled up. The engine’s dual residual maps back as r_j / (sigma * d_j), so d_j >= 1 never amplifies it and tol*sigma still holds the caller’s tol in caller units. A two-sided scale was measured and rejected: on QPLIB_8938 a d_j = 2**-11 column sent the caller-unit dual residual to 5.0e-8 (3.9e-12 unscaled), and the re-solve that would restore it must ask the engine for tol*sigma*min(d) = 4.9e-12 – below the 1.1e-11 row-residual floor, which the engine’s single tol also has to meet – and ran 18 -> 219 iterations to the limit.

Because the engine’s residuals are then in scaled coordinates, they are not used: caller_residuals() recomputes POUNCE’s own residual definitions on the caller’s data at the point the route returns, and caller_unit_converged(), the solve report and the downstream #1384 / #1596 guards (which recompute from x and the multipliers on the caller’s matrices anyway) all judge that. A QP whose columns are already unit-scaled gets d = None and is handed over unchanged. DISCOPT_POUNCE_QP_COLSCALE=0 (COLSCALE_ENV) restores the unscaled hand-off; it is an opt-out for A/B, not a pending graduation.

Graduation panel (CLAUDE.md §5, flag OFF vs ON, interleaved, every incumbent re-checked with check_feasibility at 1e-6). Plain convex QPLIB, 9 instances that route to qp-ipm: bit-identical in both arms (their columns are already balanced, d = None). The same QPLIB instances with variables re-expressed in random units (x' = 10**U(lo, hi) * x, 16 pairs that reach qp-ipm): total IPM iterations 1141 -> 751, 10 better / 2 same / 4 worse; QPLIB_8785 went from time_limit to certified optimal, and on QPLIB_8495 the OFF arm stopped at a feasible point of objective 210935 while the ON arm reached the true 42857. No instance lost certification, every certified bound stayed at or below its incumbent, and objectives that both arms returned agree to <= 1.5e-8 relative.

Attributes#

CONVEX_OPTION_KEYS

QP_ENGINE_OPTION_KEYS

COLSCALE_ENV

Exceptions#

IndefiniteQPError

The QP's Hessian was not PROVED positive semidefinite.

Functions#

certify_psd(→ bool)

True only when Q (symmetrized) is PROVED positive semidefinite (#1533 review).

convex_engine_options(→ dict)

Validate options for the convex engine and return a copy.

objective_scale(→ float)

The power of two sigma <= 1 that brings max(|P|, |c|) to at most 1.

engine_tol(→ Optional[float])

The engine tolerance that keeps tol in the caller's units under sigma.

caller_unit_converged(→ bool)

Whether the caller-unit result res meets tol (#1658).

colscale_enabled(→ bool)

Whether the qp-ipm hand-off is column-equilibrated (DISCOPT_POUNCE_QP_COLSCALE).

column_scale(→ Optional[numpy.ndarray])

Power-of-two Ruiz column equilibration of [P; A; G] (#1683).

caller_residuals(→ dict)

POUNCE's QpResiduals recomputed on the caller's (unscaled) data (#1683).

solve_lp(→ discopt.solvers.LPResult)

Solve min c'x s.t. A_ub x <= b_ub, A_eq x = b_eq with POUNCE's lp-ipm.

solve_qp(→ discopt.solvers.QPResult)

Solve min ½x'Qx + c'x s.t. linear rows with POUNCE's qp-ipm.

Module Contents#

discopt.solvers.convex_ipm_pounce.CONVEX_OPTION_KEYS: frozenset[str]#
discopt.solvers.convex_ipm_pounce.QP_ENGINE_OPTION_KEYS: frozenset[str]#
exception discopt.solvers.convex_ipm_pounce.IndefiniteQPError#

Bases: ValueError

The QP’s Hessian was not PROVED positive semidefinite.

Raised by solve_qp() when certify_psd() cannot prove Q PSD, and when POUNCE’s own PSD check refuses it. The solver="pounce" route catches it and solves the QP with the NLP engine instead, reporting a local result.

discopt.solvers.convex_ipm_pounce.certify_psd(Q: numpy.ndarray) → bool#

True only when Q (symmetrized) is PROVED positive semidefinite (#1533 review).

POUNCE’s own check accepts lam_min >= -1e-8 * max|Q_ij|, which admits an indefinite Hessian whose negative curvature is small against its largest entry but not against the box (-1e-9 x**2 + y**2 on x in [-1e4, 1e4]): the convex IPM then stops at a saddle and the route would certify it optimal. So convexity is decided here, before the QP arm can label anything optimal.

#1679: this is discopt._relax.convexity.eigenvalue.psd_certified(), the same predicate the exact QP objective certificate behind Model.convexity() and the default solver’s convex fast path use, so the two can no longer disagree on a Hessian (they did: m.convexity() accepted the float Gram matrix of a rank-deficient least-squares model – indefinite in exact arithmetic – that this proof refused). See that function for the procedure.

discopt.solvers.convex_ipm_pounce.convex_engine_options(options: dict | None) → dict#

Validate options for the convex engine and return a copy.

Raises:

ValueError – on a key outside CONVEX_OPTION_KEYS.

discopt.solvers.convex_ipm_pounce.objective_scale(P, c: numpy.ndarray) → float#

The power of two sigma <= 1 that brings max(|P|, |c|) to at most 1.

1.0 when every objective coefficient is already at most 1 in magnitude, or the objective is zero. See the module docstring (“Objective scale”).

discopt.solvers.convex_ipm_pounce.engine_tol(tol: float | None, sigma: float) → float | None#

The engine tolerance that keeps tol in the caller’s units under sigma.

None (the engine default) passes through when sigma == 1, so an unscaled solve is bit-for-bit what it was. See the module docstring.

discopt.solvers.convex_ipm_pounce.caller_unit_converged(res: Any, P, c: numpy.ndarray, tol: float) → bool#

Whether the caller-unit result res meets tol (#1658).

The primal residual is tested absolutely (sigma does not touch it); dual infeasibility and complementarity, which carry the objective’s units, against tol times _stationarity_scale_caller(). A result with no residual breakdown is not judged converged here, so it gets the tol*sigma re-solve.

discopt.solvers.convex_ipm_pounce.COLSCALE_ENV = 'DISCOPT_POUNCE_QP_COLSCALE'#
discopt.solvers.convex_ipm_pounce.colscale_enabled() → bool#

Whether the qp-ipm hand-off is column-equilibrated (DISCOPT_POUNCE_QP_COLSCALE).

discopt.solvers.convex_ipm_pounce.column_scale(P, A, G, n: int) → numpy.ndarray | None#

Power-of-two Ruiz column equilibration of [P; A; G] (#1683).

Returns d (every entry an exact power of two) such that the columns of D P D, A D and G D have infinity norm in [2**-0.5, 2**0.5) up to the coupling P introduces, or None when no column moves (d = 1), so a QP that is already unit-scaled is handed over unchanged. A column with no entry in any of the three matrices keeps d_j = 1.

discopt.solvers.convex_ipm_pounce.caller_residuals(x: numpy.ndarray, y, z, z_lb, z_ub, P, c: numpy.ndarray, A, b, G, h, lb: numpy.ndarray, ub: numpy.ndarray) → dict#

POUNCE’s QpResiduals recomputed on the caller’s (unscaled) data (#1683).

Same definitions as pounce-convex’s kkt_residuals: primal max(|Ax-b|, max(0, Gx-h), bound violation), dual ||Px + c + A'y + G'z - z_lb + z_ub||_inf and complementarity max |z_i * slack_i| over inequality rows and finite bounds. A column-equilibrated engine reports these in its own coordinates, where an infinity norm of a per-column-scaled residual has no exact inverse; the caller-unit test (caller_unit_converged()) and the solve report need them in the caller’s.

discopt.solvers.convex_ipm_pounce.solve_lp(c: numpy.ndarray, A_ub: numpy.ndarray | scipy.sparse.spmatrix | None = None, b_ub: numpy.ndarray | None = None, A_eq: numpy.ndarray | scipy.sparse.spmatrix | None = None, b_eq: numpy.ndarray | None = None, bounds: List[Tuple[float, float]] | None = None, warm_basis: object | None = None, time_limit: float | None = None, options: dict | None = None, solve_report: bool = False) → discopt.solvers.LPResult#

Solve min c'x s.t. A_ub x <= b_ub, A_eq x = b_eq with POUNCE’s lp-ipm.

bounds default to (0, +inf) (the shared LP contract). warm_basis is accepted for signature compatibility and ignored: an IPM has no basis. options must be a subset of CONVEX_OPTION_KEYS. solve_report=True attaches the solve’s pounce.solve-report/v1 document as LPResult.solve_report (#1534), built by discopt.solvers._pounce_report.convex_report().

discopt.solvers.convex_ipm_pounce.solve_qp(Q: numpy.ndarray, c: numpy.ndarray, A_ub: numpy.ndarray | scipy.sparse.spmatrix | None = None, b_ub: numpy.ndarray | None = None, A_eq: numpy.ndarray | scipy.sparse.spmatrix | None = None, b_eq: numpy.ndarray | None = None, bounds: List[Tuple[float, float]] | None = None, integrality: numpy.ndarray | None = None, time_limit: float | None = None, gap_tolerance: float = 0.0001, options: dict | None = None, solve_report: bool = False, warm_start: numpy.ndarray | collections.abc.Mapping | None = None) → discopt.solvers.QPResult#

Solve min ½x'Qx + c'x s.t. linear rows with POUNCE’s qp-ipm.

bounds default to free variables (the shared QP contract). solve_report=True attaches the solve’s pounce.solve-report/v1 document as QPResult.solve_report (#1534). warm_start seeds the interior-point iteration and never changes the answer: a primal point of length n (#1615 B-01b), or a mapping with x and any of the multiplier blocks y (A_eq rows), z (A_ub rows), z_lb, z_ub in the caller’s objective units – e.g. a previous result’s QPResult.warm_start_state (#1679).

Raises:

  • IndefiniteQPError – Q is not PSD (POUNCE’s own check).

  • ValueError – any integer-marked variable, or an unknown option.