diff --git a/docs/dev/certification-gap-plan.md b/docs/dev/certification-gap-plan.md index 8a6ebe45..37c03b20 100644 --- a/docs/dev/certification-gap-plan.md +++ b/docs/dev/certification-gap-plan.md @@ -1705,6 +1705,121 @@ convexity across zero. 3. The per-family 200-box property test stands, drawing boxes from the family's admissible sign regime. +*Re-measured and NARROWED (2026-07-26, issue #861 — the measurement wins, §0.4).* +Re-scope item 2 above is **too broad for the even powers**. The 2026-07-02 probe was +run against the *unflagged* `build_milp_relaxation`; the incremental structure +compares itself against the build with `skip_separable_floor=True, +skip_convex_lift=True`, and the sign-dependence the probe saw for `x²`/`x⁴` was +entirely the extra separable-floor row `s ≥ 0` — a row the incremental path never +regenerates and never sees. Under the flags the structure actually uses, the row +count is (re-probe: `p ∈ {2,3,4}` × boxes `[-2,3]`, `[1,3]`, `[-3,-1]`): + +| power | spanning `[-2,3]` | positive `[1,3]` | negative `[-3,-1]` | +|---|---|---|---| +| `x²` | **4** | 4 | 4 | +| `x⁴` | **4** | 4 | 4 | +| `x³` | **2** | 4 | 4 | + +So the sign-regime hazard is **odd powers only**: `f'' = p(p-1)x^(p-2) ≥ 0` for even +`p` on all of ℝ, hence one convexity, hence the same secant + 3-tangent envelope in +every regime; an odd power is S-shaped across zero and the cold build swaps to +`_emit_odd_power_hull`'s 2 facets — a facet-*count* change, which is what a fixed +sparsity pattern genuinely cannot express. The gate is therefore +`root_sign == 0 and p odd`, not `root_sign == 0`. + +Two generalizations were prerequisites, both now pinned by +`python/tests/test_861_monomial_span_zero.py`: +- **Aux enclosure.** `_monomial_aux_bounds` assumed monotonicity (endpoint + `min`/`max`), which on a straddling box floors `x²` at `min(l²,u²) > 0` and cuts + off the true point `x=0` — latent-unsound, unreachable only because of the gate. + It now reproduces `Interval.__pow__` exactly (exact square image for `p=2`, + repeated interval multiplication for `p≥3`), which is what the cold build takes + its aux bound from. Matching rather than tightening is the point: a tighter aux + bound is a *different* relaxation, and this path must be bound-neutral. +- **Validation comparison.** `_validate`'s literal row-set equality could not + validate a *pinned* box (`lb==ub`, reached whenever integer branching fixes a + variable): the cold build emits no 1-D envelope rows at zero width while the fixed + pattern must fill its four reserved rows with the tangents/secant collapsed at the + pinned point. Those rows are exactly tight, i.e. vacuous over the box, so + `_rowset` now drops rows whose maximum over the column box already satisfies them. + Dropping a row that cannot cut the box provably leaves the polytope unchanged, so + the gate stays an exact polytope-identity test — it just stops requiring that two + identical polytopes be spelled with the same redundant rows. (This also *closes* a + pre-existing blind spot: the degenerate-box case was argued sound in a comment and + never actually validated, for sign-definite roots either.) + +*Measured (2026-07-26).* In-repo corpus panel (66 `.nl`, structure admission + root +LP bound, flag OFF vs ON): **17 already-admitted instances, 0 root-bound drift, 0 +regressed to declined, 1 newly admitted (`st_miqp5`)** — which then solves to +`-333.8888902720728`, `gap_certified=True`, `node_count=1`, incumbent independently +feasibility-verified, identical to the pre-change result. On the ball_mk2_30 *class* +(n integers straddling zero, one sign-mixed row, minimize) `solve_lp_spatial_bb(..., +require_incremental=True)` went from returning `None` (n=10 and n=30) to +`status='optimal', objective=0.0, bound=0.0, gap=0.0`. Bound-neutrality was checked +directly: patched LP value vs the same-flag cold build over 62 (model, box) pairs on +straddling boxes — max |Δ| `3.8e-15`, 0 drifts — and the incremental bound never +exceeded the richer cold path. + +*Residual (revised after review of PR #873 — the real instance, not the proxy).* +Two things remain, and the first was found only by running the **real** +`ball_mk2_30` rather than the synthetic class probe (the #727 lesson again): + +1. **Admission is necessary but NOT sufficient on ball_mk2_30.** With the structure + admitted the engine now produces a sound dual bound where it previously produced + nothing — but it still returns **no incumbent**: `objective=None`, bound −24.90, + 208,067 nodes at a 300 s budget (−25.89 / 79,843 nodes at 120 s), budget honoured + in both. The issue's premise "the primal is not the problem; admitting the model + is" is therefore **falsified**: the thin shell means no node LP rounds to the + origin. #861's relaxation-coverage half is closed; its stated symptom is not. The + remainder is primal work (the #844 family). A faithful reconstruction — + `min -Σxᵢ s.t. Σ(xᵢ² − 0.995825xᵢ) ≤ 0`, `xᵢ ∈ {−1,0,1}` — is pinned in + `test_861_monomial_span_zero.py` and reproduces this exactly; the earlier synthetic + probe could not, because its optimum sits at the origin and the root LP finds it. + *Side effect, and what it turned out to be.* ball_mk2_30 previously declined in + 0.5 s under `require_incremental` (#858) and now consumes the whole fallback + reserve for no primal. The obvious response — re-tighten the guard so it excludes + the instance again — was **rejected on the evidence**, twice over: + - *No predicate can decide this in advance.* The guard would have to answer "will + a primal appear within the reserve", and the same panel that shows `ball_mk2_30` + (n=30) finding nothing shows `ball_mk2_20` finding the exact optimum `-0.0`. + Same family, same shell, opposite outcome. A "give up if no incumbent by X% of + the reserve" rule would forfeit the late incumbents this fallback exists to + catch (tln4/tln5), and tuning X until ball_mk2_30 exits early is a + single-instance fix, which §2 rejects. + - *The reserve was not producing nothing — it was producing something that got + thrown away.* `Model.solve`'s fallback merged the engine's result only when it + carried an objective, so a sound dual bound computed inside the reserve was + discarded. Measured at a 30 s budget: the fallback proved `bound = -27.88` + (≤ the 0.0 oracle) and the solve reported `bound=None`. + + So the bound merge now runs whether or not a primal was found. Sound by the same + argument the merge already relied on: `_is_in_scope` admits only MINIMIZE, both + values are valid lower bounds, and the max of two valid lower bounds is valid; no + certificate is claimed without an incumbent (`gap_certified` untouched). + Differential solve panel over the 12 fallback-reachable in-repo instances plus + `ball_mk2_{20,30}`, before vs after: **13/14 byte-identical (status, objective, + bound, gap_certified), 1 changed — `ball_mk2_30` `bound: None → -27.88` — 0 + certification regressions, 0 objective drift, 0 bounds above an incumbent or above + a known oracle, 0 bounds loosened.** A spent reserve now buys a dual bound instead + of nothing; closing the primal half remains #844 work. +2. An **odd** power whose root box straddles zero is still declined (`ok=False` → the + trusted cold build). Closing it needs a fixed row count across the 2-facet/4-row + switch, i.e. a change to the cold build's emission, not to the patch — out of + scope for #861 and not worth the churn on the trusted path. + +*Also fixed under review (non-finite endpoints).* The `p >= 3` closed form returned +`(NaN, NaN)` on an unbounded box (`0 * ±inf`), and `_monomial_rows` put NaN +coefficients into `A`/`b` there. Both are reachable on a **sign-definite** root +(`[-inf, 0]` carrying `x**4`), i.e. on models admitted before *and* after #861, and +the validation gate never sees them because all six validation boxes are finite. NaN +in an LP is worse than a loose bound — every comparison against it is `False`, so +`NaN <= incumbent` can silently disable fathoming. The enclosure now delegates to +`Interval.__pow__` when an endpoint is non-finite (the closed form cannot match it +there: `Interval` outward-rounds per step, so `[0,inf]**3` is `[-inf,inf]`, and +claiming the tighter `[0,inf]` would make the patched aux tighter than the cold +build), and the envelope emits four vacuous rows — which is exactly the cold build's +polytope, since `_emit_1d` bails under `_finite` and emits nothing. + > **STOP / ESCALATED (2026-07-02, §0.4 + §0.6) — Phase 1's bound-neutrality > premise is falsified; needs a maintainer decision before any T1.2 code lands.** > diff --git a/python/discopt/_jax/incremental_mccormick.py b/python/discopt/_jax/incremental_mccormick.py index a47bbec1..6a0bc087 100644 --- a/python/discopt/_jax/incremental_mccormick.py +++ b/python/discopt/_jax/incremental_mccormick.py @@ -18,14 +18,19 @@ only its speed. Scope: the box-dependent terms it regenerates are bilinear products ``w=x_i*x_j`` -(4 McCormick rows) and integer squares ``s=x_i**2`` (2 endpoint tangents + 1 -secant, matching an empty ``DiscretizationState``). Any other lifted term (trilinear, -univariate, fractional power, piecewise) makes validation fail -> fallback. +(4 McCormick rows), integer powers ``s=x_i**p`` (secant + 3 tangents, matching an +empty ``DiscretizationState``) and affine squares ``w=(c*x_j+d)**2``. A monomial is +mapped on any box when ``p`` is even (convex on all of R) and on a sign-definite box +when ``p`` is odd; an odd power whose ROOT box straddles zero is declined, because +there the envelope's facet *count* changes with the node's sign regime. Any other +lifted term (trilinear, univariate, fractional power, piecewise) makes validation +fail -> fallback. """ from __future__ import annotations import logging +import math import time import numpy as np @@ -93,20 +98,46 @@ def _square_aux_bounds(li, ui): def _monomial_rows(li, ui, p): - """The 4 envelope rows for ``s = x**p`` over a *sign-definite* box ``[li,ui]`` - (secant + tangents at ``li``, the midpoint, and ``ui``), each - ``(coeff_on_x, coeff_on_s, rhs)`` of an ``... <= rhs`` row. Generalizes - :func:`_square_rows` (p=2) to any integer power ``p >= 2``. - - On a sign-definite box ``x**p`` is monotone and single-convexity: convex when - ``p`` is even or ``li >= 0``; concave when ``p`` is odd and ``ui <= 0``. In the - convex case the three tangents underestimate and the secant overestimates; in - the concave case the roles flip. Matches the uniform engine's ``_emit_1d`` - envelope exactly — which underestimates a convex atom with tangents at both - endpoints AND the box midpoint (a tighter, 3-tangent hull), so the closed-form - patch must emit the same midpoint tangent to reproduce the cold build - row-for-row (validated). + """The 4 envelope rows for ``s = x**p`` over ``[li,ui]`` (secant + tangents at + ``li``, the midpoint, and ``ui``), each ``(coeff_on_x, coeff_on_s, rhs)`` of an + ``... <= rhs`` row. Generalizes :func:`_square_rows` (p=2) to any integer power + ``p >= 2``. + + Validity domain. ``x**p`` has *single* curvature — hence a secant/tangent + envelope — on: + + * any finite box when ``p`` is EVEN (``f'' = p(p-1)x**(p-2) >= 0`` everywhere, + so the atom is convex on all of R including boxes that straddle zero); or + * a sign-definite box when ``p`` is odd (convex for ``li >= 0``, concave for + ``ui <= 0``). + + In the convex case the three tangents underestimate and the secant + overestimates; in the concave case the roles flip. On a sign-STRADDLING box with + ODD ``p`` the atom is S-shaped and neither branch is valid — the caller must not + reach here (``_build_structure`` declines that monomial; the cold build emits a + different, 2-facet hull there instead). + + Matches the uniform engine's ``_emit_1d`` envelope exactly — which underestimates + a convex atom with tangents at both endpoints AND the box midpoint (a tighter, + 3-tangent hull), so the closed-form patch must emit the same midpoint tangent to + reproduce the cold build row-for-row (validated). """ + # Unbounded box: there is no secant and no finite tangent point, and evaluating + # the formulae anyway yields NaN coefficients (``inf - inf`` in the tangent + # intercept, ``0 * inf`` in the slope) that flow straight into ``A``/``b``. NaN in + # an LP row is strictly worse than a missing row — comparisons against it are all + # False, so it can silently disable fathoming rather than merely loosen a bound. + # Emit four VACUOUS rows (``0 <= 0``) instead, which is exactly what the cold + # build's polytope is here: ``_emit_1d`` bails on ``not _finite(lo, hi)`` and + # emits no envelope rows at all, leaving the aux interval bound as the whole + # relaxation. So this makes the patch AGREE with the cold build on a box where it + # previously disagreed, and ``_rowset``'s vacuity filter drops these rows on the + # patched side exactly as the cold side has none. (Latent on main too, and not + # reachable from the newly-admitted class — every validation box is finite — but + # it is the same non-finite-endpoint defect as the aux enclosure below, found by + # the regression test written for it.) + if not (math.isfinite(li) and math.isfinite(ui)): + return [(0.0, 0.0, 0.0)] * 4 mid = 0.5 * (li + ui) fl, fm, fu = li**p, mid**p, ui**p dfl, dfm, dfu = p * li ** (p - 1), p * mid ** (p - 1), p * ui ** (p - 1) @@ -135,9 +166,66 @@ def _monomial_rows(li, ui, p): def _monomial_aux_bounds(li, ui, p): - """min/max of ``x**p`` over a sign-definite ``[li,ui]`` (monotone there).""" - a, b = li**p, ui**p - return (a, b) if a <= b else (b, a) + """Enclosure of ``x**p`` over ``[li,ui]``, matching the cold build's aux column. + + The cold build takes this bound from ``evaluate_interval``, i.e. from + ``Interval.__pow__``, so this closed form must reproduce *that* — not merely a + sound enclosure — or the incremental rows describe a different polytope than the + path they must be bound-neutral against (``_validate`` compares aux bounds too). + ``Interval.__pow__`` is: + + * ``p == 2``: the exact square image ``[0, max(li^2,ui^2)]`` when the box + straddles zero, else ``[min, max]`` of the endpoint squares; + * ``p >= 3``: repeated interval MULTIPLICATION by ``[li,ui]`` (``p-1`` times), + whose corner min/max is reproduced below. That is deliberately *looser* than + the exact image on a straddling box (e.g. ``x**4`` over ``[-2,3]`` encloses + ``[-54, 81]``, not ``[0, 81]``) — matching it exactly is what keeps the fast + path bound-NEUTRAL. Do not "tighten" this: a tighter aux bound here is a + different relaxation, which is precisely what the incremental path may not be. + + A NON-FINITE endpoint delegates to ``Interval.__pow__`` itself rather than being + reproduced here, because the closed form below cannot match it there. Two ways it + diverges, both load-bearing: + + * ``0 * ±inf`` is ``NaN`` in IEEE and a bare ``min``/``max`` propagates it, + collapsing the enclosure to ``[NaN, NaN]``. That is worse than a wide bound: it + reaches the LP, and ``NaN <= incumbent`` is always ``False``, so it can + silently disable fathoming. ``Interval.__mul__`` maps that corner to 0 (C-36 / + #723). + * ``Interval`` outward-rounds after *every* step, and on an unbounded box that + changes the answer rather than the last ulp: ``[0,inf]**2`` rounds its lower + end to ``-5e-324``, which at the next multiply gives ``-inf`` instead of a + NaN-to-zero corner — so ``[0,inf]**3`` is ``[-inf, inf]``, not ``[0, inf]``. + Reporting the tighter one would make the patched aux TIGHTER than the cold + build, which this function's whole contract forbids. + + Delegating costs ~65 us versus ~2 us for the loop, so it is confined to the + non-finite case; a finite box (every real node box, and every box the corpus + exercises) keeps the fast path, where the loop reproduces ``Interval.__pow__`` + to well within the comparison tolerance (its per-step rounding is a 1-ulp + effect there, not a change of value). + + The old form (``min/max`` of the two endpoint powers) was correct only on a + sign-definite box, where ``x**p`` is monotone; it was UNSOUND on a straddling box + for even ``p`` (it would floor ``x**2`` at ``min(li^2,ui^2) > 0`` and cut off the + true point ``x=0``). That was unreachable while every monomial was gated to a + sign-definite root box, and is the reason the bound had to be generalized before + the gate could be relaxed (#861). ``test_monomial_aux_bounds_match_interval_pow`` + pins the parity with ``Interval.__pow__`` across powers and sign regimes, + half-infinite and doubly-infinite boxes included. + """ + if not (math.isfinite(li) and math.isfinite(ui)): + from discopt._jax.convexity.interval import Interval + + enc = Interval.from_bounds(np.float64(li), np.float64(ui)) ** int(p) + return float(enc.lo), float(enc.hi) + if p == 2: + return _square_aux_bounds(li, ui) + lo, hi = li, ui + for _ in range(p - 1): + corners = (lo * li, lo * ui, hi * li, hi * ui) + lo, hi = min(corners), max(corners) + return lo, hi def _affine_square_rows(coeff, const, li, ui): @@ -248,12 +336,15 @@ def _full_build(self, lb, ub): def _build_structure(self): n = len(self.model._variables) - # Per-variable ROOT sign regime (cert:T1.2). A monomial ``x**p`` has a - # box-*sign*-dependent row structure (3 rows on a sign-definite box, 4/2 - # when the box strictly spans zero), so it can be patched only when the - # variable's root box is sign-definite — which branching preserves, since - # it only shrinks boxes. ``+1`` = ``lb>=0``, ``-1`` = ``ub<=0``, ``0`` = - # spans zero (any monomial on such a var is unmappable below). + # Per-variable ROOT sign regime (cert:T1.2). ``+1`` = ``lb>=0``, ``-1`` = + # ``ub<=0``, ``0`` = spans zero. Branching only shrinks boxes, so a + # sign-definite root sign holds at every node of the subtree; a spanning root + # reaches BOTH regimes below it. An ODD-power monomial's envelope changes + # *shape* across that split — 4 secant/tangent rows on a sign-definite box vs + # the 2-facet S-hull on a straddling one — which a fixed sparsity pattern + # cannot express, so those stay unmappable (#861). An EVEN power is convex on + # all of R, so its envelope is the same 4 rows in every regime and it is + # mapped regardless of root sign. root_lb = np.array([float(np.min(v.lb)) for v in self.model._variables]) root_ub = np.array([float(np.max(v.ub)) for v in self.model._variables]) self._root_sign = np.where(root_lb >= 0.0, 1, np.where(root_ub <= 0.0, -1, 0)) @@ -326,11 +417,19 @@ def _support(k): if len(rows) != 4: raise ValueError(f"bilinear ({i},{j}) -> {len(rows)} rows, expected 4") self.bilin_rows[(i, j, a)] = rows - # monomial x_i**p, any p >= 2, gated on a sign-definite root box. + # monomial x_i**p, any p >= 2. Only an ODD power needs a sign-definite root + # box (see the ``_root_sign`` note above); an even power is convex on all of + # R and keeps its 4-row envelope across a sign change, so it is admitted on a + # straddling root too — the case that used to decline the whole structure on + # sign-mixed integer QCQPs (#861). self.mono_rows = {} for (i, p), a in self.monomial.items(): - if self._root_sign[i] == 0: - raise ValueError(f"monomial x_{i}^{p}: root box spans zero (unmappable)") + if self._root_sign[i] == 0 and p % 2 == 1: + raise ValueError( + f"monomial x_{i}^{p}: odd power on a root box spanning zero " + "(the envelope switches between the 4-row secant/tangent hull and " + "the 2-facet S-hull, which the fixed row pattern cannot express)" + ) rows = [int(k) for k in _rows_with_col(a) if _support(k) <= {i, a}] if len(rows) != 4: raise ValueError(f"monomial x_{i}^{p} -> {len(rows)} rows, expected 4") @@ -443,7 +542,7 @@ def _patch(self, lb, ub): # -- soundness gate ---------------------------------------------------- # @staticmethod - def _rowset(A, b): + def _rowset(A, b, bounds=None): """Canonical hashable representation of the polytope's rows (order-free). Sparse-native (O(nnz), never densified): each row is the sorted tuple of its @@ -451,13 +550,58 @@ def _rowset(A, b): to 0 are dropped, so an explicit structural zero (the fixed-pattern ``_patch`` can leave a zeroed target entry) compares equal to its absence in the cold build — the two matrices match iff they encode the same polytope. + + ``bounds`` (the ``(ncol,2)`` column box) additionally drops rows that are + **vacuous over that box**: ``max_{x in box} a·x <= rhs``, i.e. the row cuts + off nothing the bounds do not already exclude. Removing such a row provably + leaves the feasible set unchanged, so the comparison stays an exact + polytope-identity test — it just stops demanding that two identical polytopes + be spelled with the same redundant rows. This is what lets a *pinned* box + (``lb==ub``, reached whenever integer branching fixes a variable) validate: + the cold build emits no 1-D envelope rows at zero width (``_emit_1d`` bails + under ``_MIN_WIDTH``) and pins the aux via its ``[f(v),f(v)]`` bound, while + the fixed-pattern ``_patch`` must write *something* into its four reserved + rows and writes the tangents/secant collapsed at ``v`` — which are exactly + tight there, hence vacuous, hence dropped here. Any row that genuinely cuts + the box survives on both sides and a real mismatch is still caught. """ M = sp.csr_matrix(A) M.sort_indices() indptr, indices, data = M.indptr, M.indices, M.data b = np.asarray(b, dtype=np.float64).ravel() + vacuous = None + if bounds is not None: + # Row maximum over the box, VECTORIZED over all nonzeros at once (this + # runs on the whole lift — up to ~172k nnz on qap — twice per validation + # box, so a per-row numpy call would re-introduce the #654 pre-B&B + # overrun). Each term contributes its larger endpoint product; an exact + # zero coefficient contributes nothing (guarding ``0 * inf -> nan`` on an + # unbounded column). A non-finite maximum leaves the row in place, which + # is the conservative direction: an undroppable row can only cause a + # (sound) mismatch, never a false match. + bnd = np.asarray(bounds, dtype=np.float64) + col_lo, col_hi = bnd[:, 0], bnd[:, 1] + contrib = np.where( + data == 0.0, 0.0, np.maximum(data * col_lo[indices], data * col_hi[indices]) + ) + row_max = np.asarray( + sp.csr_matrix((contrib, indices, indptr), shape=M.shape).sum(axis=1) + ).ravel() + # The slack exists only to absorb the last-ulp difference between a row + # built by the patch and the same row built by the cold build; it is NOT + # meant to forgive a row that genuinely cuts. Measured over 22 admitted + # models x 132 validation boxes, the worst amount by which a DROPPED row + # actually cut its box was 7.1e-15 — i.e. the tolerance runs ~6 orders of + # magnitude above anything the corpus exercises. If that margin ever + # closes, this is the number to re-measure: the relative form would + # tolerate a ~1-unit cut on a row whose maximum is ~1e9, so a future + # large-coefficient lift is where it would first matter. + slack = 1e-9 * (1.0 + np.abs(b) + np.abs(row_max)) + vacuous = np.isfinite(row_max) & (row_max <= b + slack) out = [] for k in range(M.shape[0]): + if vacuous is not None and vacuous[k]: + continue entries = tuple( (int(indices[t]), rv) for t in range(indptr[k], indptr[k + 1]) @@ -488,12 +632,15 @@ def _validation_boxes(self): Every box is a *reachable* B&B sub-box of the root: branching only shrinks a box, so a var that is sign-definite at the root (``_root_sign != 0``) keeps that sign — a positive var never gets ``lb<0``, a negative var never gets - ``ub>0``. A **spanning** var (``_root_sign==0``), however, carries no - monomial (gated out in :meth:`_build_structure`) and its real nodes DO carry - negative / zero-spanning bounds, so the boxes below deliberately drive those - vars through negative-lb, zero-spanning (``lb<0=0``-only set never exercised (C-21). + ``ub>0``. A **spanning** var (``_root_sign==0``) carries only even-power + monomials (odd ones are gated out in :meth:`_build_structure`) and its real + nodes DO carry negative / zero-spanning bounds, so the boxes below + deliberately drive those vars through negative-lb, zero-spanning + (``lb<0=0``-only set + never exercised (C-21). Since #861 those same boxes are what proves an + even-power envelope on a straddling box reproduces the cold build: the + ``span``/``span_wide`` trials put ``lb<0 bool: enough that a single one runs past the top-of-loop deadline poll. On tln5 at a 21 s budget: 5 nodes in 43.8 s (2.08x, slowest node 42.4 s) versus 13158 nodes in 21.0 s (1.00x, slowest node 0.04 s) once the engine passes its own deadline. - ball_mk2_30 declines the structure for a genuine structural reason (``monomial - x_0^2: root box spans zero``); there the fallback now refuses to run at all - (``require_incremental``) rather than spend 61 s on a single root LP for nothing. + A model whose structure declines for a genuine structural reason — e.g. an odd + power over a root box straddling zero, whose envelope changes facet count across + the sign change — makes the fallback refuse to run at all + (``require_incremental``) rather than spend a whole reserve on a single root LP + for nothing; measured at 61 s against a 21 s reserve on ball_mk2_30 back when its + ``x_0**2`` monomial also declined. #861 has since admitted the even powers, so + ball_mk2_30 now runs this path and produces a sound bound but still no incumbent — + the remaining gap there is primal, not relaxation coverage. **Known cost.** The incumbents are *worse* than the ones the degraded cold path happened to report (tln4 19.6 vs 8.7, tln5 32.8 vs 15.1) — its slower, full- @@ -4169,12 +4174,29 @@ def solve( # # require_incremental=True: without the incremental McCormick # structure the engine has no cuts, no feasibility pump, and - # rebuilds the relaxation per node. On ball_mk2_30 — where the - # structure legitimately declines (``monomial x_0^2: root box - # spans zero``) — that cold path spent 61 s on the ROOT LP alone - # against a 21 s reserve: 0 nodes, no incumbent, 2.91x over - # budget. It cannot produce a primal inside a fallback-sized - # budget, so declining costs no gain and removes the overrun. + # rebuilds the relaxation per node. Measured on ball_mk2_30 while + # its monomial still declined, that cold path spent 61 s on the + # ROOT LP alone against a 21 s reserve: 0 nodes, no incumbent, + # 2.91x over budget. It cannot produce a primal inside a + # fallback-sized budget, so declining costs no gain and removes + # the overrun. ball_mk2_30 itself now MAPS (#861 narrowed the + # monomial gate to odd powers on a straddling root) and so takes + # the incremental path here — but measured, it still returns no + # incumbent, just a sound bound, having spent the reserve. + # + # The guard is deliberately NOT tightened to re-exclude it. Its + # premise ("the structure builds" ⇒ "this path can find a primal") + # is a proxy, and #861 widened the gap between the two — but no + # predicate can decide in advance whether a primal is coming, and + # a "give up if no incumbent by X% of the reserve" rule would + # forfeit precisely the late incumbents this fallback exists to + # catch (tln4/tln5 are found deep in the reserve, not early). + # Tuning such a rule until ball_mk2_30 exits early would be a + # single-instance fix, which this repo rejects. Instead the + # reserve's OUTPUT is no longer discarded: the bound merge below + # runs whether or not a primal was found, so a spent reserve now + # buys a tighter dual bound rather than nothing. Closing the primal + # half remains #844 work. _fb = solve_lp_spatial_bb( self, time_limit=_fb_reserve, @@ -4182,21 +4204,39 @@ def solve( use_obbt=False, require_incremental=True, ) + if _fb is not None: + # Keep the TIGHTER of the two dual bounds, and do it whether or not + # the fallback found a primal. The engine's frontier bound is a + # valid LOWER bound for a minimize (``_is_in_scope`` admits only + # minimize), so the max of two valid lower bounds is a valid lower + # bound — the same merge this branch has always done, just no longer + # conditioned on an incumbent it does not depend on. + # + # Why it moved (#861 review): once a model is admitted but hard, the + # fallback spends its whole reserve, computes a sound bound, finds no + # primal, and everything was thrown away — the reserve produced + # nothing *even though it produced something valid*. Measured on + # ball_mk2_30 at a 30 s budget: the fallback returned bound -27.88 + # (<= the 0.0 oracle) and the solve reported ``bound=None``. That is + # the real cost the review flagged as "spends the budget for + # nothing", and it is a reporting gap, not a reason to re-tighten + # ``require_incremental``: no predicate can decide in advance whether + # a primal is coming, and a budget-fraction early exit would forfeit + # exactly the late incumbents this fallback exists to catch. + if _fb.bound is not None: + result.bound = ( + _fb.bound if result.bound is None else max(result.bound, _fb.bound) + ) + result.node_count = (result.node_count or 0) + _fb.node_count if _fb is not None and _fb.objective is not None: from discopt.solver import _unpack_solution result.status = _fb.status result.objective = _fb.objective result.x = _unpack_solution(self, np.asarray(_fb.x)) - # Keep the TIGHTER of the two dual bounds and never claim a - # certificate the fallback did not actually prove. - if _fb.bound is not None: - result.bound = ( - _fb.bound if result.bound is None else max(result.bound, _fb.bound) - ) result.gap = _fb.gap + # Never claim a certificate the fallback did not actually prove. result.gap_certified = _fb.status == "optimal" - result.node_count = (result.node_count or 0) + _fb.node_count except Exception as _fb_exc: # Never break a solve -- but never swallow silently either. A bare # ``except`` here previously turned a hard TypeError into an invisible diff --git a/python/tests/test_861_monomial_span_zero.py b/python/tests/test_861_monomial_span_zero.py new file mode 100644 index 00000000..006af64f --- /dev/null +++ b/python/tests/test_861_monomial_span_zero.py @@ -0,0 +1,451 @@ +"""#861 — the monomial envelope no longer declines when the root box spans zero. + +``IncrementalMcCormickLP`` gated EVERY monomial ``x_i**p`` on a sign-definite root +box, so a model whose integers straddle zero (ball_mk2_30: 30 integers, one +sign-mixed "thin shell" row) declined the whole incremental structure — and with it +the cuts, the feasibility pump and, under ``require_incremental=True``, the entire +LP-per-node solve, which then returned no incumbent at all. + +The gate was wider than the mathematics. ``x**p`` for EVEN ``p`` has +``f'' = p(p-1)x**(p-2) >= 0`` on all of R, so it is convex across a sign change and +the cold build emits the *same* 4-row secant/tangent envelope in every sign regime — +measured on ``build_milp_relaxation`` before this change: + + p=2 box=[-2,3] -> 4 rows p=2 box=[1,3] -> 4 rows p=2 box=[-3,-1] -> 4 rows + p=4 box=[-2,3] -> 4 rows p=4 box=[1,3] -> 4 rows p=4 box=[-3,-1] -> 4 rows + p=3 box=[-2,3] -> 2 rows p=3 box=[1,3] -> 4 rows p=3 box=[-3,-1] -> 4 rows + +Only the ODD powers change facet COUNT across the sign change (the S-shaped atom's +2-facet hull vs the 4-row envelope), and only those are still unmappable by a fixed +sparsity pattern. So even powers are admitted on any root box; odd powers keep the +sign-definite requirement. + +Two things had to be generalized before the gate could move, and both are pinned +here: the aux-column enclosure (the old endpoint ``min``/``max`` assumed monotonicity +and would have FLOORED ``x**2`` above zero on a straddling box — cutting off the true +point ``x=0``), and the validation gate's row comparison (a *pinned* box, reached +whenever integer branching fixes a variable, gets no envelope rows from the cold +build but four exactly-tight — hence vacuous — rows from the fixed-pattern patch). +""" + +from __future__ import annotations + +import os + +os.environ.setdefault("JAX_PLATFORMS", "cpu") +os.environ.setdefault("JAX_ENABLE_X64", "1") + +import itertools + +import discopt.modeling as dm +import numpy as np +import pytest +from discopt._jax.incremental_mccormick import IncrementalMcCormickLP, _monomial_aux_bounds +from discopt._jax.lp_spatial_bb import solve_lp_spatial_bb +from discopt._jax.term_classifier import classify_nonlinear_terms + + +def _monomial_model(p: int, lo: float, hi: float, n: int = 3): + """``n`` integers over ``[lo,hi]`` carrying bare ``x_i**p`` monomials.""" + m = dm.Model(f"mono_p{p}") + xs = [m.integer(f"x{i}", lb=lo, ub=hi) for i in range(n)] + m.minimize(sum(x**p for x in xs)) + m.subject_to(sum(x**p for x in xs) >= 1) + return m + + +def _ball_mk2_class(n: int = 30): + """The ball_mk2_30 class: ``n`` integers whose root boxes straddle zero, a single + sign-mixed row over all of them, MINIMIZE, optimum 0 at the origin. Named + instances are gate probes only — this is the *shape* that #861 declined.""" + m = dm.Model("ball_mk2_class") + xs = [m.integer(f"x{i}", lb=-1, ub=1) for i in range(n)] + m.minimize(sum(x * x for x in xs)) + shell = sum((-1.0) ** i * x for i, x in enumerate(xs)) + m.subject_to(shell <= 1) + m.subject_to(shell >= -1) + return m + + +def _ball_mk2_real(n: int = 30): + """Faithful reconstruction of MINLPLib's ``ball_mk2_30`` — the instance #861 was + actually filed against: + + min -Σ xᵢ s.t. Σ (xᵢ² - 0.995825·xᵢ) ≤ 0, xᵢ ∈ {-1,0,1} + + This is NOT the same problem as :func:`_ball_mk2_class`, and the difference is the + point. Every term ``x² - 0.995825x`` is ≥ 0 at an integer ``x`` with equality only + at ``x=0``, so the origin is the *only* feasible integer point and the optimum is + 0.0 (matching ``minlplib.solu``) — but the objective ``-Σ xᵢ`` pulls the + relaxation the other way, to fractional ``xᵢ ≈ 0.995825`` where the shell is + slack. That is the "thin shell": the LP optimum is nowhere near the integer one, + which is what makes the primal hard here and trivial in ``_ball_mk2_class`` + (whose optimum sits at the origin and is found by the root LP). + + Reconstructing it matters because ``_ball_mk2_class`` alone cannot detect the real + failure — the #727 RLT lesson in CLAUDE.md: a mechanism validated only on a + synthetic proxy can be a no-op on the real class. + """ + m = dm.Model("ball_mk2_real") + xs = [m.integer(f"x{i}", lb=-1, ub=1) for i in range(n)] + m.minimize(-sum(xs)) + m.subject_to(sum(x * x - 0.995825 * x for x in xs) <= 0.0) + return m + + +def _structure(model): + return IncrementalMcCormickLP(model, classify_nonlinear_terms(model)) + + +# --------------------------------------------------------------------------- # +# The gate itself +# --------------------------------------------------------------------------- # + + +@pytest.mark.parametrize("p", [2, 4, 6]) +def test_even_power_monomial_maps_on_a_root_box_spanning_zero(p): + """The #861 regression: before the fix this raised ``monomial x_0^{p}: root box + spans zero (unmappable)`` and left ``ok=False``.""" + assert _structure(_monomial_model(p, -2, 2)).ok + + +@pytest.mark.parametrize("p", [2, 3, 4, 5]) +@pytest.mark.parametrize("box", [(0, 3), (-3, 0)]) +def test_sign_definite_root_still_maps_for_every_power(p, box): + """Unchanged behaviour on the regime that already worked.""" + assert _structure(_monomial_model(p, *box)).ok + + +@pytest.mark.parametrize("p", [3, 5]) +def test_odd_power_monomial_still_declines_on_a_root_box_spanning_zero(p): + """An odd power's envelope switches between the 4-row secant/tangent hull and the + 2-facet S-hull across zero — a facet-COUNT change the fixed sparsity pattern + cannot express, so it must keep declining (soundly: the caller cold-builds).""" + assert not _structure(_monomial_model(p, -2, 2)).ok + + +# --------------------------------------------------------------------------- # +# Aux-column enclosure parity (the soundness prerequisite) +# --------------------------------------------------------------------------- # + + +@pytest.mark.parametrize("p", [2, 3, 4, 5, 6, 7, 8]) +@pytest.mark.parametrize( + "box", + [ + (-2.0, 3.0), + (1.0, 3.0), + (-3.0, -1.0), + (0.0, 2.5), + (-2.5, 0.0), + (-1.5, -1.5), + # Half-infinite and doubly-infinite boxes. These are where a hand-rolled + # closed form diverges from Interval: `0 * ±inf` is NaN in IEEE, and + # Interval's per-step outward rounding turns an exact 0 endpoint into a + # denormal that then multiplies to ±inf rather than to a NaN-to-zero corner. + # `_monomial_aux_bounds` delegates here rather than reproducing it. + (-np.inf, 0.0), + (0.0, np.inf), + (-np.inf, np.inf), + (-np.inf, -1.0), + (1.0, np.inf), + (-2.0, np.inf), + (-np.inf, 2.0), + ], +) +def test_monomial_aux_bounds_match_interval_pow(p, box): + """``_monomial_aux_bounds`` must reproduce the enclosure the COLD build takes from + ``Interval.__pow__`` — not merely a sound one — or the two paths describe + different polytopes. Pins the parity so a change to the interval arithmetic + surfaces here rather than as a silent bound difference.""" + from discopt._jax.convexity.interval import Interval + + lo, hi = box + ref = Interval.from_bounds(np.array([lo]), np.array([hi])) ** p + got = _monomial_aux_bounds(lo, hi, p) + assert not any(g != g for g in got), f"NaN enclosure on {box}, p={p}: {got}" + for g, t in zip(got, (float(ref.lo[0]), float(ref.hi[0]))): + if g == t: # covers the ±inf endpoints, which approx() cannot compare + continue + assert g == pytest.approx(t, rel=1e-9, abs=1e-9) + + +@pytest.mark.parametrize("p", [3, 4, 5, 6]) +def test_unbounded_box_puts_no_nan_anywhere_in_the_node_lp(p): + """Regression for the NaN family on an unbounded box (caught in review of this + PR): ``(nan, nan)`` from the ``p >= 3`` aux enclosure, and NaN coefficients from + the envelope rows. + + The variable below is sign-DEFINITE (``ub <= 0`` → ``_root_sign = -1``), so the + structure was admitted before *and* after #861 and all six validation boxes are + finite — the gate never sees the unbounded node box. NaN reaching the LP is worse + than a loose bound: every comparison against it is ``False``, so ``NaN <= + incumbent`` can silently disable fathoming. Assert on the WHOLE assembled LP, not + just the aux bounds — the aux fix alone left 4 NaNs in ``A``/``b``. + """ + m = dm.Model("halfinf") + x = m.integer("x", lb=-np.inf, ub=0) + y = m.integer("y", lb=-5, ub=0) + m.minimize(x**p + y**p) + m.subject_to(x + y >= -8) + inc = _structure(m) + assert inc.ok + A, b, bounds = inc.assemble(np.array([-np.inf, -5.0]), np.array([0.0, 0.0])) + assert not np.isnan(bounds).any(), f"NaN in aux bounds for p={p}: {bounds}" + assert not np.isnan(A.data).any(), f"NaN in constraint matrix for p={p}" + assert not np.isnan(b).any(), f"NaN in rhs for p={p}" + + +@pytest.mark.parametrize("p", [2, 3, 4]) +def test_unbounded_box_envelope_matches_the_cold_build_emptiness(p): + """On a non-finite box the cold build emits NO envelope rows (``_emit_1d`` bails + under ``_finite``), leaving the aux interval bound as the entire relaxation. The + fixed-pattern patch cannot delete rows, so it must fill them with VACUOUS ones — + which describes the same polytope. Pins that they carry no coefficients at all.""" + from discopt._jax.incremental_mccormick import _monomial_rows + + for box in [(-np.inf, 0.0), (0.0, np.inf), (-np.inf, np.inf), (1.0, np.inf)]: + rows = _monomial_rows(box[0], box[1], p) + assert len(rows) == 4 + for cx, cs, rhs in rows: + assert (cx, cs, rhs) == (0.0, 0.0, 0.0), f"non-vacuous row on {box}: {rows}" + + +def test_even_power_aux_floor_admits_the_origin_on_a_straddling_box(): + """The old endpoint-``min`` form returned ``min(li^2, ui^2) = 1`` on ``[-1, 3]``, + which floors ``x**2`` above its true value at ``x=0`` and cuts off a feasible + point. The enclosure must reach 0 wherever the box contains 0.""" + lo, hi = _monomial_aux_bounds(-1.0, 3.0, 2) + assert lo <= 0.0 <= hi + assert hi == pytest.approx(9.0) + + +# --------------------------------------------------------------------------- # +# Soundness of the admitted envelope +# --------------------------------------------------------------------------- # + + +@pytest.mark.parametrize("p", [2, 4]) +def test_spanning_monomial_envelope_cuts_no_feasible_point(p): + """Feasible-point sampling: every ``(x_i, x_i**p)`` with ``x_i`` in the node box + must satisfy every patched envelope row. Run over sub-boxes that straddle zero, + sit on either side of it, and pin the variable.""" + inc = _structure(_monomial_model(p, -2, 2)) + assert inc.ok + boxes = [(-2.0, 2.0), (-2.0, 0.0), (0.0, 2.0), (-1.0, 1.0), (1.0, 1.0), (-2.0, -2.0)] + checked = 0 + for lo, hi in boxes: + lb = np.full(inc.n, lo) + ub = np.full(inc.n, hi) + A, b, bounds = inc.assemble(lb, ub) + A = A.tocsr() + for (i, a, pw), rows in inc.mono_rows.items(): + aux_lo, aux_hi = bounds[a] + for v in np.linspace(lo, hi, 9): + s = v**pw + assert aux_lo - 1e-9 <= s <= aux_hi + 1e-9, f"aux bound cuts x={v}" + for k in rows: + lhs = 0.0 + for t in range(A.indptr[k], A.indptr[k + 1]): + col = int(A.indices[t]) + if col == i: + lhs += float(A.data[t]) * v + elif col == a: + lhs += float(A.data[t]) * s + assert lhs <= float(b[k]) + 1e-9, f"row {k} cuts x={v} (p={pw})" + checked += 1 + # Guard against the loops silently emptying (a changed mono_rows key, a box list + # that stops admitting): 6 boxes x 3 monomials x 9 samples x 4 rows. + assert checked == 648, f"expected 648 row assertions, executed {checked}" + + +@pytest.mark.parametrize("p", [2, 4]) +def test_spanning_monomial_patch_is_bound_neutral_against_the_cold_build(p): + """The incremental path may change speed, never the bound: the patched LP's + optimal value must equal the value of the same-flag cold build on the same box, + including boxes that straddle zero.""" + inc = _structure(_monomial_model(p, -2, 2)) + assert inc.ok + rng = np.random.default_rng(861) + compared = 0 + for _ in range(10): + lb = rng.integers(-2, 3, size=inc.n).astype(float) + ub = np.array([rng.integers(int(v), 3) for v in lb], dtype=float) + patched = inc.solve_assembled(*inc.assemble(lb, ub))[0] + Af, bf, bdf, _, _, _ = inc._full_build(lb, ub) + cold = inc.solve_assembled(Af, bf, bdf)[0] + if patched is None or cold is None: + continue + assert patched == pytest.approx(cold, rel=1e-9, abs=1e-9) + compared += 1 + assert compared >= 5, "too few comparable boxes to call this a check" + + +@pytest.mark.parametrize("p", [2, 4]) +def test_spanning_monomial_bound_never_exceeds_the_box_optimum(p): + """The relaxation is a valid lower bound: on a fixed box its LP value must not + exceed the true (brute-forced) integer optimum over that box.""" + lo, hi = -2, 2 + n = 3 + model = _monomial_model(p, lo, hi, n=n) + inc = _structure(model) + assert inc.ok + compared = 0 + for box in [(-2, 2), (-2, 0), (0, 2), (-1, 1)]: + lb = np.full(n, float(box[0])) + ub = np.full(n, float(box[1])) + bound = inc.solve_assembled(*inc.assemble(lb, ub))[0] + if bound is None: + continue + true = None + for pt in itertools.product(range(box[0], box[1] + 1), repeat=n): + v = np.array(pt, dtype=float) + if float(np.sum(v**p)) >= 1.0 - 1e-9: + obj = float(np.sum(v**p)) + true = obj if true is None else min(true, obj) + if true is not None: + assert bound <= true + 1e-6, f"bound {bound} above box optimum {true}" + compared += 1 + # Without this the test degrades to a no-op if every LP returns None or every box + # turns out infeasible (the sibling bound-neutrality test guards the same way). + assert compared == 4, ( + f"expected 4 boxes compared against a brute-forced optimum, got {compared}" + ) + + +# --------------------------------------------------------------------------- # +# The validation gate's vacuous-row filter +# --------------------------------------------------------------------------- # + + +def test_rowset_drops_only_rows_that_cannot_cut_the_box(): + """The filter that lets a *pinned* box validate must drop exactly the rows whose + maximum over the box already satisfies them, and keep every row that bites.""" + import scipy.sparse as sp + + A = sp.csr_matrix(np.array([[1.0, 1.0], [1.0, 0.0]])) + b = np.array([10.0, 0.5]) # row 0 vacuous over the box below, row 1 cuts it + bounds = np.array([[0.0, 1.0], [0.0, 2.0]]) + unfiltered = IncrementalMcCormickLP._rowset(A, b) + filtered = IncrementalMcCormickLP._rowset(A, b, bounds) + assert len(unfiltered) == 2 + assert filtered == [(((0, 1.0),), 0.5)] + + +def test_pinned_variable_box_validates_for_a_spanning_root(): + """``_validation_boxes`` drives every spanning var through a degenerate + (``lb==ub``) trial — reachable whenever integer branching fixes a variable. The + cold build emits no envelope rows at zero width while the fixed pattern must fill + its four reserved rows, so the structure only validates because those rows are + exactly tight (vacuous) there. Assert the regime is actually exercised.""" + inc = _structure(_monomial_model(2, -2, 2)) + assert inc.ok + assert {"span", "degen", "neg"} <= inc._validated_regimes + + +# --------------------------------------------------------------------------- # +# Gate probe: the class the issue was filed against +# --------------------------------------------------------------------------- # + + +@pytest.mark.parametrize("n", [10, 30]) +def test_ball_mk2_class_is_admitted_and_solved_under_require_incremental(n): + """Before #861 this returned ``None``: the structure declined on + ``monomial x_0^2: root box spans zero`` and ``require_incremental=True`` (PR #858) + turned that into "no incumbent". It must now certify the optimum, 0.0.""" + result = solve_lp_spatial_bb(_ball_mk2_class(n), time_limit=60.0, require_incremental=True) + assert result is not None, "engine still declines the ball_mk2_30 class" + assert result.status == "optimal" + assert result.objective == pytest.approx(0.0, abs=1e-6) + assert result.bound <= result.objective + 1e-6 # certificate invariant + + +def test_real_ball_mk2_30_relaxation_is_admitted(): + """The REAL instance's relaxation coverage — the gap #861 reports — is closed. + + Fails before the change (``monomial x_0^2: root box spans zero (unmappable)`` + → ``ok=False``), passes after, with all 30 monomials mapped. + """ + inc = _structure(_ball_mk2_real(30)) + assert inc.ok, "the real ball_mk2_30 relaxation still cannot be built" + assert len(inc.monomial) == 30 + assert all(inc._root_sign[i] == 0 for (i, _p) in inc.monomial), "expected straddling roots" + + +def test_real_ball_mk2_30_bound_is_sound_and_no_false_certificate(): + """On the real instance the engine now produces a *bound* where it previously + produced nothing — and that bound must never cross the 0.0 oracle + (``minlplib.solu``). + + Deliberately NOT asserted: that an incumbent is found. It is not, at any budget + tried (see ``test_real_ball_mk2_30_still_finds_no_incumbent``), so #861's stated + symptom is only partly addressed — this PR closes the relaxation-coverage half. + The assertions below are written to keep holding if a later primal fix lands. + """ + result = solve_lp_spatial_bb(_ball_mk2_real(30), time_limit=20.0, require_incremental=True) + assert result is not None, "engine declines the real ball_mk2_30" + assert result.bound is not None + assert result.bound <= 0.0 + 1e-6, f"dual bound {result.bound} crossed the 0.0 oracle" + if result.objective is not None: + # If a primal ever appears it must be feasible-valued and certificate-consistent. + assert result.objective >= 0.0 - 1e-6, "incumbent below the true optimum" + assert result.bound <= result.objective + 1e-6 + + +@pytest.mark.slow +def test_real_ball_mk2_30_still_finds_no_incumbent(): + """Pins the RESIDUAL so it cannot be mistaken for fixed, and so the day it stops + being true someone is told. + + Admitting the model was necessary but not sufficient: the relaxation now builds + and the engine explores thousands of nodes inside its budget, but the thin shell + means no node LP rounds to the origin, so ``objective`` stays ``None`` — #861's + "returns no incumbent" persists. Measured on this reconstruction: 791 nodes / + bound -27.88 at a 20 s budget, 12236 nodes / bound -26.89 at 60 s, both budgets + honoured. That is primal work (the #844 family), not relaxation coverage. + + XFAIL-shaped on purpose: it PASSES while the residual exists and XPASSes loudly + if a primal fix lands, at which point delete it and assert the optimum instead. + """ + result = solve_lp_spatial_bb(_ball_mk2_real(30), time_limit=20.0, require_incremental=True) + assert result is not None + if result.objective is not None: # pragma: no cover - the day this fires, celebrate + pytest.fail( + f"ball_mk2_30 now yields an incumbent ({result.objective}) — the #861 " + "residual is closed; replace this test with an optimality assertion." + ) + assert result.status == "time_limit" + assert result.node_count > 0, "admitted but explored no nodes — budget spent elsewhere" + + +# --------------------------------------------------------------------------- # +# The #844 fallback's dual bound is no longer discarded (raised in review of #873) +# --------------------------------------------------------------------------- # + + +@pytest.mark.slow +def test_no_incumbent_fallback_still_reports_its_dual_bound(): + """A spent fallback reserve must buy a BOUND even when it buys no primal. + + #861 admits ball_mk2_30 into the LP-per-node engine, so ``Model.solve``'s + no-incumbent fallback (#844) now runs it instead of declining under + ``require_incremental``. It spends the whole reserve and finds no incumbent — and + before this fix the fallback's result was merged into ``result`` only when it + carried an objective, so a *sound dual bound* was computed and then thrown away: + the solve reported ``bound=None`` while the fallback had proved ``bound=-27.88``. + That, not the admission itself, was the real cost behind the review's + "spends the budget for nothing". + + Asserts properties rather than values, since the bound tightens with the budget: + a bound exists, it never crosses the 0.0 oracle (``minlplib.solu``), and no + certificate is claimed without an incumbent. Fails before the fix (``bound is + None``), passes after. + """ + result = _ball_mk2_real(30).solve(time_limit=6.0) + assert result.bound is not None, "the fallback's dual bound was discarded again" + assert result.bound <= 0.0 + 1e-6, f"dual bound {result.bound} crossed the 0.0 oracle" + # No incumbent -> no certificate, whatever the bound says. + if result.objective is None: + assert not getattr(result, "gap_certified", False) + else: # a primal fix landed: the certificate invariant must still hold + assert result.bound <= result.objective + 1e-6