pounce_algorithm/ipopt_alg.rs
1//! Main optimization loop — port of
2//! `Algorithm/IpIpoptAlg.{hpp,cpp}`.
3//!
4//! Phase 7 ships the loop scaffold matching `Optimize()` lines
5//! 292-563 in upstream. The body invokes:
6//!
7//! 1. `IterateInitializer::set_initial_iterates`
8//! 2. (loop) `OutputIteration` → `CheckConvergence` →
9//! `UpdateBarrierParameter` → `UpdateHessian` →
10//! `ComputeSearchDirection` → `ComputeAcceptableTrialPoint` →
11//! `AcceptTrialPoint`
12//! 3. `correct_bound_multiplier` (kappa_sigma) per `MAIN_LOOP.md`
13//! §"Bound multiplier reset" lines 1055-1134
14//! 4. exception → `SolverReturn` mapping per the table in
15//! `MAIN_LOOP.md`.
16//!
17//! The NLP handle and search-direction calculator are optional:
18//! when both are present, `iterate()` computes a real Newton step and
19//! drives the line search. Without them, `iterate()` runs the bookkeeping
20//! pieces (mu update, hessian update, conv check, kappa_sigma reset)
21//! and is exercised by structural unit tests. The full path lights up
22//! once `pounce-nlp::OrigIpoptNLP` lands.
23
24use crate::alg_builder::AlgorithmBundle;
25use crate::conv_check::r#trait::ConvergenceStatus;
26use crate::intermediate::{CtxGuard, IntermediateContext};
27use crate::ipopt_cq::IpoptCqHandle;
28use crate::ipopt_data::IpoptDataHandle;
29use crate::ipopt_nlp::IpoptNlp;
30use crate::iter_dump::IterDumper;
31use crate::iterate_dump::emit_record as emit_iterate_record;
32use crate::kkt::pd_search_dir_calc::PdSearchDirCalc;
33use crate::line_search::backtracking::Outcome;
34use crate::restoration::{RestorationOutcome, RestorationPhase};
35use pounce_common::diagnostics::DiagnosticsState;
36use pounce_common::types::{Index, Number};
37use pounce_linalg::Vector;
38use pounce_nlp::alg_types::SolverReturn;
39use pounce_nlp::return_codes::AlgorithmMode;
40use pounce_nlp::tnlp::{IpoptCq as TnlpIpoptCq, IpoptData as TnlpIpoptData, IterStats, TNLP};
41use std::cell::RefCell;
42use std::rc::Rc;
43
44/// Dual-divergence guard (pounce#246): only dual-infeasibility growth in the
45/// *elevated* regime (`inf_du` above this) counts toward the streak, so the
46/// noisy early iterations of a normal solve never build one.
47const DUAL_DIV_COUNT_FLOOR: Number = 1e2;
48/// Dual-divergence guard: the guard fires only once `inf_du` is this large in
49/// absolute terms — well above the transient peaks a converging solve reaches
50/// (e.g. the least-square-init path recovers from ~6e9 on emfl050), and below
51/// the ~1e10 where the KKT factorizations begin to choke, so the diversion to
52/// restoration happens *before* the seconds-long factorizations start.
53const DUAL_DIV_FIRE_TOL: Number = 1e8;
54
55/// gh #534 — how many consecutive outer NLP errors the progress test reads.
56/// Four samples give three ratios: enough that a single lucky step cannot pass
57/// the test, short enough to still be inside the endgame it is meant to
58/// recognise. `eigena2`'s quoted tail is exactly four iterations long
59/// (`1.19e-5 → 2.96e-6 → 7.38e-7 → 1.84e-7`).
60const DECLINE_PROGRESS_SAMPLES: usize = 4;
61/// gh #534 — default `resto_decline_progress_ratio`: every one of those ratios
62/// must be at least this contraction for the decline to be deferred. `eigena2`
63/// quarters (ratio `0.249`) and passes; `eigenb2`'s tail *rises*
64/// (`1.88e-7, 2.69e-7, 2.89e-7, 2.93e-7`) and fails, which is the intended
65/// split — the issue calls `eigenb2` a plausible genuine stall and the guard
66/// plausibly right there.
67const DEFAULT_DECLINE_PROGRESS_RATIO: Number = 0.5;
68/// gh #534 — outer iterations a deferred continuation gets to produce a strict
69/// certificate before it is cut and the floor reported. `eigena2`'s
70/// extrapolation needs three; ten leaves room for a slower but still genuine
71/// endgame while keeping the cost of a lost bet bounded and small.
72const DECLINE_CONTINUATION_BUDGET: Index = 10;
73/// gh #534 — default for `resto_decline_deferrals`. One deferral is enough for
74/// the reported case (the continuation either converges within the budget or it
75/// does not); more entries would mostly re-bet on a point the first bet already
76/// failed to improve.
77const DEFAULT_RESTO_DECLINE_DEFERRALS: usize = 1;
78
79pub struct IpoptAlgorithm {
80 pub data: IpoptDataHandle,
81 pub cq: IpoptCqHandle,
82 pub bundle: AlgorithmBundle,
83 /// Optional NLP handle. Required for any step that evaluates
84 /// problem functions or pulls bound expansion matrices (init,
85 /// search direction, line-search trial-point evaluation). Absent
86 /// in the structural unit tests of Phases 5-6.
87 pub nlp: Option<Rc<RefCell<dyn IpoptNlp>>>,
88 /// Optional TNLP handle — the user-facing problem. When present,
89 /// `iterate()` fires `TNLP::intermediate_callback` once per outer
90 /// iteration so callers can monitor progress or request early
91 /// termination (returning `false` from the callback surfaces as
92 /// `SolverReturn::UserRequestedStop`). Kept separate from `nlp`
93 /// because the algorithm-side NLP is the *compressed* `OrigIpoptNlp`
94 /// view (fixed-variable elimination, c/d split) while the callback
95 /// payload needs to expose the original-coordinate iterate.
96 pub tnlp: Option<Rc<RefCell<dyn TNLP>>>,
97 /// Search-direction calculator (`PdSearchDirCalc`). Lands once a
98 /// concrete `SymLinearSolver` backend (MUMPS / FERAL) is wired
99 /// through `AlgBuilder` in Phase 7's tail.
100 pub search_dir: Option<PdSearchDirCalc>,
101 /// Restoration-phase strategy. Invoked when the line search
102 /// returns [`Outcome::Failed`] (port of upstream
103 /// `IpBacktrackingLineSearch::ActivateLineSearch`'s resto
104 /// fallback). Optional: in its absence, line-search failure maps
105 /// directly to [`SolverReturn::RestorationFailure`] so the main
106 /// loop's exit-code semantics match upstream's "no resto built"
107 /// case.
108 pub restoration: Option<Box<dyn RestorationPhase>>,
109
110 /// `kappa_sigma` for the post-AcceptTrialPoint multiplier reset
111 /// (`IpIpoptAlg.cpp:correct_bound_multiplier`, line 1055-1134).
112 pub kappa_sigma: Number,
113 pub max_iter: Index,
114 /// Initial primal step length offered to the line search at the
115 /// top of each iteration. Mirrors `IpBacktrackingLineSearch`'s
116 /// fraction-to-the-boundary primal step (with τ = `data.curr_tau`).
117 /// In v1.0 the structural value here is 1.0 and the FTB cap is
118 /// applied per-component when the line-search driver computes
119 /// trial slacks; the simplification holds for non-degenerate runs.
120 pub alpha_init: Number,
121 /// Tiny-step relative tolerance — port of upstream
122 /// `IpBacktrackingLineSearch::tiny_step_tol_` (default `10·EPSILON`).
123 /// Step is "tiny" when `max_i |δx_i|/(1+|x_i|) ≤ tiny_step_tol`
124 /// (and same for s, and `c_viol ≤ 1e-4`).
125 pub tiny_step_tol: Number,
126 /// Port of upstream `IpIpoptAlg.cpp` divergence guard: when
127 /// `max_i |x_i|` exceeds this threshold the optimization aborts with
128 /// `SolverReturn::DivergingIterates`. Default `1e20` matches the
129 /// registered `diverging_iterates_tol` option. Catches MESH and
130 /// similar cases where the normal-mode IPM heads off to infinity
131 /// (orig `f` to ±1e33 by iter 90) before line-search failure forces
132 /// a degenerate restoration entry.
133 pub diverging_iterates_tol: Number,
134 /// #248 divergence persistence — consecutive iterations the primal
135 /// iterate has kept *growing* while past `diverging_iterates_tol` on a
136 /// structurally unbounded side. A genuine recession ray sustains this;
137 /// a transient ill-scaling excursion on a bounded-below problem peaks
138 /// and recedes (MINLPLib `jit1`: `|x|` climbs to ~16 then falls back to
139 /// ~2.9 at the finite optimum). Reset to zero whenever the iterate is
140 /// within the threshold or is not growing.
141 divergence_streak: u32,
142 /// Largest `|x|` seen in the current growth run (companion to
143 /// [`Self::divergence_streak`]). Zero when no run is active.
144 divergence_prev_amax: Number,
145 /// #252 objective at the previous over-threshold iterate of the current
146 /// growth run (companion to [`Self::divergence_streak`]). A genuine
147 /// recession ray drives the (minimized) objective toward `−∞`, so a
148 /// diverging iterate only counts toward the streak when the objective is
149 /// *still descending* against this reference. A transient ill-scaling
150 /// excursion past a finite optimum grows `|x|` while the objective
151 /// *worsens* (the linear tail dominates), so it never accumulates the
152 /// streak — this is the fix for the unbounded-box (`ub = +∞`) B&B node
153 /// subproblems of jit1 that #248's growth-only check still mislabelled
154 /// `UNBOUNDED`. `+∞` when no run is active.
155 divergence_prev_f: Number,
156 /// #252 objective *decrease* at the previous step of the current growth
157 /// run (`prev_prev_f − prev_f`), the companion that lets the streak
158 /// require the descent to be *non-decelerating*. A recession ray's
159 /// per-step objective drop keeps up or accelerates as `|x|` grows
160 /// geometrically (`f` is at least linear along the ray); an excursion
161 /// converging to a finite optimum has a per-step drop that shrinks
162 /// toward zero. `NaN` (non-finite) until the run has a first finite
163 /// decrease to compare against, which bootstraps the check.
164 divergence_prev_decrease: Number,
165 /// #285 recession-ray persistence — consecutive iterations for which the
166 /// *checked recession-ray proof* ([`Self::curr_is_recession_ray`]) held
167 /// while the primal iterate kept growing. This is a second, independent
168 /// unboundedness path that catches a genuine recession ray in
169 /// `null(A_eq)` over free variables whose `|x|` grows only *linearly*
170 /// (the regularized zero-Hessian step in an equality null space marches
171 /// out at a bounded rate), so it never crosses `diverging_iterates_tol`
172 /// (`1e20`) within `max_iter` and the geometric-growth
173 /// [`Self::divergence_streak`] never accumulates. Reset to zero whenever
174 /// the proof fails or the iterate stops growing.
175 recession_streak: u32,
176 /// Largest `|x|` seen in the current recession-ray run (companion to
177 /// [`Self::recession_streak`]). Zero when no run is active.
178 recession_prev_amax: Number,
179 /// Companion threshold on the dual step — when both primal and dual
180 /// steps are tiny in two consecutive iterations the algorithm
181 /// declares convergence at the best attainable accuracy. Default
182 /// `1e-2` matches upstream.
183 pub tiny_step_y_tol: Number,
184 /// `dual_diverging_streak` (pounce#246) — number of consecutive
185 /// iterations of *growing* dual infeasibility (in the elevated regime,
186 /// `inf_du > `[`DUAL_DIV_COUNT_FLOOR`]) that must accumulate before the
187 /// dual-divergence guard fires. When the streak reaches the limit and
188 /// `inf_du > `[`DUAL_DIV_FIRE_TOL`], the outer routes to restoration.
189 ///
190 /// **`0` (off) is the default**, set from the option of the same name
191 /// (`application.rs`). It defaulted to `15` when introduced; see the option
192 /// help in `upstream_options.rs` for why that changed, and
193 /// [`Self::honour_best_acceptable_after_dual_guard`] for what protects a
194 /// solve when it is enabled. See the guard itself in [`Self::iterate`].
195 pub dual_diverging_streak: usize,
196 dual_inf_prev: Number,
197 dual_growth_streak: usize,
198 /// Set true when the previous iterate was tagged tiny; on the
199 /// second consecutive tiny step the loop sets `data.tiny_step_flag`
200 /// so the mu update can attempt to terminate. Mirrors
201 /// `IpBacktrackingLineSearch::tiny_step_last_iteration_`.
202 pub tiny_step_last_iteration: bool,
203 /// Cycle-detection state for [`Self::invoke_restoration`]: the
204 /// outer `(x, s)` snapshot from the previous restoration entry,
205 /// cleared on any iteration that exits via a normal line-search
206 /// accept. When restoration is invoked twice in a row and the
207 /// outer iterate has not moved between entries (relative
208 /// 2-norm < 1e-10 on both `x` and `s`), the inner resto-IPM is
209 /// returning Recovered points indistinguishable from `curr` — a
210 /// cycle. Surfaces as `ErrorInStepComputation`. Mirrors the
211 /// *intent* of upstream `IpBacktrackingLineSearch.cpp:580-600`'s
212 /// almost-feasible-resto guard while staying robust against the
213 /// `inf_pr` micro-drift seen on ACOPR14 (delta ~3e-12 per entry,
214 /// inf_du essentially constant) where a scalar-`inf_pr` heuristic
215 /// fails. Productive single-restoration sequences (BT8, HIMMELBJ,
216 /// LINSPANH, LSNNODOC, ODFITS, OET3) clear the snapshot via
217 /// `Outcome::Accepted` between entries and are unaffected.
218 last_resto_entry_x: Option<Box<dyn Vector>>,
219 last_resto_entry_s: Option<Box<dyn Vector>>,
220 /// Snapshot of the *recovery* iterate from the previous
221 /// restoration. Compared against the next entry's `(x, s)` to
222 /// detect "outer made no progress between consecutive resto
223 /// invocations". When this distance is below threshold for
224 /// several consecutive entries, terminate — catching
225 /// slow-non-convergence cycles (ACOPR14, TRO3X3, ACOPR30) where
226 /// resto's *inner* moves substantively each call but the *outer*
227 /// makes no progress between calls. Cleared on any LS-accepted
228 /// step.
229 last_resto_recovery_x: Option<Box<dyn Vector>>,
230 last_resto_recovery_s: Option<Box<dyn Vector>>,
231 /// Count of consecutive restoration entries on which the outer
232 /// step (recovery → next-entry) was below the iterate-distance
233 /// threshold. Cleared on any LS-accepted step. Limit chosen to
234 /// let MAKELA3, HAIFAM, HALDMADS, ROBOT, TENBARS2 — which need
235 /// 2-3 consecutive resto entries to recover — pass through.
236 resto_no_outer_progress_count: usize,
237 /// `resto_decline_deferrals` (gh #534) — how many times the
238 /// acceptable-point restoration decline in [`Self::invoke_restoration`] may
239 /// be *deferred* on a solve whose NLP error is still contracting. `0`
240 /// restores the pre-#534 behaviour (decline immediately, always).
241 ///
242 /// See [`Self::may_defer_acceptable_decline`] for the progress test and
243 /// [`Self::honour_decline_floor`] for what makes a spent deferral harmless.
244 pub resto_decline_deferrals: usize,
245 /// `resto_decline_progress_ratio` (gh #534) — the contraction each of the
246 /// last [`DECLINE_PROGRESS_SAMPLES`]` - 1` iterations must have achieved for
247 /// the decline to be deferred. Default
248 /// [`DEFAULT_DECLINE_PROGRESS_RATIO`]. A value of `1` admits any
249 /// non-increasing window and a large one drops the progress requirement
250 /// altogether, which is the "patch the guard and see" experiment the issue
251 /// asks for, available without patching.
252 pub resto_decline_progress_ratio: Number,
253 /// The most recent outer-iteration NLP errors, oldest first, with
254 /// [`Self::nlp_err_recent_len`] entries live. Feeds the gh #534 progress
255 /// test and nothing else.
256 nlp_err_recent: [Number; DECLINE_PROGRESS_SAMPLES],
257 nlp_err_recent_len: usize,
258 /// gh #534 — deferrals of the acceptable-point decline spent so far.
259 decline_deferrals_used: usize,
260 /// gh #534 — the iterate the guard would have returned had it not been
261 /// deferred. Captured at the *first* deferral only, because that point is
262 /// precisely the answer the pre-#534 build reports; it is the floor the
263 /// continuation is never allowed to fall below.
264 decline_floor: Option<VetoSnapshot>,
265 /// gh #534 — outer iteration by which the deferred continuation must have
266 /// produced a strict certificate. Past it the continuation is cut and the
267 /// floor is reported, so the bet costs a bounded number of iterations.
268 decline_deadline_iter: Option<Index>,
269 /// Count of consecutive restoration entries on which the outer
270 /// constraint violation at entry was already below `tol` (the
271 /// outer optimality tolerance). Matches the *intent* of upstream
272 /// `IpBacktrackingLineSearch.cpp:580-600`'s almost-feasible-resto
273 /// guard while using a looser cv threshold (`tol` vs `1e-2·tol`)
274 /// — catches DECONVBNE's resto-thrash where each cycle re-enters
275 /// at cv ≈ 3e-10 < tol with bound multipliers reset to 1, the
276 /// outer's σ-blowup explodes inf_du to 1.9e7, alpha-min triggers
277 /// resto re-entry, and the (inf_pr, inf_du) post-recovery state
278 /// is essentially identical across cycles but `x` drifts enough
279 /// that [`Self::last_resto_recovery_x`]-based detection misses.
280 /// Cumulative (never cleared on LS-accept), since DECONVBNE's
281 /// cycle interleaves R-recoveries with sub-tol accepts that
282 /// accomplish no real outer progress. Fires after 3 near-feasible
283 /// entries — surfaces as `StopAtAcceptablePoint` since the
284 /// recovered point already satisfies constraint feasibility
285 /// within `tol`.
286 resto_near_feasible_count: usize,
287 /// Snapshot of the most recent iterate that the convergence check
288 /// flagged "acceptable" (NLP error ≤ `acceptable_tol`). Mirrors
289 /// upstream `IpBacktrackingLineSearch::acceptable_iterate_`
290 /// (`IpBacktrackingLineSearch.cpp:1286-1310`). Used by
291 /// [`Self::restore_acceptable_point`] to roll back when restoration
292 /// fails — if such an iterate exists, the algorithm exits with
293 /// `SolverReturn::StopAtAcceptablePoint` rather than
294 /// `RestorationFailure`. Cleared/refreshed on every iteration that
295 /// satisfies the acceptable predicate.
296 acceptable_iterate: Option<crate::iterates_vector::IteratesVector>,
297 /// The first iterate whose *strict* certificate the masked-scale veto
298 /// refused (gh #200), kept so the refusal can be undone verbatim if the
299 /// continued run does not do better. Deliberately not the acceptable
300 /// snapshot: that one is overwritten unconditionally and drifts.
301 vetoed: Option<VetoSnapshot>,
302 /// The iterate at which a refused *acceptable-level* termination would have
303 /// fired. Held separately from `vetoed` because it restores under a weaker
304 /// status, and claiming `Success` for it would over-report.
305 vetoed_acceptable: Option<VetoSnapshot>,
306 /// Whether a strict refusal has already been *seen*, independent of whether
307 /// a snapshot was successfully captured for it.
308 ///
309 /// This is the first-only latch, held apart from `vetoed` on purpose.
310 /// Testing `vetoed.is_none()` instead would let a refusal whose capture
311 /// failed be "completed" at a later iterate — the veto flag on the
312 /// convergence check is sticky, so it still reads true next pass, and the
313 /// fallback would then restore a point that never passed the strict test.
314 /// With the latch, a failed capture stays failed and the fallback declines.
315 ///
316 /// Declining is *not* the baseline outcome — the baseline stopped and
317 /// reported a certificate at the uncaptured iterate, and declining fails to
318 /// reproduce it. It is the least-bad handling of an unidentifiable baseline,
319 /// not a faithful one.
320 vetoed_seen: bool,
321 /// Same latch for the acceptable-level refusal.
322 vetoed_acceptable_seen: bool,
323 /// Whether the dual-divergence guard (pounce#246) actually fired this
324 /// solve. Gates the *use* of [`Self::best_acceptable`], so solves the guard
325 /// never touches behave identically — see
326 /// [`Self::honour_best_acceptable_after_dual_guard`].
327 dual_guard_fired: bool,
328 /// Best (lowest scaled objective) acceptable-quality iterate seen anywhere
329 /// in this solve. Recorded unconditionally — including *before* any
330 /// diversion, which is the point: the guard returns to the driver before
331 /// the recording site on the iteration it fires, so gating the recording on
332 /// `dual_guard_fired` would miss everything up to and including the
333 /// diversion. Only read when `dual_guard_fired`.
334 best_acceptable: Option<VetoSnapshot>,
335 /// `kkt_fidelity_tol` (pounce#173), needed here — not just at termination —
336 /// because the fallback's tiebreak has to predict the post-solve status
337 /// gate. See [`Self::honour_refused_certificate`]. Zero (the default)
338 /// disables the gate, and with it every tiebreak effect it has.
339 pub kkt_fidelity_tol: Number,
340 acceptable_iter_number: Index,
341 /// Shared per-solve diagnostics state. `None` unless the CLI
342 /// requested `--dump <cat>:<spec>`. When set, the outer loop
343 /// advances the state's iter counter and the augmented-system
344 /// solver consults it to gate KKT dumps.
345 diagnostics: Option<Rc<DiagnosticsState>>,
346 /// Optional interactive debugger. Shared (`Rc<RefCell<…>>`) so the
347 /// same debugger instance also drives the restoration inner IPM —
348 /// one debugger sees both levels. Fired at every
349 /// [`crate::debug::Checkpoint`]. See `crate::debug`.
350 debug: Option<Rc<RefCell<dyn crate::debug::DebugHook>>>,
351
352 // ---- Restoration-phase audit counters (pounce#12). ----
353 //
354 // Drained into `SolveStatistics` by `IpoptApplication::optimize_constrained`
355 // after the solve completes. Counts are cumulative across the run.
356 /// Number of `invoke_restoration` entries.
357 pub resto_calls: Index,
358 /// Sum of inner-IPM iter counts across every restoration call.
359 pub resto_inner_iters: Index,
360 /// Number of outer iters that ran in restoration mode (R-line
361 /// equivalents in `print_level=5` output).
362 pub resto_outer_iters: Index,
363 /// Cumulative wall-clock seconds spent inside `perform_restoration`.
364 pub resto_wall_secs: Number,
365
366 // ---- Per-iteration history capture (pounce#8, pounce#71). ----
367 //
368 // The per-iteration trajectory is no longer accumulated on the
369 // algorithm: `iterate()` emits a structured `pounce::iteration`
370 // event each step, and `pounce_observability::IterCollectorLayer`
371 // rebuilds the `IterRecord`s into the active `IterCaptureGuard`
372 // that `IpoptApplication` installs around the solve.
373 /// When `false`, the per-iteration table that `iterate()` writes
374 /// straight to stdout is suppressed. Wired from
375 /// `IpoptApplication`'s `print_level` option: level 0 turns this
376 /// off (matches upstream's "no console output" contract). Default
377 /// `true` so CLI / direct-driver users keep the familiar trace.
378 pub print_iter_output: bool,
379}
380
381impl IpoptAlgorithm {
382 pub fn new(data: IpoptDataHandle, cq: IpoptCqHandle, mut bundle: AlgorithmBundle) -> Self {
383 // The builder may pre-populate `bundle.search_dir` when given a
384 // `LinearBackendFactory`; lift it onto the algorithm so the
385 // iterate body can call into it directly.
386 let search_dir = bundle.search_dir.take();
387 Self {
388 data,
389 cq,
390 bundle,
391 nlp: None,
392 tnlp: None,
393 search_dir,
394 restoration: None,
395 kappa_sigma: 1e10,
396 max_iter: 3000,
397 alpha_init: 1.0,
398 tiny_step_tol: 10.0 * Number::EPSILON,
399 diverging_iterates_tol: 1e20,
400 divergence_streak: 0,
401 divergence_prev_amax: 0.0,
402 divergence_prev_f: Number::INFINITY,
403 divergence_prev_decrease: Number::NAN,
404 recession_streak: 0,
405 recession_prev_amax: 0.0,
406 tiny_step_y_tol: 1e-2,
407 dual_diverging_streak: 15,
408 dual_inf_prev: 0.0,
409 dual_growth_streak: 0,
410 tiny_step_last_iteration: false,
411 last_resto_entry_x: None,
412 last_resto_entry_s: None,
413 last_resto_recovery_x: None,
414 last_resto_recovery_s: None,
415 resto_no_outer_progress_count: 0,
416 resto_decline_deferrals: DEFAULT_RESTO_DECLINE_DEFERRALS,
417 resto_decline_progress_ratio: DEFAULT_DECLINE_PROGRESS_RATIO,
418 nlp_err_recent: [Number::NAN; DECLINE_PROGRESS_SAMPLES],
419 nlp_err_recent_len: 0,
420 decline_deferrals_used: 0,
421 decline_floor: None,
422 decline_deadline_iter: None,
423 resto_near_feasible_count: 0,
424 acceptable_iterate: None,
425 vetoed: None,
426 vetoed_acceptable: None,
427 dual_guard_fired: false,
428 best_acceptable: None,
429 vetoed_seen: false,
430 vetoed_acceptable_seen: false,
431 kkt_fidelity_tol: 0.0,
432 acceptable_iter_number: 0,
433 diagnostics: None,
434 debug: None,
435 resto_calls: 0,
436 resto_inner_iters: 0,
437 resto_outer_iters: 0,
438 resto_wall_secs: 0.0,
439 print_iter_output: true,
440 }
441 }
442
443 /// Stash the current iterate as the "last acceptable" backup —
444 /// port of `IpBacktrackingLineSearch::StoreAcceptablePoint`
445 /// (`IpBacktrackingLineSearch.cpp:1286-1293`).
446 fn store_acceptable_point(&mut self) {
447 let d = self.data.borrow();
448 if let Some(curr) = d.curr.as_ref() {
449 self.acceptable_iterate = Some(curr.clone());
450 self.acceptable_iter_number = d.iter_count;
451 }
452 }
453
454 /// Record this outer iteration's NLP error for the gh #534 progress test.
455 ///
456 /// One push per `iterate()` call, so the samples are consecutive outer
457 /// iterations by construction. Deliberately *not* cleared when restoration
458 /// recovers: a recovery that helped shows up as continued contraction and a
459 /// recovery that hurt shows up as a jump, and the ratio test reads both
460 /// correctly without needing to know which happened.
461 fn note_nlp_err(&mut self, nlp_err: Number) {
462 push_sample(
463 &mut self.nlp_err_recent,
464 &mut self.nlp_err_recent_len,
465 nlp_err,
466 );
467 }
468
469 /// Whether the last [`DECLINE_PROGRESS_SAMPLES`] outer iterations each cut
470 /// the NLP error by at least `resto_decline_progress_ratio` (gh #534).
471 ///
472 /// The question the restoration-decline guard never asked: *is this solve
473 /// still converging?* A full window is required, so the test cannot pass on
474 /// a short history — the early iterations of every solve included.
475 ///
476 /// The test itself lives in the pure [`window_is_contracting`], for the
477 /// reason [`ranks_better_within_band`] does: what it must and must not fire
478 /// on is stated in the issue as two recorded traces, and those are provable
479 /// by deterministic unit test rather than inferable from a solve.
480 fn nlp_err_contracting(&self) -> bool {
481 if self.nlp_err_recent_len < DECLINE_PROGRESS_SAMPLES {
482 return false;
483 }
484 window_is_contracting(&self.nlp_err_recent, self.resto_decline_progress_ratio)
485 }
486
487 /// The live progress window, oldest first, for the gh #534 trace lines.
488 fn nlp_err_window_str(&self) -> String {
489 let live = &self.nlp_err_recent[..self.nlp_err_recent_len];
490 let parts: Vec<String> = live.iter().map(|e| format!("{e:.3e}")).collect();
491 format!("[{}]", parts.join(" -> "))
492 }
493
494 /// Roll the iterate back to the last acceptable snapshot — port of
495 /// `IpBacktrackingLineSearch::RestoreAcceptablePoint`
496 /// (`IpBacktrackingLineSearch.cpp:1295-1310`). Returns `true` if a
497 /// snapshot was available and applied; `false` otherwise (caller
498 /// then surfaces the original failure status).
499 fn restore_acceptable_point(&mut self) -> bool {
500 let Some(prev) = self.acceptable_iterate.clone() else {
501 return false;
502 };
503 let mut d = self.data.borrow_mut();
504 d.set_trial(prev);
505 // `accept_trial_point` promotes `trial → curr`, mirroring the
506 // upstream sequence `set_trial(...); AcceptTrialPoint();`.
507 d.accept_trial_point();
508 true
509 }
510
511 /// Whether a diverging primal iterate is consistent with the feasible
512 /// region actually being *unbounded* (issue #248).
513 ///
514 /// `DivergingIterates` is Ipopt's unboundedness verdict, but a large
515 /// `|x_i|` only proves unboundedness if variable `i` is free to escape
516 /// to infinity in the direction it is heading — i.e. it has no finite
517 /// bound on that side. This lifts a vector of ones from the compressed
518 /// lower/upper bound spaces through the `Px_L` / `Px_U` expansion
519 /// matrices to obtain full-length indicators of which variables carry a
520 /// finite bound, then returns `true` only when some component whose
521 /// magnitude exceeds `diverging_iterates_tol` is heading toward a side
522 /// with no finite bound.
523 ///
524 /// When every large component is pinned by a finite bound — in
525 /// particular when all variables are boxed, so the feasible region is a
526 /// bounded box and unboundedness is structurally impossible — this
527 /// returns `false`, and the caller reports the best iterate via the
528 /// normal convergence / restoration path instead of a spurious
529 /// `Unbounded`.
530 /// #248: consecutive growing, over-threshold iterations required before
531 /// a structurally-free divergence is reported as `DivergingIterates`.
532 /// `jit1`'s transient excursion lasts ~2 growing steps and then
533 /// recedes, so a small persistence requirement clears it without
534 /// materially delaying a genuine ray.
535 const DIVERGENCE_PERSIST_ITERS: u32 = 4;
536 /// #248: an iterate counts as "still growing" toward divergence when it
537 /// grows at least this factor over the previous over-threshold iterate.
538 /// A recession ray in an interior-point method grows geometrically; an
539 /// iterate settling onto a finite optimum above the threshold does not.
540 const DIVERGENCE_GROWTH_FACTOR: Number = 2.0;
541 /// #252: the objective descent must *keep up* — each step's drop must be
542 /// at least this fraction of the previous step's drop for the iterate to
543 /// count toward the divergence streak. A recession ray descends `f` to
544 /// `−∞` with per-step drops that grow (ratio ≥ 1) as `|x|` grows
545 /// geometrically; an excursion converging to a finite optimum decelerates
546 /// (ratio → 0). The slack below 1 tolerates ordinary interior-point noise
547 /// on a genuine ray without admitting a decelerating excursion — jit1's
548 /// node subproblems shrink the drop by 3–15× per step, far past this bar.
549 const DIVERGENCE_DESCENT_KEEPUP: Number = 0.9;
550 /// #248: absolute runaway backstop. An iterate this large is reported
551 /// unbounded regardless of persistence. It sits at or below the default
552 /// `diverging_iterates_tol = 1e20`, so the default behaviour (fire the
553 /// instant `|x|` crosses the threshold) is preserved, while a low
554 /// user threshold no longer fires on the way to a finite optimum.
555 const DIVERGENCE_ABS_RUNAWAY: Number = 1e18;
556
557 /// #285: magnitude floor for the checked recession-ray unboundedness path.
558 /// Below this the (slightly more expensive) recession proof is not even
559 /// attempted, so it is inert on every normal, well-scaled solve. Above it,
560 /// unboundedness is only ever concluded through the full checked proof in
561 /// [`Self::curr_is_recession_ray`] — a genuinely *feasible* iterate of this
562 /// magnitude already witnesses an unbounded feasible region, and the proof
563 /// additionally certifies the escape direction. Sits far below the
564 /// `diverging_iterates_tol` (`1e20`) magnitude guard so a linearly-growing
565 /// ray (which never reaches `1e20` within `max_iter`) is still caught.
566 const RECESSION_MIN_NORM: Number = 1e10;
567 /// #285: consecutive growing, proof-passing iterations required before the
568 /// recession-ray path reports `DivergingIterates`. A bounded feasible
569 /// region cannot supply a *growing* sequence of feasible over-floor
570 /// iterates, so persistence is defense-in-depth against a lone numerical
571 /// fluke rather than a soundness requirement.
572 const RECESSION_PERSIST_ITERS: u32 = 4;
573 /// #285: relative feasibility bar for the recession proof. The current
574 /// iterate counts as feasible (hence a witness that the feasible region
575 /// reaches its magnitude) when its unscaled max-norm primal infeasibility
576 /// is at most this fraction of `|x|_∞`. The check is *relative* on purpose:
577 /// evaluating `A_eq x − b` at `|x| ~ 1e17` carries floating-point roundoff
578 /// that scales with `|x|`, while a genuinely infeasible excursion (e.g.
579 /// mid-restoration) has a residual comparable to `|x|` itself.
580 const RECESSION_FEAS_REL: Number = 1e-6;
581 /// #285: relative bar for "the escape direction lies in `null(A_eq)`".
582 /// `‖J_c x‖_∞ ≤ this · |x|_∞` certifies that moving along `d ≈ x` preserves
583 /// the (linearized) equality constraints — `A_eq d ≈ 0`.
584 const RECESSION_DIR_TOL: Number = 1e-6;
585 /// #285: relative descent bar. The objective must strictly decrease along
586 /// the escape direction with a real margin — `∇f·x ≤ −this · ‖∇f‖ ‖x‖` —
587 /// so a variable drifting orthogonally to the objective (`∇f·x ≈ 0`) can
588 /// never be mistaken for a recession ray driving `f → −∞`.
589 const RECESSION_DESC_REL: Number = 1e-6;
590
591 /// Update the divergence-persistence state for the current iterate and
592 /// return whether `DivergingIterates` should be reported now (issues
593 /// #248 / #252). `amax` is `max_i |x_i|`; `structural_free` is the result
594 /// of [`Self::divergence_is_true_unboundedness`] (already gated on
595 /// `amax > diverging_iterates_tol`); `f` is the (minimized, internally
596 /// scaled) objective at the current iterate, supplied only while
597 /// `structural_free` holds.
598 ///
599 /// A large `|x|` is reported as unbounded only when it is heading to an
600 /// unbounded side (`structural_free`) *and* the divergence looks like a
601 /// genuine recession ray: the iterate keeps *growing* while the objective
602 /// keeps *descending toward `−∞` without decelerating* — the per-step drop
603 /// holds up as `|x|` grows geometrically — for
604 /// [`Self::DIVERGENCE_PERSIST_ITERS`] consecutive iterations (or it has
605 /// blown past the absolute runaway backstop). Two failure modes are thereby
606 /// left to the normal convergence machinery instead of being mislabelled
607 /// `UNBOUNDED`:
608 ///
609 /// * #248 — a transient ill-scaling excursion that peaks in `|x|` and
610 /// recedes never sustains the growth streak.
611 /// * #252 — an excursion that *keeps* growing in `|x|` toward an unbounded
612 /// box side (a jit1 B&B node subproblem with `ub = +∞`), lowering `f` as
613 /// it goes, but with a per-step objective drop that *decelerates* toward
614 /// zero: it is settling onto a finite optimum, not riding a recession
615 /// ray. The descent must keep up (not merely exist), so this no longer
616 /// accumulates the streak.
617 fn update_divergence_verdict(
618 &mut self,
619 amax: Option<Number>,
620 structural_free: bool,
621 f: Option<Number>,
622 ) -> bool {
623 let over = matches!(amax, Some(a) if a > self.diverging_iterates_tol) && structural_free;
624 if !over {
625 self.divergence_streak = 0;
626 self.divergence_prev_amax = 0.0;
627 self.divergence_prev_f = Number::INFINITY;
628 self.divergence_prev_decrease = Number::NAN;
629 return false;
630 }
631 let a = amax.expect("over implies amax is Some");
632 // A recession ray in an interior-point method grows the iterate
633 // geometrically *and* drives the objective down without bound, with a
634 // per-step drop that keeps up as `|x|` grows. A finite-optimum
635 // excursion may grow `|x|` and even lower `f` for a few steps, but its
636 // per-step objective drop decelerates toward zero as it settles onto
637 // the finite floor. Require all three — growth, descent, and
638 // non-decelerating descent — before a step counts toward the streak.
639 let growing = a >= self.divergence_prev_amax * Self::DIVERGENCE_GROWTH_FACTOR;
640 // `f` is `None` only when `structural_free` is false, already handled
641 // by the `!over` branch; treat a missing value as non-descending so a
642 // run can never accumulate without objective evidence.
643 let fv = f.unwrap_or(Number::INFINITY);
644 let decrease = self.divergence_prev_f - fv;
645 let descending = decrease > 0.0;
646 // Non-decelerating: the drop must be at least a fixed fraction of the
647 // previous step's drop. Bootstrapped `true` until a first finite
648 // decrease has been recorded (`divergence_prev_decrease` non-finite),
649 // so the run's opening steps are admitted on growth + descent alone.
650 let keeping_up = !self.divergence_prev_decrease.is_finite()
651 || decrease >= self.divergence_prev_decrease * Self::DIVERGENCE_DESCENT_KEEPUP;
652 if growing && descending && keeping_up {
653 self.divergence_streak += 1;
654 } else {
655 // Over the threshold on an unbounded side, but the divergence is
656 // not sustaining a recession ray's growth-and-accelerating-descent
657 // profile — the hallmark of a scaling excursion toward a finite
658 // optimum. Drop the streak; a genuine ray re-accumulates it on its
659 // next qualifying step (or trips the absolute runaway backstop).
660 self.divergence_streak = 0;
661 }
662 self.divergence_prev_amax = a;
663 self.divergence_prev_f = fv;
664 // Record the baseline for the next step's keep-up comparison only from
665 // a finite, real decrease; skip the `+∞` opening step and reset the
666 // baseline whenever the objective stops descending.
667 self.divergence_prev_decrease = if decrease.is_finite() && descending {
668 decrease
669 } else {
670 Number::NAN
671 };
672 a >= Self::DIVERGENCE_ABS_RUNAWAY
673 || self.divergence_streak >= Self::DIVERGENCE_PERSIST_ITERS
674 }
675
676 fn divergence_is_true_unboundedness(&self, x: &dyn Vector) -> bool {
677 self.free_to_escape_over(x, self.diverging_iterates_tol)
678 }
679
680 /// Shared core of the free-variable structural check: returns `true` when
681 /// some component of `x` with magnitude exceeding `thresh` is heading
682 /// toward a side (positive → upper, negative → lower) that carries *no*
683 /// finite bound, so it is free to escape to infinity. Parameterized on the
684 /// magnitude threshold so both the `diverging_iterates_tol` (`1e20`)
685 /// divergence guard and the lower `RECESSION_MIN_NORM` recession-ray path
686 /// (#285) share one implementation.
687 fn free_to_escape_over(&self, x: &dyn Vector, thresh: Number) -> bool {
688 use pounce_linalg::DenseVector;
689
690 let cq = self.cq.borrow();
691 let nlp = cq.nlp().borrow();
692
693 // Full-length 0/1 indicators of finite lower / upper bounds,
694 // built by scattering ones through the bound expansion matrices.
695 let mut ones_l = nlp.x_l().make_new();
696 ones_l.set(1.0);
697 let mut has_lb = x.make_new();
698 nlp.px_l().mult_vector(1.0, &*ones_l, 0.0, &mut *has_lb);
699
700 let mut ones_u = nlp.x_u().make_new();
701 ones_u.set(1.0);
702 let mut has_ub = x.make_new();
703 nlp.px_u().mult_vector(1.0, &*ones_u, 0.0, &mut *has_ub);
704
705 let downcast = |v: &dyn Vector| -> Option<Vec<Number>> {
706 v.as_any()
707 .downcast_ref::<DenseVector>()
708 .map(|d| d.expanded_values())
709 };
710
711 // POUNCE is dense-only; if a backing is unexpectedly non-dense we
712 // cannot prove the divergence is spurious, so fall back to the
713 // original (magnitude-only) verdict to avoid changing behaviour.
714 let (Some(xv), Some(lb), Some(ub)) = (downcast(x), downcast(&*has_lb), downcast(&*has_ub))
715 else {
716 return true;
717 };
718
719 for i in 0..xv.len() {
720 if xv[i].abs() > thresh {
721 let free_to_diverge = if xv[i] > 0.0 {
722 ub[i] == 0.0
723 } else {
724 lb[i] == 0.0
725 };
726 if free_to_diverge {
727 return true;
728 }
729 }
730 }
731 false
732 }
733
734 /// #285: checked recession-ray unboundedness proof at the current iterate.
735 ///
736 /// Returns `true` only when the current iterate `x` (with `|x|_∞ = amax`,
737 /// already known `> RECESSION_MIN_NORM` by the caller) *proves* the
738 /// problem is unbounded below via a genuine recession ray — the same
739 /// standard the LP/symmetric path holds itself to, not a magnitude
740 /// heuristic. All of the following must hold:
741 ///
742 /// 1. **Feasible witness.** The iterate's unscaled primal infeasibility is
743 /// at most `RECESSION_FEAS_REL · amax`. A genuinely feasible iterate of
744 /// norm `≥ 1e10` witnesses that the feasible region reaches that far —
745 /// a *bounded* region cannot contain it. (Relative bar: the residual of
746 /// `A_eq x − b` carries roundoff that scales with `|x|`.)
747 /// 2. **Free to escape.** Some over-floor component heads toward a side
748 /// with no finite variable bound ([`Self::free_to_escape_over`] at
749 /// `RECESSION_MIN_NORM`).
750 /// 3. **Direction in `null(A_eq)`.** `‖J_c x‖_∞ ≤ RECESSION_DIR_TOL · amax`
751 /// — moving along `d ≈ x` preserves the equality constraints.
752 /// 4. **Inequalities not blocking.** No finitely-bounded inequality row is
753 /// driven toward its bound along `d ≈ x`
754 /// ([`Self::recession_blocked_by_inequality`]).
755 /// 5. **Objective descending.** `∇f·x ≤ −RECESSION_DESC_REL · ‖∇f‖ ‖x‖` —
756 /// the objective strictly decreases along the escape direction with a
757 /// real (non-orthogonal) margin, so `f → −∞` along the ray.
758 ///
759 /// On a *bounded* problem at least one of (1)/(2)/(3)/(4)/(5) fails, so
760 /// this can never manufacture a spurious `DivergingIterates`.
761 fn curr_is_recession_ray(&self, x: &dyn Vector, amax: Number) -> bool {
762 // (1) Feasible witness (relative bar).
763 let primal_inf = self.cq.borrow().curr_unscaled_primal_infeasibility_max();
764 if !(primal_inf.is_finite() && primal_inf <= Self::RECESSION_FEAS_REL * amax) {
765 return false;
766 }
767 // (2) Some over-floor component free to escape to infinity.
768 if !self.free_to_escape_over(x, Self::RECESSION_MIN_NORM) {
769 return false;
770 }
771 // (3) Escape direction lies in the equality null space. A non-finite
772 // (NaN) residual is treated as failing, so the direction is only
773 // accepted on a genuinely small, finite `‖J_c x‖∞`.
774 let jc_x_amax = self.cq.borrow().curr_jac_c_times_vec(x).amax();
775 if !jc_x_amax.is_finite() || jc_x_amax > Self::RECESSION_DIR_TOL * amax {
776 return false;
777 }
778 // (4) No finitely-bounded inequality blocks the direction.
779 if self.recession_blocked_by_inequality(x, amax) {
780 return false;
781 }
782 // (5) Objective strictly descending along the escape direction.
783 let (dot, gnorm, xnorm) = {
784 let cq = self.cq.borrow();
785 let g = cq.curr_grad_f();
786 (g.dot(x), g.nrm2(), x.nrm2())
787 };
788 if !(dot < 0.0 && dot <= -Self::RECESSION_DESC_REL * gnorm * xnorm) {
789 return false;
790 }
791 true
792 }
793
794 /// #285: does any *finitely-bounded* inequality constraint block motion
795 /// along the escape direction `d ≈ x`? For each inequality row the
796 /// constraint value `d(x)` changes at rate `(J_d x)_j` per unit of the
797 /// direction; if that row has a finite upper bound and the rate is
798 /// positive (or a finite lower bound and the rate is negative) beyond a
799 /// relative tolerance, moving out along `d` would eventually violate it,
800 /// so it is not a feasible recession direction. Bounds are detected via
801 /// the `Pd_L / Pd_U` expansion matrices exactly as the variable-bound
802 /// check uses `Px_L / Px_U`.
803 fn recession_blocked_by_inequality(&self, x: &dyn Vector, amax: Number) -> bool {
804 use pounce_linalg::DenseVector;
805
806 let cq = self.cq.borrow();
807 // Rate of change of each inequality value along d ≈ x (length m_ineq).
808 // Compute first so the internal `nlp.borrow_mut()` is released before
809 // the immutable borrow below.
810 let jd_x = cq.curr_jac_d_times_vec(x);
811 let (has_dlb, has_dub) = {
812 let nlp = cq.nlp().borrow();
813 let mut ones_dl = nlp.d_l().make_new();
814 ones_dl.set(1.0);
815 let mut has_dlb = jd_x.make_new();
816 nlp.pd_l().mult_vector(1.0, &*ones_dl, 0.0, &mut *has_dlb);
817
818 let mut ones_du = nlp.d_u().make_new();
819 ones_du.set(1.0);
820 let mut has_dub = jd_x.make_new();
821 nlp.pd_u().mult_vector(1.0, &*ones_du, 0.0, &mut *has_dub);
822 // Order matters: bind `(has_dlb, has_dub)` in that exact order so
823 // the finite-lower / finite-upper indicators are not transposed.
824 // #314: this pair was returned swapped, inverting the bound
825 // semantics below — a ray *increasing* a lower-bounded row (moving
826 // deeper into the feasible set, slack growing) was wrongly treated
827 // as blocked, so a genuine inequality-slack recession ray was never
828 // proven unbounded.
829 (has_dlb, has_dub)
830 };
831
832 let downcast = |v: &dyn Vector| -> Option<Vec<Number>> {
833 v.as_any()
834 .downcast_ref::<DenseVector>()
835 .map(|d| d.expanded_values())
836 };
837 // Dense-only fallback: if we cannot inspect the rows, conservatively
838 // treat the direction as blocked (no spurious unbounded verdict).
839 let (Some(jd), Some(dlb), Some(dub)) =
840 (downcast(&*jd_x), downcast(&*has_dlb), downcast(&*has_dub))
841 else {
842 return true;
843 };
844 let tol = Self::RECESSION_DIR_TOL * amax;
845 for j in 0..jd.len() {
846 // Increasing a row that has a finite upper bound, or decreasing a
847 // row that has a finite lower bound, would leave the feasible set.
848 if (jd[j] > tol && dub[j] != 0.0) || (jd[j] < -tol && dlb[j] != 0.0) {
849 return true;
850 }
851 }
852 false
853 }
854
855 /// #285: update the recession-ray persistence state and return whether
856 /// `DivergingIterates` should be reported now. `amax` is `|x|_∞`;
857 /// `is_ray` is the result of [`Self::curr_is_recession_ray`]. The verdict
858 /// fires once the checked proof has held for
859 /// [`Self::RECESSION_PERSIST_ITERS`] consecutive *growing* iterations — a
860 /// bounded region cannot supply a growing sequence of feasible over-floor
861 /// iterates, so this is impossible to satisfy on a bounded problem.
862 fn update_recession_verdict(&mut self, amax: Number, is_ray: bool) -> bool {
863 if !is_ray {
864 self.recession_streak = 0;
865 self.recession_prev_amax = 0.0;
866 return false;
867 }
868 if amax > self.recession_prev_amax {
869 self.recession_streak += 1;
870 } else {
871 // Proof holds but the iterate is not growing (a stalled or rejected
872 // step). Restart the run at the current witness rather than firing
873 // on a plateau; a genuine ray resumes growing next step.
874 self.recession_streak = 1;
875 }
876 self.recession_prev_amax = amax;
877 self.recession_streak >= Self::RECESSION_PERSIST_ITERS
878 }
879
880 /// Honour a certificate the masked-scale veto refused, when the run that
881 /// was allowed to continue did not end in one of its own (gh #200).
882 ///
883 /// The veto's bargain is "never worse off": it refuses a point that had
884 /// *already passed the strict test*, betting that continuing reaches a
885 /// better one. This is the losing side of that bet — so hand back exactly
886 /// what would have been returned without the veto, point and status both.
887 ///
888 /// Two details make that guarantee real rather than approximate:
889 ///
890 /// - It runs on **every** non-success exit, applied once where the driver
891 /// loop's result is finalized. Wiring individual termination sites was
892 /// tried and is not safe: there are sixteen, and the ones easiest to
893 /// overlook are the ones most likely to fire here — the veto's extra
894 /// iterations are exactly what pushes a run past `max_cpu_time`.
895 /// - It restores the **refused iterate itself** (`vetoed`), not the last
896 /// acceptable snapshot. `store_acceptable_point` overwrites
897 /// unconditionally, so after the veto the stored point drifts to whatever
898 /// the continued run last touched — which may be worse than the point
899 /// that was refused.
900 ///
901 /// "Better" is **status-dominant lexicographic**: the reported status first,
902 /// and the objective only to break a tie *within equal status*. Both halves
903 /// matter and the order between them is not cosmetic — see the `Success`
904 /// branch, where reading it as a plain objective comparison costs a status.
905 fn honour_refused_certificate(&mut self, result: SolverReturn) -> SolverReturn {
906 if matches!(result, SolverReturn::Success) {
907 // The continued run produced a certificate of its own — but not
908 // necessarily a better *outcome*.
909 //
910 // This is what makes "never worse" hold even when the bet loses in a
911 // way that still converges: on a non-convex problem the extra travel
912 // can reach a different, worse stationary point, and the budget cap
913 // (`VETO_MAX_EXTRA_ITERS`) can also hand back a late-but-converged
914 // one. Neither may silently replace a better answer the solver
915 // already had in hand.
916 //
917 // The comparison is NOT objective-only. That was the original bug
918 // here: both points passed `passes_component_tols`, which looked
919 // like a licence to treat them as equally valid certificates and
920 // just take the lower objective. They are not equally valid when
921 // `kkt_fidelity_tol` is set — `apply_kkt_fidelity_gate` re-grades a
922 // `Success` on the unscaled KKT error afterwards, on a strictly
923 // finer criterion than the convergence test. Taking a 3-ulp
924 // objective win at a point whose unscaled error is 5x worse traded
925 // `Solve_Succeeded` for `Solved_To_Acceptable_Level`: a status
926 // regression against baseline, which is the strongest form of the
927 // guarantee breaking. So rank by the status each point will actually
928 // be *reported* under, and only then by objective.
929 let Some((refused, refused_status)) = self.baseline_outcome() else {
930 return result;
931 };
932 self.assert_comparable_scale(&refused);
933 let (curr_f, curr_kkt) = self.curr_obj_and_unscaled_kkt();
934 // Rank each candidate by the status it will actually be *reported*
935 // under, which for a `Success` means after the fidelity gate has had
936 // its say.
937 let continued_success = self.survives_fidelity_gate(curr_kkt);
938 let refused_success = matches!(refused_status, SolverReturn::Success)
939 && self.survives_fidelity_gate(refused.unscaled_kkt);
940 let keep_refused = match (continued_success, refused_success) {
941 // Equal reported status: the objective breaks the tie, which is
942 // legitimate because both points are feasible to tolerance.
943 //
944 // Negated `<=`, not `>`: they differ at NaN, and the difference
945 // matters. A `Converged` exit at an iterate whose objective is
946 // NaN but whose residuals are finite and tiny is reachable (the
947 // convergence test never inspects `f`), and `NaN > x` is false,
948 // which would keep the NaN point over a finite refused one.
949 // Phrased as a negated `<=`, an incomparable objective fails to
950 // justify keeping the continued point and the refused one wins.
951 (true, true) | (false, false) => !(curr_f <= refused.obj),
952 // The refused point keeps a status the continued one loses.
953 (false, true) => true,
954 (true, false) => false,
955 };
956 if !keep_refused {
957 return result;
958 }
959 self.restore_snapshot(&refused);
960 // The restored point's own status, which is what the baseline
961 // reported for it. For a strict refusal that is `Success` even when
962 // it fails the fidelity gate — the gate re-grades the restored point
963 // downstream, exactly as it would have re-graded the baseline's. For
964 // an acceptable-level refusal it is `StopAtAcceptablePoint`, since
965 // claiming `Success` for a point that only ever qualified at the
966 // acceptable level would over-report.
967 return refused_status;
968 }
969 // The continued run did not certify — but its final point can still be a
970 // *better* would-be certificate than the one the baseline stopped at, and
971 // restoring the chronologically-first refusal unconditionally throws it
972 // away (gh #327). The masking veto keeps refusing at the true optimum too
973 // (its unscaled error stays above `acceptable_tol` under an extreme
974 // objective scale), so a run that actually reaches the optimum never gets
975 // to certify there and instead exits non-`Success` — typically on a tiny
976 // step once it settles. Rolling straight back to the first refusal then
977 // hands back the point the baseline stopped at, which can be far worse:
978 // on `min 1/x` over `[1e-12, 10]` the solve reaches x≈10 (f≈0.1) but was
979 // rolled back to the first refusal at x≈2.84 (f≈0.35) and reported
980 // success there.
981 //
982 // The extra candidate is admitted *narrowly*, and the gate is
983 // load-bearing: the continued point may displace the refused snapshot
984 // only if it itself passes the strict per-component tolerances — i.e. it
985 // is a would-be strict certificate the veto refused solely because of
986 // masking. That is precisely what tells the settled true optimum apart
987 // from a merely lower objective reached on an unbounded ray (e.g.
988 // `A(x−a)⁴ − K·√(1+y²)`, unbounded below in y): the diverging iterate
989 // never passes the strict test, so it can never win here, and those runs
990 // stay bit-for-bit as before. When the gate does open, keep whichever
991 // point ranks better under the same feasibility-aware key the dual-guard
992 // fallback uses, and report the baseline's restored status either way —
993 // never worse than baseline on status, never worse (often better) on the
994 // point.
995 let Some((refused, restored_status)) = self.baseline_outcome() else {
996 return result;
997 };
998 self.assert_comparable_scale(&refused);
999 let curr_nlp_err = self.cq.borrow().curr_nlp_error();
1000 let curr_passes_strict =
1001 self.bundle
1002 .conv_check
1003 .current_passes_strict(curr_nlp_err, &self.data, &self.cq);
1004 let curr_f = self.cq.borrow().curr_f();
1005 let curr_viol = self.cq.borrow().curr_unscaled_primal_infeasibility_max();
1006 // The second admissible candidate: a continued run that ends *at the
1007 // acceptable level itself* (gh #533). The `curr_passes_strict` gate
1008 // exists to tell a settled optimum from a diverging ray, and on this
1009 // exit the exit itself already answers that — `StopAtAcceptablePoint` is
1010 // only reachable at a point that passed the acceptable per-component
1011 // tolerances, either by qualifying here or by being the stashed
1012 // acceptable iterate a rollback restored. A diverging iterate cannot
1013 // produce it.
1014 //
1015 // This matters because the gh #533 progress refusal is frequently paid
1016 // off by a *better acceptable point* rather than by a strict
1017 // certificate: the streak refuses while the solve is still descending,
1018 // the solve descends, and then settles somewhere better but still short
1019 // of `tol`. Without this the refused point is restored and the entire
1020 // continuation is discarded — never worse than baseline, but never
1021 // better either, which for that whole class is pure cost.
1022 //
1023 // Gated on the *restored* status also being `StopAtAcceptablePoint`, so
1024 // a strict refusal's `Success` is never reported at a point that only
1025 // ever qualified at the acceptable level.
1026 let continued_is_acceptable_exit = matches!(result, SolverReturn::StopAtAcceptablePoint)
1027 && matches!(restored_status, SolverReturn::StopAtAcceptablePoint);
1028 // Keep the continued point in place only when it is an admissible
1029 // candidate that also ranks strictly better; otherwise restore the
1030 // refused snapshot exactly as before. `ranks_better` treats a non-finite
1031 // continued objective as worst, so a NaN-objective continuation never
1032 // displaces a finite refused point (the NaN-loses convention the
1033 // `Success` branch relies on).
1034 let keep_continued = (curr_passes_strict || continued_is_acceptable_exit)
1035 && self.ranks_better(curr_f, curr_viol, refused.obj, refused.constr_viol);
1036 if !keep_continued {
1037 self.restore_snapshot(&refused);
1038 }
1039 if self.cq.borrow().curr_f().is_finite() {
1040 restored_status
1041 } else {
1042 result
1043 }
1044 }
1045
1046 /// What the baseline — the same solve with the veto disabled — would have
1047 /// returned, as (point, status), or `None` if nothing was ever refused.
1048 ///
1049 /// The **chronologically first** refusal, not the strictest one. Both arms
1050 /// follow the same trajectory until the first refusal, so that iterate is
1051 /// where the baseline stopped and what it reported. A refusal recorded later
1052 /// sits on the continued trajectory, which the baseline never walked — its
1053 /// point was never on offer, and restoring it would neither reproduce the
1054 /// baseline nor be comparable to it.
1055 ///
1056 /// Both kinds do occur, and in either order: an acceptable-level refusal
1057 /// needs `acceptable_iter` consecutive qualifying iterates, so a strict
1058 /// refusal can precede it, while a run that first drifts through the
1059 /// acceptable band can refuse there and only later pass the strict test.
1060 /// Preferring `Success` unconditionally was wrong for exactly the second
1061 /// case — it compared against a strict point from iteration 50-odd when the
1062 /// baseline had already stopped and reported acceptable at iteration 43.
1063 fn baseline_outcome(&self) -> Option<(VetoSnapshot, SolverReturn)> {
1064 // A refusal that was seen but not captured makes the baseline
1065 // unidentifiable, so decline rather than guess. Without this, a failed
1066 // strict capture alongside a successful acceptable one would present the
1067 // acceptable snapshot as the baseline outcome — but that snapshot sits
1068 // on the continued trajectory, so this would silently reintroduce the
1069 // very misidentification the chronological rule exists to prevent.
1070 // Declining loses the restore; misidentifying reports a wrong point
1071 // under a confident status.
1072 //
1073 // Unreachable today (`data.curr` is always `Some` inside `iterate()`, so
1074 // `snapshot_current` cannot fail), but the latches make the state
1075 // representable, and it must not be handled by accident.
1076 if (self.vetoed_seen && self.vetoed.is_none())
1077 || (self.vetoed_acceptable_seen && self.vetoed_acceptable.is_none())
1078 {
1079 return None;
1080 }
1081 match (&self.vetoed, &self.vetoed_acceptable) {
1082 // Ties go to the strict refusal, and the tie is reachable: both can
1083 // arm in the same call when the acceptable streak crosses on the
1084 // same iterate a strict certificate is refused. Strict is correct
1085 // there because of the baseline's own branch order — the `Converged`
1086 // gate (`opt_error.rs`, in `check_convergence_with_state`) precedes
1087 // `note_acceptable`, so the baseline returned `Converged` at that
1088 // iterate. Reordering those two branches would invert this.
1089 (Some(strict), Some(acc)) => Some(if strict.iter <= acc.iter {
1090 (strict.clone(), SolverReturn::Success)
1091 } else {
1092 (acc.clone(), SolverReturn::StopAtAcceptablePoint)
1093 }),
1094 (Some(strict), None) => Some((strict.clone(), SolverReturn::Success)),
1095 (None, Some(acc)) => Some((acc.clone(), SolverReturn::StopAtAcceptablePoint)),
1096 (None, None) => None,
1097 }
1098 }
1099
1100 /// Capture the current iterate as a veto snapshot, or `None` if there is no
1101 /// current iterate to capture.
1102 ///
1103 /// All-or-nothing by construction — see [`VetoSnapshot`].
1104 fn snapshot_current(&self, iter: Index) -> Option<VetoSnapshot> {
1105 let iterate = self.data.borrow().curr.as_ref().cloned()?;
1106 let cq = self.cq.borrow();
1107 Some(VetoSnapshot {
1108 iterate,
1109 iter,
1110 obj: cq.curr_f(),
1111 mu: self.data.borrow().curr_mu,
1112 unscaled_kkt: cq.curr_unscaled_nlp_error(),
1113 constr_viol: cq.curr_unscaled_primal_infeasibility_max(),
1114 obj_scale: cq.obj_scaling_factor(),
1115 })
1116 }
1117
1118 /// Current objective and max-norm unscaled KKT error, read together so the
1119 /// pair cannot describe different iterates.
1120 fn curr_obj_and_unscaled_kkt(&self) -> (Number, Number) {
1121 let cq = self.cq.borrow();
1122 (cq.curr_f(), cq.curr_unscaled_nlp_error())
1123 }
1124
1125 /// Guard the precondition of every scaled-objective comparison in
1126 /// [`Self::honour_refused_certificate`]: the factor must not have moved
1127 /// between the refusal and now, or the two numbers are not comparable.
1128 fn assert_comparable_scale(&self, snap: &VetoSnapshot) {
1129 debug_assert_eq!(
1130 snap.obj_scale,
1131 self.cq.borrow().obj_scaling_factor(),
1132 "objective scaling factor moved during the solve; the refused and \
1133 continued objectives are scaled differently and cannot be compared \
1134 (gh #200)"
1135 );
1136 }
1137
1138 /// Whether a point with this unscaled KKT error would keep `Solve_Succeeded`
1139 /// through [`IpoptApplication::apply_kkt_fidelity_gate`].
1140 ///
1141 /// Mirrors that gate rather than approximating it: same quantity
1142 /// (`final_unscaled_kkt_error`), same strict comparison, same "non-positive
1143 /// tolerance disables". With the default `kkt_fidelity_tol = 0` this is
1144 /// always `true`, so every caller collapses to the plain objective
1145 /// comparison and the mechanism's behaviour is unchanged.
1146 fn survives_fidelity_gate(&self, unscaled_kkt: Number) -> bool {
1147 // Phrased as the negation of the gate's own `> tol` test rather than as
1148 // `<= tol`, because the two disagree at NaN and the gate is the
1149 // authority: it demotes only on `> tol`, so a NaN error keeps `Success`
1150 // there and must keep it here. Written as `<= tol` this mirror said the
1151 // opposite, which would rank a NaN-error continued point below a refused
1152 // one. Benign in that direction — it restores the baseline point — but a
1153 // mirror that disagrees with the thing it mirrors is a latent trap.
1154 !(self.kkt_fidelity_tol > 0.0) || !(unscaled_kkt > self.kkt_fidelity_tol)
1155 }
1156
1157 /// Record the current iterate as the best acceptable-quality point seen so
1158 /// far in this solve (pounce#250 follow-up).
1159 ///
1160 /// Recording runs on **every** acceptable iterate, not only after the
1161 /// dual-divergence guard has fired. Gating it on the guard was the first
1162 /// attempt and left a hole: the guard fires and returns to the driver
1163 /// *before* this site is reached on that iteration (see the guard block in
1164 /// [`Self::iterate`]), so nothing at or before the diversion was ever
1165 /// captured. A diversion that wrecks the solve immediately — reaching no
1166 /// acceptable point afterwards — therefore had nothing to hand back, which
1167 /// is precisely the case the fallback exists for. `autocorr_bern55-06` hid
1168 /// this, because its better point happens to arrive at iteration 86, well
1169 /// after the guard fires at 23.
1170 ///
1171 /// Recording always is still behaviour-neutral, because the record is only
1172 /// ever *read* under `dual_guard_fired` — see
1173 /// [`Self::honour_best_acceptable_after_dual_guard`]. A solve the guard
1174 /// never touches computes a comparison per acceptable iterate and nothing
1175 /// else.
1176 ///
1177 /// The cost is one `f64` comparison per acceptable iterate; the iterate is
1178 /// cloned only on an actual improvement, so this does not double the
1179 /// per-iteration clone `store_acceptable_point` already pays.
1180 ///
1181 /// "Best" is a feasibility-aware ranking, **not** the lowest objective:
1182 /// candidates are ordered by [`Self::ranks_better`]'s `(feasible_enough,
1183 /// objective)` key, so objective only decides among points already inside a
1184 /// capped feasibility band. Being *bounded* by `acceptable_constr_viol_tol`
1185 /// is not the same as *not trading* feasibility within it — that band is a
1186 /// user option and can be widened to `1e1` or beyond. A pure-objective argmax
1187 /// over it has no lower bound on the feasibility it will spend, and one
1188 /// option-value away it returns a point `pounce verify` rejects under a
1189 /// `Solved_To_Acceptable_Level` status (gh #267). Whether an early
1190 /// low-objective iterate is even a candidate is the user's
1191 /// `acceptable_constr_viol_tol`; the capped feasibility key is what keeps a
1192 /// grossly-infeasible one from winning even when the band admits it.
1193 fn record_best_acceptable(&mut self, curr_f: Number) {
1194 if !curr_f.is_finite() {
1195 return;
1196 }
1197 // Same quantity the acceptable-point gate keys on, so the recorded
1198 // feasibility matches the band the candidate just passed.
1199 let curr_viol = self.cq.borrow().curr_unscaled_primal_infeasibility_max();
1200 // Reject before cloning: only a strictly better candidate — by the
1201 // feasibility-aware key, not objective alone — is worth a snapshot.
1202 if let Some(best) = self.best_acceptable.as_ref() {
1203 let (b_obj, b_viol) = (best.obj, best.constr_viol);
1204 if !self.ranks_better(curr_f, curr_viol, b_obj, b_viol) {
1205 return;
1206 }
1207 }
1208 let iter = self.data.borrow().iter_count;
1209 let Some(snap) = self.snapshot_current(iter) else {
1210 return;
1211 };
1212 // Scaled objectives are only comparable under an unchanged factor; if it
1213 // ever moved, keep the earlier point rather than compare noise.
1214 if let Some(best) = self.best_acceptable.as_ref() {
1215 if snap.obj_scale != best.obj_scale {
1216 return;
1217 }
1218 }
1219 self.best_acceptable = Some(snap);
1220 }
1221
1222 /// Cap on the feasibility band [`Self::ranks_better`] admits, matching the
1223 /// upstream default `acceptable_constr_viol_tol`. The fallback treats a point
1224 /// as "feasible enough to win on objective" only within this band, *however
1225 /// loose the user made `acceptable_constr_viol_tol`*, so widening that option
1226 /// cannot let the fallback trade feasibility for objective (gh #267).
1227 const FEASIBLE_ENOUGH_CAP: Number = 1e-2;
1228
1229 /// Whether candidate `(a_obj, a_viol)` ranks strictly better than
1230 /// `(b_obj, b_viol)` for the best-acceptable fallback (gh #267, gh #280).
1231 ///
1232 /// The key is `(band_clamped_viol, objective)` compared lexicographically,
1233 /// where each violation is clamped *up* to
1234 /// `band = min(acceptable_constr_viol_tol, FEASIBLE_ENOUGH_CAP)` before it is
1235 /// compared. Inside the band every point clamps to `band`, so they tie on
1236 /// feasibility and objective decides — objective still rules *only among
1237 /// points already feasible-enough*. Outside the band the actual violation
1238 /// decides, so the less-infeasible point always wins and a
1239 /// strictly-more-infeasible point can never rank better (gh #280 — the
1240 /// earlier `feasible_enough` partition fell through to objective-only once
1241 /// both points were outside the band). The cap keeps the band no looser than
1242 /// the upstream default: `acceptable_constr_viol_tol` is user-widenable, and
1243 /// admitting a wide band into the *objective-decides* region would let a
1244 /// grossly-infeasible low-objective iterate win. Capping the band bounds that.
1245 ///
1246 /// At default (or tighter) tolerances this is behaviour-neutral: every
1247 /// recorded point already passed the `acceptable_constr_viol_tol` gate, so
1248 /// with that band at or below the cap every candidate clamps to `band` and
1249 /// objective alone decides, exactly as before. The feasibility ordering only
1250 /// bites once the user loosens `acceptable_constr_viol_tol` past its default
1251 /// and two candidates both sit outside the cap.
1252 ///
1253 /// A non-finite objective ranks worst and can never win — feasibility never
1254 /// rescues a `NaN`/`Inf` `f`. This mirrors the `NaN`-loses convention the
1255 /// gh #200 comparisons already rely on, and it keeps a `NaN`-objective
1256 /// returned point losing to a finite recorded one in
1257 /// [`Self::honour_best_acceptable_after_dual_guard`].
1258 ///
1259 /// The ranking itself lives in the pure [`ranks_better_within_band`] so its
1260 /// never-worse-off guarantee can be proven by deterministic unit tests rather
1261 /// than inferred from a host-dependent end-to-end objective comparison (see
1262 /// gh #267, which flagged an earlier CLI test for measuring the wrong,
1263 /// host-varying property). This method only resolves the admitted band.
1264 fn ranks_better(&self, a_obj: Number, a_viol: Number, b_obj: Number, b_viol: Number) -> bool {
1265 let band = self
1266 .bundle
1267 .conv_check
1268 .acceptable_constr_viol_tol_or_default()
1269 .min(Self::FEASIBLE_ENOUGH_CAP);
1270 ranks_better_within_band(a_obj, a_viol, b_obj, b_viol, band)
1271 }
1272
1273 /// Make the dual-divergence guard's diversion non-destructive (pounce#250
1274 /// follow-up).
1275 ///
1276 /// The guard bets that routing to restoration beats grinding on, and nothing
1277 /// made losing that bet safe: it could return a materially worse point than
1278 /// the solve already had, under a status that does not admit it.
1279 ///
1280 /// WHAT THIS DOES AND DOES NOT GUARANTEE. It guarantees the diverted run
1281 /// never returns worse than the best acceptable-quality point **that same
1282 /// run visited**. It does *not* guarantee the diverted run is no worse than
1283 /// not diverting at all — that counterfactual solve never happened, and its
1284 /// points were never on offer to compare against. The distinction is not
1285 /// academic: on the Linux CI host `deb7` returns 97.56 with the guard off and
1286 /// 127.87 with it on at streak 15, and this fallback cannot close that gap,
1287 /// because 127.87 is the best acceptable point the diverted run ever reached.
1288 /// Bounding the diversion's damage is a weaker property than making the
1289 /// diversion harmless, and only the weaker one is available from inside a
1290 /// single solve. It is a large part of why the guard is off by default.
1291 ///
1292 /// The observed case is `autocorr_bern55-06`. The guard fires at iteration
1293 /// 23, the diverted run reaches the true optimum (-2304.0000278, matching
1294 /// Ipopt to 12 significant figures) and holds it from iteration 57 to 86 —
1295 /// but the dual residual sawtooths between 1e-8 and 2e-1 there, so it never
1296 /// strings together the `acceptable_iter` consecutive qualifying iterates
1297 /// that would stop the solve. It then enters restoration a second time,
1298 /// wanders into a worse basin, and terminates `StopAtAcceptablePoint` at
1299 /// -2263.46 — 1.8 % worse, with an overall NLP error of 1.0. The better
1300 /// point was *visited and passed the acceptable test*; it was simply
1301 /// overwritten, because `store_acceptable_point` keeps the latest rather
1302 /// than the best.
1303 ///
1304 /// So: on a non-`Success` exit, if the best acceptable-quality iterate seen
1305 /// anywhere in the solve beats the point being returned, hand that back
1306 /// instead. This is the same "never worse off" bargain the gh #200 veto
1307 /// makes, applied to the other bet in the algorithm.
1308 ///
1309 /// "Beats" is the feasibility-aware ranking in [`Self::ranks_better`], not a
1310 /// bare objective comparison: the recorded point wins only if it is
1311 /// feasible-enough while the returned point is not, or both are in the same
1312 /// feasibility class and it has a lower objective. Ranking by objective alone
1313 /// let a widened `acceptable_constr_viol_tol` band trade feasibility for
1314 /// objective here — restoring a lower-objective point that `pounce verify`
1315 /// rejects, under a success-mapped status (gh #267). The key prevents that:
1316 /// objective can only win among points already inside the capped acceptable
1317 /// feasibility band.
1318 ///
1319 /// Note "anywhere in the solve", not "since the guard fired":
1320 /// [`Self::record_best_acceptable`] runs unconditionally and explains why —
1321 /// points at or before the diversion have to be on offer, or a diversion that
1322 /// wrecks the solve immediately has nothing to hand back. Only this *read* is
1323 /// gated on `dual_guard_fired`.
1324 ///
1325 /// A strict `Success` is never overridden — that point carries a real
1326 /// certificate, and a lower objective at a merely-acceptable point must not
1327 /// displace it.
1328 ///
1329 /// Tuning the guard's firing threshold was tried first and rejected: no
1330 /// setting separates the models it helps from the ones it harms, and the
1331 /// effect turned out to differ by host anyway (see the option help in
1332 /// `upstream_options.rs`). Fixing the consequence is what remained available.
1333 fn honour_best_acceptable_after_dual_guard(&mut self, result: SolverReturn) -> SolverReturn {
1334 if !self.dual_guard_fired || matches!(result, SolverReturn::Success) {
1335 return result;
1336 }
1337 let Some(best) = self.best_acceptable.clone() else {
1338 return result;
1339 };
1340 let (curr_f, _) = self.curr_obj_and_unscaled_kkt();
1341 let curr_viol = self.cq.borrow().curr_unscaled_primal_infeasibility_max();
1342 let curr_scale = self.cq.borrow().obj_scaling_factor();
1343 // Only comparable under the same factor, sign included.
1344 if curr_scale != best.obj_scale {
1345 return result;
1346 }
1347 // Restore only when the recorded point ranks strictly better under the
1348 // feasibility-aware key — more feasible, or equally feasible at a lower
1349 // objective. `ranks_better` also handles the `NaN` case the previous
1350 // bare `!(curr_f <= best.obj)` did: a non-finite returned objective
1351 // ranks worst, so a finite recorded point wins and is restored.
1352 if self.ranks_better(best.obj, best.constr_viol, curr_f, curr_viol) {
1353 tracing::debug!(target: "pounce::algorithm",
1354 "[POUNCE] dual-divergence diversion ended worse than a point already \
1355 in hand (obj {:.10e} viol {:.3e} -> obj {:.10e} viol {:.3e}, iter {}); \
1356 restoring it (pounce#250, gh#267).",
1357 curr_f, curr_viol, best.obj, best.constr_viol, best.iter,
1358 );
1359 self.restore_snapshot(&best);
1360 // Swap the *point*, but never let the swap erase why the solve
1361 // stopped. A budget that was exhausted stays reported as exhausted:
1362 // a caller polling for "did I run out of time" must not be told
1363 // "solved to acceptable level" merely because a better point was
1364 // recoverable. Only the outcomes that carry no such fact of their
1365 // own are relabelled to describe what is now being returned.
1366 return match result {
1367 SolverReturn::MaxiterExceeded
1368 | SolverReturn::CpuTimeExceeded
1369 | SolverReturn::WallTimeExceeded
1370 | SolverReturn::UserRequestedStop => result,
1371 _ => SolverReturn::StopAtAcceptablePoint,
1372 };
1373 }
1374 result
1375 }
1376
1377 /// Make a refused snapshot the current iterate again.
1378 fn restore_snapshot(&mut self, snap: &VetoSnapshot) {
1379 let mut d = self.data.borrow_mut();
1380 d.set_trial(snap.iterate.clone());
1381 d.accept_trial_point();
1382 // The restored point's own barrier parameter, not the continued run's —
1383 // see `VetoSnapshot::mu`.
1384 d.curr_mu = snap.mu;
1385 }
1386
1387 /// Decide whether the acceptable-point restoration decline may be deferred
1388 /// this once (gh #534), and arm the bookkeeping that bounds the bet.
1389 ///
1390 /// Four conditions, all required:
1391 ///
1392 /// * the option leaves deferrals available at all
1393 /// (`resto_decline_deferrals`, `0` = pre-#534 behaviour);
1394 /// * the budget is not already spent;
1395 /// * the NLP error has contracted on every one of the last
1396 /// [`DECLINE_PROGRESS_SAMPLES`]` - 1` iterations
1397 /// ([`Self::nlp_err_contracting`]) — the progress test the guard lacked;
1398 /// * the iteration budget has room for a continuation, and the entry point
1399 /// can actually be captured. Without a floor there is nothing to fall
1400 /// back to, and a bet with no floor is exactly what must not be placed.
1401 ///
1402 /// The deadline is clamped below `max_iter` so a lost bet can never turn a
1403 /// reportable `StopAtAcceptablePoint` into `Maximum_Iterations_Exceeded`:
1404 /// the continuation is always cut before the iteration budget runs out. A
1405 /// *time* budget is not clamped the same way — elapsed time is an external
1406 /// fact and the deadline cannot predict it — so a solve that expires inside
1407 /// the continuation window still reports the time limit, at the floor
1408 /// iterate rather than at whatever the continuation last touched.
1409 fn may_defer_acceptable_decline(&mut self) -> bool {
1410 if self.decline_deferrals_used >= self.resto_decline_deferrals {
1411 return false;
1412 }
1413 if !self.nlp_err_contracting() {
1414 return false;
1415 }
1416 let iter = self.data.borrow().iter_count;
1417 // No room to continue: the deadline below would fire on the very next
1418 // iteration, so the deferral would buy nothing and cost a restoration.
1419 if iter.saturating_add(1) >= self.max_iter {
1420 return false;
1421 }
1422 if self.decline_floor.is_none() {
1423 let Some(snap) = self.snapshot_current(iter) else {
1424 return false;
1425 };
1426 self.decline_floor = Some(snap);
1427 }
1428 self.decline_deferrals_used += 1;
1429 self.decline_deadline_iter = Some(
1430 iter.saturating_add(DECLINE_CONTINUATION_BUDGET)
1431 .min(self.max_iter.saturating_sub(1)),
1432 );
1433 true
1434 }
1435
1436 /// The deferred continuation ran out of budget without a strict certificate
1437 /// (gh #534). Report the floor — the point the pre-#534 guard would have
1438 /// returned — unless the continuation is standing somewhere at least as
1439 /// good.
1440 fn terminate_at_decline_floor(&mut self) -> IterateOutcome {
1441 let Some(floor) = self.decline_floor.clone() else {
1442 // Unreachable in practice: the deadline is only ever set after a
1443 // floor is captured. Stopping at the current point is still the
1444 // right thing if it somehow is not — the point passed the
1445 // acceptable-level triplet when the deferral was taken.
1446 return IterateOutcome::Terminate(SolverReturn::StopAtAcceptablePoint);
1447 };
1448 tracing::debug!(target: "pounce::algorithm",
1449 "[POUNCE] deferred restoration decline expired at iter {} without a strict \
1450 certificate; falling back to the floor from iter {} (gh #534).",
1451 self.data.borrow().iter_count, floor.iter,
1452 );
1453 if !self.continuation_outranks(&floor) {
1454 self.restore_snapshot(&floor);
1455 }
1456 IterateOutcome::Terminate(SolverReturn::StopAtAcceptablePoint)
1457 }
1458
1459 /// Whether the current iterate is a *better* answer than the gh #534 floor.
1460 ///
1461 /// Two gates, in order. The current point must itself pass the
1462 /// acceptable-level triplet — the floor is going to be reported under
1463 /// `Solved_To_Acceptable_Level`, and a continuation that wandered off is not
1464 /// entitled to that status however attractive its objective looks. Only then
1465 /// does [`Self::ranks_better`]'s feasibility-first key decide, and only under
1466 /// an unmoved objective scaling factor, since the two objectives are
1467 /// otherwise not comparable.
1468 fn continuation_outranks(&self, floor: &VetoSnapshot) -> bool {
1469 let (curr_f, _) = self.curr_obj_and_unscaled_kkt();
1470 if !curr_f.is_finite() {
1471 return false;
1472 }
1473 let nlp_err = self.cq.borrow().curr_nlp_error();
1474 if !self
1475 .bundle
1476 .conv_check
1477 .current_is_acceptable_with_state(nlp_err, &self.data, &self.cq)
1478 {
1479 return false;
1480 }
1481 let (curr_viol, curr_scale) = {
1482 let cq = self.cq.borrow();
1483 (
1484 cq.curr_unscaled_primal_infeasibility_max(),
1485 cq.obj_scaling_factor(),
1486 )
1487 };
1488 if curr_scale != floor.obj_scale {
1489 return false;
1490 }
1491 !self.ranks_better(floor.obj, floor.constr_viol, curr_f, curr_viol)
1492 }
1493
1494 /// Make a deferred restoration decline non-destructive (gh #534).
1495 ///
1496 /// The deferral is a bet that a contracting endgame is three iterations from
1497 /// a certificate. This is what makes losing it cost only those iterations:
1498 /// whatever the continued run ends up returning, if it is not at least as
1499 /// good an answer as the floor — the point the pre-#534 guard would have
1500 /// reported — the floor is restored and reported instead.
1501 ///
1502 /// Applied once, in [`Self::optimize`], for the same reason the gh #200 and
1503 /// pounce#250 hooks are: the driver loop has many `return`s and this must
1504 /// see all of them.
1505 ///
1506 /// A strict `Success` is never overridden — that is the bet paying off, and
1507 /// a real certificate outranks any acceptable-level point by construction.
1508 /// The budget statuses keep their own status, as they do in
1509 /// [`Self::honour_best_acceptable_after_dual_guard`]: a caller polling for
1510 /// "did I run out of time" must be told so, even while the point it gets
1511 /// back is swapped for the better one.
1512 fn honour_decline_floor(&mut self, result: SolverReturn) -> SolverReturn {
1513 if matches!(result, SolverReturn::Success) {
1514 if self.decline_floor.is_some() {
1515 tracing::debug!(target: "pounce::algorithm",
1516 "[POUNCE] the deferred restoration decline paid off: the continuation \
1517 reached a strict certificate (gh #534).",
1518 );
1519 }
1520 return result;
1521 }
1522 let Some(floor) = self.decline_floor.clone() else {
1523 return result;
1524 };
1525 if self.continuation_outranks(&floor) {
1526 return result;
1527 }
1528 tracing::debug!(target: "pounce::algorithm",
1529 "[POUNCE] the deferred restoration decline did not pay off; restoring the \
1530 floor from iter {} (obj {:.10e} viol {:.3e}) and reporting it (gh #534).",
1531 floor.iter, floor.obj, floor.constr_viol,
1532 );
1533 self.restore_snapshot(&floor);
1534 match result {
1535 SolverReturn::MaxiterExceeded
1536 | SolverReturn::CpuTimeExceeded
1537 | SolverReturn::WallTimeExceeded
1538 | SolverReturn::UserRequestedStop => result,
1539 _ => SolverReturn::StopAtAcceptablePoint,
1540 }
1541 }
1542
1543 /// Terminal fallback for a near-feasible numerical breakdown (a
1544 /// restoration cycle or a failed step computation). If a finite
1545 /// acceptable iterate was recorded earlier in the solve, roll back
1546 /// to it and stop at [`SolverReturn::StopAtAcceptablePoint`] (mapped
1547 /// by the application layer to `Solved_To_Acceptable_Level`) rather
1548 /// than surfacing the hard `fallback` error. This mirrors upstream
1549 /// `IpBacktrackingLineSearch`'s `ACCEPTABLE_POINT_REACHED`
1550 /// precedence: when the line search exhausts but an acceptable point
1551 /// was stored, that point is returned instead of the failure. With
1552 /// no snapshot — or if the restored objective is non-finite — the
1553 /// original `fallback` status is surfaced unchanged, so genuinely
1554 /// failed/infeasible solves keep their honest status. Catches
1555 /// degenerate LPs (kleemin8, nsir2) whose μ-endgame reaches the
1556 /// optimum, then destabilizes on the ill-conditioned vertex and
1557 /// cycles in restoration instead of stopping at the acceptable
1558 /// iterate it already passed through.
1559 fn terminate_acceptable_or(&mut self, fallback: SolverReturn) -> IterateOutcome {
1560 if self.restore_acceptable_point() && self.cq.borrow().curr_f().is_finite() {
1561 IterateOutcome::Terminate(SolverReturn::StopAtAcceptablePoint)
1562 } else {
1563 IterateOutcome::Terminate(fallback)
1564 }
1565 }
1566
1567 /// The single place this module turns a local-infeasibility conclusion
1568 /// into a returned status (gh #505).
1569 ///
1570 /// Three separate routes reach that verdict — the conv-check's rapid
1571 /// detection, restoration layer 2, and the slow-cycle exits — and the same
1572 /// defect was found in two of them independently: returning the hard
1573 /// verdict without consulting the acceptable-point stash, so a solve that
1574 /// had already passed through an acceptable iterate discarded it. Only the
1575 /// cycle exits got it right, and nothing structural said the other two were
1576 /// wrong.
1577 ///
1578 /// That is the shape of a defect that comes back. The route a solve takes
1579 /// to the verdict is an internal detail — the user sees one status either
1580 /// way — so the *decision* about what that status means must not live at
1581 /// each route. It lives here, and
1582 /// [`no_route_concludes_local_infeasibility_alone`] is a tripwire against
1583 /// a new site rebuilding the outcome inline.
1584 ///
1585 /// The cycle exits are not routed through here because their fallback is
1586 /// chosen between `LocalInfeasibility` and `ErrorInStepComputation` at the
1587 /// call site; they already reach `terminate_acceptable_or`, which is the
1588 /// behaviour this guarantees.
1589 ///
1590 /// Scope: this governs how *this module* returns the verdict. Other layers
1591 /// name `SolverReturn::LocalInfeasibility` for their own reasons — the SQP
1592 /// status map and the ℓ₁ elastic path in `application.rs`, for instance —
1593 /// and are outside both this funnel and its tripwire.
1594 fn terminate_local_infeasibility(&mut self) -> IterateOutcome {
1595 self.terminate_acceptable_or(SolverReturn::LocalInfeasibility)
1596 }
1597
1598 pub fn with_nlp(mut self, nlp: Rc<RefCell<dyn IpoptNlp>>) -> Self {
1599 self.nlp = Some(nlp);
1600 self
1601 }
1602
1603 /// Install a user-facing TNLP handle. Enables per-iteration
1604 /// `TNLP::intermediate_callback` invocation from `optimize()`.
1605 pub fn with_tnlp(mut self, tnlp: Rc<RefCell<dyn TNLP>>) -> Self {
1606 self.tnlp = Some(tnlp);
1607 self
1608 }
1609
1610 /// Build an [`IterStats`] payload from the current `IpoptData` /
1611 /// `IpoptCq` state. Mirrors the field set the upstream Ipopt main
1612 /// loop hands to `IntermediateCallback` after each `AcceptTrialPoint`.
1613 fn build_iter_stats(&self) -> IterStats {
1614 let d = self.data.borrow();
1615 let c = self.cq.borrow();
1616 let dnrm = match d.delta.as_ref() {
1617 Some(delta) => delta.x.amax().max(delta.s.amax()),
1618 None => 0.0,
1619 };
1620 IterStats {
1621 // alg_mod tracking (regular vs restoration) is a follow-up;
1622 // every callback fire from the outer loop reports RegularMode.
1623 mode: AlgorithmMode::RegularMode,
1624 iter: d.iter_count,
1625 obj_value: c.curr_f(),
1626 inf_pr: c.curr_primal_infeasibility_max(),
1627 inf_du: c.curr_dual_infeasibility_max(),
1628 mu: d.curr_mu,
1629 d_norm: dnrm,
1630 regularization_size: d.info_regu_x,
1631 alpha_du: d.info_alpha_dual,
1632 alpha_pr: d.info_alpha_primal,
1633 ls_trials: d.info_ls_count,
1634 }
1635 }
1636
1637 /// Fire `TNLP::intermediate_callback` if a TNLP handle and NLP
1638 /// handle are installed. Wraps the call in an [`IntermediateContext`]
1639 /// guard so downstream inspector entry points (the C API's
1640 /// `GetIpoptCurrent*`) can read live state for the duration. Returns
1641 /// `true` to continue, `false` if the user requested termination.
1642 fn fire_intermediate(&self) -> bool {
1643 let Some(tnlp) = self.tnlp.as_ref() else {
1644 return true;
1645 };
1646 let Some(nlp) = self.nlp.as_ref() else {
1647 return true;
1648 };
1649 let stats = self.build_iter_stats();
1650 let _guard = CtxGuard::install(IntermediateContext {
1651 data: Rc::clone(&self.data),
1652 cq: Rc::clone(&self.cq),
1653 nlp: Rc::clone(nlp),
1654 });
1655 tnlp.borrow_mut().intermediate_callback(
1656 stats,
1657 &TnlpIpoptData::default(),
1658 &TnlpIpoptCq::default(),
1659 )
1660 }
1661
1662 pub fn with_search_dir(mut self, sd: PdSearchDirCalc) -> Self {
1663 self.search_dir = Some(sd);
1664 self
1665 }
1666
1667 pub fn with_restoration(mut self, resto: Box<dyn RestorationPhase>) -> Self {
1668 self.restoration = Some(resto);
1669 self
1670 }
1671
1672 /// Install the shared diagnostics state. The state is propagated
1673 /// to the augmented-system solver at the top of [`Self::optimize`]
1674 /// so dump sites can consult per-iter gating.
1675 pub fn with_diagnostics(mut self, diag: Rc<DiagnosticsState>) -> Self {
1676 self.diagnostics = Some(diag);
1677 self
1678 }
1679
1680 /// Install an interactive debugger hook. Fired at each checkpoint
1681 /// in [`Self::optimize`]; returning [`crate::debug::DebugAction::Stop`]
1682 /// ends the solve with `SolverReturn::UserRequestedStop`.
1683 pub fn with_debug_hook(mut self, hook: Rc<RefCell<dyn crate::debug::DebugHook>>) -> Self {
1684 self.debug = Some(hook);
1685 self
1686 }
1687
1688 /// Shared handle to the installed debugger, if any — used to forward
1689 /// it into the restoration inner IPM.
1690 pub fn debug_hook(&self) -> Option<Rc<RefCell<dyn crate::debug::DebugHook>>> {
1691 self.debug.as_ref().map(Rc::clone)
1692 }
1693
1694 /// Fire the debugger hook (if installed) at `cp`, building a live
1695 /// [`crate::debug::DebugCtx`] over cheap handle clones. Returns the
1696 /// requested action, defaulting to `Resume` when no hook is set.
1697 fn fire_debug(&mut self, cp: crate::debug::Checkpoint) -> crate::debug::DebugAction {
1698 use crate::debug::{DebugAction, DebugCtx};
1699 // Clone the Rc so the hook borrow is released before we touch
1700 // `self.bundle` to apply any live option changes below.
1701 let Some(hook) = self.debug.as_ref().map(Rc::clone) else {
1702 return DebugAction::Resume;
1703 };
1704 let mut ctx = DebugCtx::new(Rc::clone(&self.data), Rc::clone(&self.cq), cp);
1705 let action = hook.borrow_mut().at_checkpoint(&mut ctx);
1706 // Drain any tolerances the hook asked to hot-swap and write them
1707 // into the live convergence-check policy, so the next iteration's
1708 // termination test uses the new value (no `resolve` needed).
1709 for (name, value) in ctx.take_live_tolerances() {
1710 self.bundle.conv_check.set_tolerance(&name, value);
1711 }
1712 action
1713 }
1714
1715 /// Run the restoration phase, bracketed by the `PreRestoration` /
1716 /// `PostRestoration` debug checkpoints so a debugger can inspect the
1717 /// iterate just before entry and just after exit. With no debugger
1718 /// installed this is exactly `invoke_restoration()`.
1719 fn invoke_restoration_debugged(&mut self) -> IterateOutcome {
1720 if let Some(o) = self.debug_stop(crate::debug::Checkpoint::PreRestoration) {
1721 return o;
1722 }
1723 let outcome = self.invoke_restoration();
1724 if let Some(o) = self.debug_stop(crate::debug::Checkpoint::PostRestoration) {
1725 return o;
1726 }
1727 outcome
1728 }
1729
1730 /// Fire a sub-iteration checkpoint from inside [`Self::iterate`].
1731 /// Returns `Some(Terminate(UserRequestedStop))` if the debugger asked
1732 /// to stop, so the caller can `return` it; `None` to continue.
1733 fn debug_stop(&mut self, cp: crate::debug::Checkpoint) -> Option<IterateOutcome> {
1734 if self.debug.is_none() {
1735 return None;
1736 }
1737 if self.fire_debug(cp) == crate::debug::DebugAction::Stop {
1738 Some(IterateOutcome::Terminate(SolverReturn::UserRequestedStop))
1739 } else {
1740 None
1741 }
1742 }
1743
1744 /// Fire the terminal post-mortem checkpoint (if a debugger is set),
1745 /// carrying the solve outcome so the hook can decide whether to pause
1746 /// at the final iterate. The action is advisory — the loop returns
1747 /// `result` regardless — so the hook just gets a last look.
1748 fn fire_debug_terminal(&mut self, result: SolverReturn) {
1749 use crate::debug::{Checkpoint, DebugCtx};
1750 let Some(hook) = self.debug.as_ref() else {
1751 return;
1752 };
1753 let mut ctx = DebugCtx::new(
1754 Rc::clone(&self.data),
1755 Rc::clone(&self.cq),
1756 Checkpoint::Terminated,
1757 )
1758 .with_status(format!("{result:?}"));
1759 let _ = hook.borrow_mut().at_checkpoint(&mut ctx);
1760 }
1761
1762 /// Cheap mid-iteration time-budget check (pounce#242). Returns the
1763 /// terminal [`SolverReturn`] when the shared [`Deadline`] has been
1764 /// crossed, so the caller can bail *within* an iteration — after the
1765 /// KKT factorization, before the line search — rather than only at the
1766 /// next outer-iteration convergence check. Returns `None` (never
1767 /// terminating) when no deadline is installed, keeping the
1768 /// direct-driver / unit-test paths on their `overall_alg`-based gate.
1769 /// A `None` here is not "no budget" but "check it at the coarse site".
1770 fn deadline_status(&self) -> Option<SolverReturn> {
1771 let d = self.data.borrow();
1772 let kind = d.deadline.as_ref()?.exceeded()?;
1773 Some(match kind {
1774 pounce_common::timing::DeadlineKind::Cpu => SolverReturn::CpuTimeExceeded,
1775 pounce_common::timing::DeadlineKind::Wall => SolverReturn::WallTimeExceeded,
1776 })
1777 }
1778
1779 /// One iteration body — port of `Optimize()`'s inner loop.
1780 /// Returns either `Continue` to keep iterating or a terminal
1781 /// [`SolverReturn`] mirroring upstream's exception → return-code
1782 /// translation table (see `MAIN_LOOP.md` §"Exception mapping").
1783 fn iterate(&mut self) -> IterateOutcome {
1784 // Shared timing accumulator — cheap Rc clone so each phase can
1785 // bump its own counter without re-borrowing `data`.
1786 let timing = self.data.borrow().timing.clone();
1787
1788 // Per-iteration span so every event emitted in this body (the
1789 // structured iteration record, restoration/linear-solve spans)
1790 // is tagged with the iteration index.
1791 let _iter_span =
1792 tracing::info_span!("iteration", iter = self.data.borrow().iter_count).entered();
1793
1794 // 1. Output iteration row. Header every 10 iters; the row itself
1795 // is built plain by the strategy (so the column widths stay
1796 // exact and unit-testable) and wrapped in a tiger/rust style
1797 // at the print site (pounce#71). `anstream::stdout()` strips
1798 // the escapes automatically when stdout is redirected or
1799 // `NO_COLOR` is set, so non-TTY output is plain text.
1800 //
1801 // Print BEFORE `reset_info` so the row reflects the accepted
1802 // step from the previous iteration (alphas, ls count,
1803 // alpha_char), matching upstream's `IpIpoptAlgorithm::Optimize`
1804 // ordering.
1805 timing.output_iteration.start();
1806 self.bundle.iter_output.write_output();
1807 if self.print_iter_output {
1808 use std::io::Write as _;
1809 let (iter_count, alpha_pr, alpha_char) = {
1810 let d = self.data.borrow();
1811 (d.iter_count, d.info_alpha_primal, d.info_alpha_primal_char)
1812 };
1813 let row = self.bundle.iter_output.format_row(&self.data, &self.cq);
1814 // Iteration 0 is the initial point — no step has been taken
1815 // yet, so `alpha_primal` is 0; treat it as a full step
1816 // (neutral black) rather than a stalling alarm (red).
1817 let style_alpha = if iter_count == 0 { 1.0 } else { alpha_pr };
1818 let style = pounce_common::style::iteration_row_style(style_alpha, alpha_char);
1819 let mut out = anstream::stdout();
1820 // Write errors (e.g. a closed pipe / `head` on the output)
1821 // are deliberately ignored: a vanished terminal must not
1822 // panic the solver, unlike the old `println!`.
1823 if iter_count % 10 == 0 {
1824 let _ = write!(out, "{}", crate::output::orig::OrigIterationOutput::HEADER);
1825 }
1826 let _ = writeln!(out, "{}{}{}", style.render(), row, style.render_reset());
1827 }
1828 timing.output_iteration.end();
1829
1830 // Structured per-iteration event (pounce#71) — the single source
1831 // of truth for the per-iteration trajectory. The JSON log sink
1832 // and the solve-report collector
1833 // (`pounce_observability::IterCollectorLayer`) both derive from
1834 // it. The text console layer filters this target out (its human
1835 // form is the colored table above).
1836 //
1837 // Skipped entirely when nothing consumes it (no iter-history
1838 // capture active and JSON logging off) so the default run pays
1839 // no per-iteration field-evaluation / allocation cost.
1840 if pounce_observability::iteration_event_wanted() {
1841 let d = self.data.borrow();
1842 let c = self.cq.borrow();
1843 let alpha_char = d.info_alpha_primal_char;
1844 let alpha_char_s = alpha_char.to_string();
1845 let d_norm = match &d.delta {
1846 Some(delta) => delta.x.amax().max(delta.s.amax()),
1847 None => 0.0,
1848 };
1849 tracing::info!(
1850 target: pounce_observability::ITER_TARGET,
1851 iter = d.iter_count,
1852 objective = c.unscaled_curr_f(),
1853 inf_pr = c.curr_primal_infeasibility_max(),
1854 inf_du = c.curr_dual_infeasibility_max(),
1855 mu = d.curr_mu,
1856 d_norm = d_norm,
1857 regularization = d.info_regu_x,
1858 alpha_dual = d.info_alpha_dual,
1859 alpha_primal = d.info_alpha_primal,
1860 ls_trials = d.info_ls_count,
1861 alpha_char = alpha_char_s.as_str(),
1862 resto_kind = pounce_common::style::resto_kind_str(alpha_char),
1863 );
1864 }
1865
1866 // Reset per-iteration info on data (after printing previous
1867 // iter's accepted-step info; before the next line search).
1868 self.data.borrow_mut().reset_info();
1869
1870 // 2. Convergence check.
1871 timing.check_convergence.start();
1872 let nlp_err = self.cq.borrow().curr_nlp_error();
1873 let iter_count = self.data.borrow().iter_count;
1874 if !nlp_err.is_finite() {
1875 timing.check_convergence.end();
1876 return IterateOutcome::Terminate(SolverReturn::InvalidNumberDetected);
1877 }
1878 // gh #534 progress history. One sample per outer iteration, recorded
1879 // before any of the guards below can divert, so the samples the
1880 // restoration-decline test reads are consecutive by construction.
1881 self.note_nlp_err(nlp_err);
1882 // Divergence guard — port of upstream `IpIpoptAlg.cpp` post-
1883 // AcceptTrialPoint check. When `max_i |x_i|` exceeds the
1884 // registered `diverging_iterates_tol` (default `1e20`), exit
1885 // cleanly with `DivergingIterates` rather than spiralling into
1886 // a degenerate restoration whose inner sub-NLP can't recover
1887 // (MESH: orig `f` already at -3.6e33 by iter 90, restoration
1888 // entered too late to bound `x`).
1889 //
1890 // A large `|x|` alone does not prove unboundedness, though:
1891 // `DivergingIterates` is Ipopt's *unboundedness* signal (it maps
1892 // to the AMPL 300 "unbounded" range), and under severe objective
1893 // ill-scaling the normal-mode IPM can take a large but transient
1894 // excursion on a problem that is bounded below with a finite
1895 // optimum (issue #248: MINLPLib `jit1`). Only conclude divergence
1896 // when the growth is *structurally* consistent with an unbounded
1897 // feasible region — some over-threshold component heading toward a
1898 // side with no finite bound. If every large component is pinned by
1899 // a finite bound (in particular, all variables boxed), the growth
1900 // is a scaling artifact, so let the normal convergence / iteration
1901 // machinery return the best iterate instead of a spurious
1902 // `Unbounded`.
1903 // Evaluate the structural check under an immutable borrow, then
1904 // update the persistence state and take the verdict separately so
1905 // the mutable field updates don't clash with the `data` borrow.
1906 // Two independent unboundedness paths share this block:
1907 // * the `diverging_iterates_tol` (`1e20`) magnitude guard, gated on
1908 // the free-variable structural check + geometric-growth streak
1909 // (issues #248 / #252); and
1910 // * the #285 recession-ray path — a *checked proof*, active from a
1911 // far lower magnitude floor, that catches a genuine recession ray
1912 // in `null(A_eq)` over free variables whose `|x|` grows only
1913 // linearly and so never reaches `1e20` within `max_iter`.
1914 let (amax, structural_free, is_ray) = {
1915 let data = self.data.borrow();
1916 match data.curr.as_ref() {
1917 Some(curr) => {
1918 let amax = curr.x.amax();
1919 let structural = amax > self.diverging_iterates_tol
1920 && self.divergence_is_true_unboundedness(&*curr.x);
1921 let is_ray = amax > Self::RECESSION_MIN_NORM
1922 && self.curr_is_recession_ray(&*curr.x, amax);
1923 (Some(amax), structural, is_ray)
1924 }
1925 None => (None, false, false),
1926 }
1927 };
1928 // Evaluate the (scaled) objective only while a structural divergence
1929 // is live — the streak's descent gate needs it, and it costs an
1930 // objective evaluation, so skip it on the common non-diverging path.
1931 let curr_f = structural_free.then(|| self.cq.borrow().curr_f());
1932 // Evaluate both streak updates (no short-circuit) so each keeps its
1933 // state current, then fire if either concludes divergence.
1934 let fire_magnitude = self.update_divergence_verdict(amax, structural_free, curr_f);
1935 let fire_recession = self.update_recession_verdict(amax.unwrap_or(0.0), is_ray);
1936 if fire_magnitude || fire_recession {
1937 if fire_recession && !fire_magnitude {
1938 tracing::debug!(target: "pounce::algorithm",
1939 "[POUNCE] recession-ray guard fired at iter {} (|x|_inf={:.2e}); \
1940 reporting DivergingIterates (pounce#285).",
1941 self.data.borrow().iter_count,
1942 amax.unwrap_or(f64::NAN),
1943 );
1944 }
1945 timing.check_convergence.end();
1946 return IterateOutcome::Terminate(SolverReturn::DivergingIterates);
1947 }
1948 // Dual-divergence guard (pounce#246). The primal guard above only
1949 // catches `|x|` blowing up; a bad warm start can instead send the
1950 // *dual* infeasibility diverging — `inf_du` 1 -> 1e14, the inertia
1951 // regularization -> 1e14, the barrier parameter frozen, full steps
1952 // still accepted by the filter because primal feasibility inches
1953 // down — while `|x|` stays bounded. `diverging_iterates_tol` never
1954 // trips, restoration is never entered, and the solve grinds in
1955 // ever-more-ill-conditioned KKT factorizations that each take
1956 // seconds (the emfl050 warm-start overshoot: one 3.8 s factorization
1957 // per iteration, forever). Detect a sustained streak of growing dual
1958 // infeasibility in the elevated regime and route to restoration —
1959 // the same recovery the least-square-multiplier init path reaches on
1960 // its own — before the factorizations start choking. Gated on a
1961 // large absolute `inf_du` so a well-behaved solve whose dual
1962 // residual transiently rises (then falls) is never diverted:
1963 // restoration is a heavier hammer than the guard should swing at a
1964 // merely-bumpy-but-converging iterate.
1965 //
1966 // OFF BY DEFAULT (pounce#250 follow-up). The emfl050 overshoot above is
1967 // how this was justified, and it did not reproduce: that measurement was
1968 // caller-side JAX compilation, and the build predating the guard solves
1969 // both emfl050 instances to the same optimum in the same time. What is
1970 // left is an effect on four of 1284 MINLPLib models that is knife-edge
1971 // and non-monotone in `dual_diverging_streak` — a better local optimum on
1972 // deb7/deb9 at exactly 15, and pooling_rt2stp turning Solve_Succeeded
1973 // into Maximum_Iterations_Exceeded at 10 and 15 only. Kept because it
1974 // does help when it helps, but not imposed. Full account in the option
1975 // help (`upstream_options.rs`).
1976 //
1977 // Two things to know before changing this:
1978 //
1979 // * `curr_dual_infeasibility_max` is the RAW ‖∇L‖∞, not divided by the
1980 // `s_d` optimality scaling the convergence check applies, and this runs
1981 // *before* `conv_check`. So the thresholds below are not on the same
1982 // quantity the solver's own tolerances are on, and the claim that they
1983 // are scale-robust holds only while `nlp_scaling_method != none`. No
1984 // exploit is known; the margin is thinner than it looks.
1985 // * The `DivergingIterates` fallback at the end is unreachable from every
1986 // shipped front end — CLI, pounce-py and cinterface all wire a
1987 // restoration provider, so the guard can only ever route to
1988 // restoration. Do not assume it is dead code and delete the provider
1989 // check; do not assume it is live and rely on the status either.
1990 if self.dual_diverging_streak > 0 {
1991 let inf_du = self.cq.borrow().curr_dual_infeasibility_max();
1992 if inf_du.is_finite() && inf_du > self.dual_inf_prev && inf_du > DUAL_DIV_COUNT_FLOOR {
1993 self.dual_growth_streak += 1;
1994 } else {
1995 self.dual_growth_streak = 0;
1996 }
1997 self.dual_inf_prev = inf_du;
1998 if self.dual_growth_streak >= self.dual_diverging_streak && inf_du > DUAL_DIV_FIRE_TOL {
1999 self.dual_growth_streak = 0;
2000 self.dual_inf_prev = 0.0;
2001 // Arm the "never worse off" bookkeeping for the bet about to be
2002 // placed (pounce#250 follow-up).
2003 self.dual_guard_fired = true;
2004 timing.check_convergence.end();
2005 tracing::debug!(target: "pounce::algorithm",
2006 "[POUNCE] dual-divergence guard fired at iter {} (inf_du={:.2e}); \
2007 routing to restoration (pounce#246).",
2008 self.data.borrow().iter_count, inf_du,
2009 );
2010 if self.restoration.is_some() {
2011 return self.invoke_restoration_debugged();
2012 }
2013 return IterateOutcome::Terminate(SolverReturn::DivergingIterates);
2014 }
2015 }
2016 let conv_status = self
2017 .bundle
2018 .conv_check
2019 .check_convergence_with_state(nlp_err, iter_count, &self.data, &self.cq);
2020 // Snapshot the *first* refused certificate. Baseline would have stopped
2021 // and returned exactly this point, so keeping it — and only it — is what
2022 // makes the "never worse" guarantee exact rather than approximate. A
2023 // later refusal is also a valid certificate but not necessarily a
2024 // better one, so it must not overwrite this.
2025 if !self.vetoed_seen && self.bundle.conv_check.certificate_vetoed() {
2026 // Latch on *seeing* the refusal, not on the snapshot being present:
2027 // the veto flag is sticky, so keying off `vetoed.is_none()` would
2028 // let a failed capture be completed at a later, arbitrary iterate.
2029 // See `IpoptAlgorithm::vetoed_seen`.
2030 self.vetoed_seen = true;
2031 self.vetoed = self.snapshot_current(iter_count);
2032 }
2033 if !self.vetoed_acceptable_seen && self.bundle.conv_check.acceptable_certificate_vetoed() {
2034 self.vetoed_acceptable_seen = true;
2035 self.vetoed_acceptable = self.snapshot_current(iter_count);
2036 }
2037 match conv_status {
2038 ConvergenceStatus::Continue => {}
2039 ConvergenceStatus::Converged => {
2040 timing.check_convergence.end();
2041 return IterateOutcome::Terminate(SolverReturn::Success);
2042 }
2043 ConvergenceStatus::ConvergedToAcceptable => {
2044 timing.check_convergence.end();
2045 return IterateOutcome::Terminate(SolverReturn::StopAtAcceptablePoint);
2046 }
2047 ConvergenceStatus::MaxIterExceeded => {
2048 timing.check_convergence.end();
2049 return IterateOutcome::Terminate(SolverReturn::MaxiterExceeded);
2050 }
2051 ConvergenceStatus::CpuTimeExceeded => {
2052 timing.check_convergence.end();
2053 return IterateOutcome::Terminate(SolverReturn::CpuTimeExceeded);
2054 }
2055 ConvergenceStatus::WallTimeExceeded => {
2056 timing.check_convergence.end();
2057 return IterateOutcome::Terminate(SolverReturn::WallTimeExceeded);
2058 }
2059 ConvergenceStatus::LocallyInfeasible => {
2060 timing.check_convergence.end();
2061 // gh #505: consult the acceptable-point stash, as the
2062 // restoration-cycle exits below already do (`:2686`, `:2716`,
2063 // both via `terminate_acceptable_or`). This arm used to return
2064 // without it, so a solve that had passed through an acceptable
2065 // iterate — stashed, un-vetoed, sitting there as a rollback
2066 // target — discarded it and surfaced the hard verdict instead.
2067 // The stashing code sits *after* this match, so the firing
2068 // iteration returns before it would even consider stashing;
2069 // only iterates from earlier in the solve are on offer, which
2070 // is exactly what a rollback target is.
2071 //
2072 // This is about what to *return* once the verdict has fired,
2073 // not about when it fires. Whether the rapid detector should
2074 // have convicted this point at all is a separate question,
2075 // answered by the violation floor in `OptErrorConvCheck`
2076 // (gh #519).
2077 //
2078 // Inert on genuinely infeasible models: `store_acceptable_point`
2079 // is gated on `current_is_acceptable_with_state`, which requires
2080 // `acceptable_tol` *and* the unscaled violation against
2081 // `acceptable_constr_viol_tol`, and the scale-relative veto
2082 // blocks the stash outright for a row violated relative to its
2083 // own magnitude. Nothing is stashed on such a model, so
2084 // `terminate_acceptable_or` falls through to the verdict
2085 // unchanged. `infeasible_models_are_never_reported_solved`
2086 // (`infeasible_status_tol_invariance.rs`) is the standing guard.
2087 return self.terminate_local_infeasibility();
2088 }
2089 ConvergenceStatus::Failed => {
2090 timing.check_convergence.end();
2091 return IterateOutcome::Terminate(SolverReturn::InternalError);
2092 }
2093 }
2094
2095 // Stash the iterate if it satisfies the per-component
2096 // `acceptable_*_tol` triplet. Mirrors upstream
2097 // `IpBacktrackingLineSearch.cpp:282-289` — checked at the top
2098 // of every line-search call so the most recent acceptable
2099 // iterate is always available as a rollback target if
2100 // restoration later fails. The recorder feeds
2101 // `acceptable_obj_change_tol`'s stability cross-check on
2102 // subsequent iterates.
2103 if self
2104 .bundle
2105 .conv_check
2106 .current_is_acceptable_with_state(nlp_err, &self.data, &self.cq)
2107 {
2108 self.store_acceptable_point();
2109 let curr_f = self.cq.borrow().curr_f();
2110 self.bundle.conv_check.set_curr_acceptable_obj(curr_f);
2111 // pounce#250 follow-up: keep the *best* acceptable iterate, not just
2112 // the latest. `store_acceptable_point` overwrites unconditionally,
2113 // so once the dual-divergence guard diverts a solve the rollback
2114 // target drifts to whatever the diverted run last touched — which
2115 // may be far worse than a point already in hand. Recorded on every
2116 // acceptable iterate (including before any diversion) and read only
2117 // when the guard fired; see `honour_best_acceptable_after_dual_guard`.
2118 self.record_best_acceptable(curr_f);
2119 }
2120 timing.check_convergence.end();
2121
2122 // gh #534: a deferred restoration decline is a bet with a deadline.
2123 // Checked *after* the convergence check, so a strict certificate the
2124 // continuation reached in the meantime wins the bet rather than being
2125 // pre-empted by its own expiry; and after the acceptable stash, so a
2126 // continuation that ended somewhere better has been recorded before the
2127 // floor comparison reads it.
2128 if self.decline_deadline_iter.is_some_and(|d| iter_count > d) {
2129 return self.terminate_at_decline_floor();
2130 }
2131
2132 // 3. Hessian update. Must run BEFORE `update_barrier_parameter`
2133 // so the adaptive-μ oracles (probing, quality-function) drive
2134 // their affine/centering solves against `W(curr_N)`, not the
2135 // stale `W(curr_{N-1})` left in `data.w` by the previous iter's
2136 // tail-end Hessian update. Upstream calls `UpdateHessian()`
2137 // first in every main-loop body (`IpIpoptAlg.cpp:386`); pounce
2138 // previously reordered this to the tail, which made iters 1+
2139 // pick μ from the prior iterate's Hessian on adaptive-mu +
2140 // quality-function — visible on CRESC50 as a catastrophic
2141 // early-iter divergence (theta=5.8e5 by iter 61 vs upstream
2142 // never entering restoration).
2143 timing.update_hessian.start();
2144 let _ = self.bundle.hess.update_hessian(&self.data, &self.cq);
2145 timing.update_hessian.end();
2146
2147 // 4. Barrier parameter. Pass nlp + search_dir through so the
2148 // adaptive μ oracles (probing, quality-function) can drive
2149 // their own affine-step solves; monotone ignores them.
2150 // Snapshot the tiny-step flag (set by the previous iteration's
2151 // tiny-step branch) and the entry mu — if μ can't reduce while
2152 // the flag is on, upstream `IpMonotoneMuUpdate.cpp:158-161`
2153 // throws TINY_STEP_DETECTED → STOP_AT_TINY_STEP, which we
2154 // realise as a clean termination here.
2155 //
2156 // Both updates terminate, by different routes (pounce#512).
2157 // Monotone has one throw site covering its whole update, so the
2158 // μ-unchanged comparison below reconstructs it exactly, gated on
2159 // `terminates_on_tiny_step()`. `IpAdaptiveMuUpdate.cpp` throws at
2160 // two specific sites (`:330-333`, `:377-380`) and merely fixes μ
2161 // and keeps iterating elsewhere, so the comparison would over-fire
2162 // there — on the no-bounds short-circuit, which returns before
2163 // upstream even reads the flag, and on a free-mode oracle that
2164 // re-picks the current μ. The adaptive update therefore raises
2165 // `request_tiny_step_stop` at its own two sites and opts out of
2166 // the comparison. (An earlier comment here claimed the adaptive
2167 // update never self-terminates; it does — `force_no_progress` is
2168 // what happens on the iterations that do *not* throw.)
2169 timing.update_barrier_parameter.start();
2170 let tiny_at_entry = self.data.borrow().tiny_step_flag;
2171 let mu_before = self.data.borrow().curr_mu;
2172 let mu_terminates_on_tiny = self.bundle.mu_update.terminates_on_tiny_step();
2173 let next_mu = self.bundle.mu_update.update_barrier_parameter(
2174 &self.data,
2175 &self.cq,
2176 self.nlp.as_ref(),
2177 self.search_dir.as_mut(),
2178 );
2179 self.data.borrow_mut().curr_mu = next_mu;
2180 timing.update_barrier_parameter.end();
2181
2182 // pounce#510 — line-search reset. Upstream's μ updates own a
2183 // `linesearch_` handle and call `linesearch_->Reset()` (which
2184 // clears the filter via `FilterLSAcceptor::Reset`,
2185 // `IpFilterLSAcceptor.cpp:524-532`) at four fixed points:
2186 // `IpAdaptiveMuUpdate.cpp:339` (fixed-mode decrease), `:386`
2187 // (free→fixed switch), `:431` (**unconditionally** on every
2188 // free-mode iteration, μ moved or not), and
2189 // `IpMonotoneMuUpdate.cpp:165` (after a monotone reduction).
2190 // Pounce's `MuUpdate` trait has no line-search handle, so each
2191 // update raises `request_ls_reset` at exactly those points and
2192 // we honour it here — the same plumbing `request_resto` uses
2193 // below.
2194 //
2195 // This used to be inferred from `next_mu != mu_before`. That
2196 // proxy is right for the monotone update but wrong for the
2197 // adaptive one, which resets every free-mode iteration
2198 // regardless of μ: whenever μ stayed numerically put (the
2199 // free-mode endgame, and any iteration after a restoration that
2200 // returns at the same μ) the filter kept entries computed
2201 // against a barrier parameter and an iterate the algorithm had
2202 // already left. On #505's reproducer that rejected every trial
2203 // step from α=2.4e-6 down to 1e-12 on the filter alone and
2204 // forced a spurious restoration.
2205 //
2206 // Both flags are consumed here, but the tiny-step stop is
2207 // answered first: at each of the two adaptive sites that raise
2208 // it, upstream's `TINY_STEP_DETECTED` throw sits *above* the
2209 // reset it would otherwise reach (`cpp:330-333` before `:339`,
2210 // `:377-380` before `:386`), so a terminating iteration never
2211 // resets the line search.
2212 let (tiny_step_stop_requested, ls_reset) = {
2213 let mut d = self.data.borrow_mut();
2214 let flags = (d.request_tiny_step_stop, d.request_ls_reset);
2215 d.request_tiny_step_stop = false;
2216 d.request_ls_reset = false;
2217 flags
2218 };
2219 if tiny_step_stop_requested
2220 || (tiny_at_entry
2221 && mu_terminates_on_tiny
2222 && (next_mu - mu_before).abs() < Number::EPSILON)
2223 {
2224 return IterateOutcome::Terminate(SolverReturn::StopAtTinyStep);
2225 }
2226 if ls_reset {
2227 self.bundle.line_search.reset();
2228 }
2229
2230 // pounce#58 — iterate-quality guard for the probing oracle.
2231 // The μ-update layer sets `request_resto` when the input
2232 // iterate is too corrupted for the probing rule to produce a
2233 // sane μ (see `mu/adaptive.rs` Probing dispatch). Restoration
2234 // re-initialises the multipliers and gives the outer loop a
2235 // clean iterate to continue from. When no restoration phase
2236 // is configured (embedded callers, tests), emit a one-line
2237 // notice and continue with the current μ — the guard has
2238 // already prevented the destabilising 4-order μ jump.
2239 let request_resto = {
2240 let mut d = self.data.borrow_mut();
2241 let f = d.request_resto;
2242 d.request_resto = false;
2243 f
2244 };
2245 if request_resto {
2246 if self.restoration.is_some() {
2247 return self.invoke_restoration_debugged();
2248 } else {
2249 tracing::warn!(target: "pounce::algorithm",
2250 "[POUNCE] probing-oracle iterate-quality guard fired \
2251 at iter {}, but no restoration phase is configured; \
2252 continuing with μ={:.3e}.",
2253 self.data.borrow().iter_count,
2254 next_mu,
2255 );
2256 }
2257 }
2258
2259 // Sub-iteration checkpoint: μ has been updated for this iteration.
2260 if let Some(o) = self.debug_stop(crate::debug::Checkpoint::AfterBarrierUpdate) {
2261 return o;
2262 }
2263
2264 // 5. Search direction. Skipped without an NLP + search_dir.
2265 // (Hessian was updated in step 3 above before the barrier-μ
2266 // oracle so that adaptive-μ uses W(curr_N), not stale W.)
2267 if let (Some(nlp), Some(sd)) = (self.nlp.as_ref(), self.search_dir.as_mut()) {
2268 timing.compute_search_direction.start();
2269 // Fields are declared `Empty` and filled by the linear
2270 // solver (matrix size, factor nnz, inertia, ordering — see
2271 // `pounce_feral::record_factor_stats`) and below
2272 // (regularization), so the `linear_solve` span carries the
2273 // KKT-solve characteristics for the JSON sink (pounce#71).
2274 let ls_span = tracing::info_span!(
2275 target: "pounce::linsol",
2276 "linear_solve",
2277 n = tracing::field::Empty,
2278 matrix_nnz = tracing::field::Empty,
2279 factor_nnz = tracing::field::Empty,
2280 inertia_neg = tracing::field::Empty,
2281 fill_ratio = tracing::field::Empty,
2282 ordering = tracing::field::Empty,
2283 regularization = tracing::field::Empty,
2284 );
2285 let ls_enter = ls_span.enter();
2286 let ok = sd.compute_search_direction(&self.data, &self.cq, nlp);
2287 ls_span.record("regularization", self.data.borrow().info_regu_x);
2288 // Within-span marker so the enriched `linear_solve` fields
2289 // (filled by the solver above) surface to the JSON sink at
2290 // debug level; off at the default `info` level.
2291 tracing::debug!(target: "pounce::linsol", "kkt solve complete");
2292 drop(ls_enter);
2293 timing.compute_search_direction.end();
2294 // Fine-grained time-budget gate (pounce#244). The KKT solve now
2295 // checks the shared deadline *between* its major factorization
2296 // steps (inertia correction / iterative refinement) and aborts
2297 // cooperatively when the budget is crossed — bounding the
2298 // overshoot to roughly one factorization instead of the whole
2299 // multi-factorization sweep that #242's post-solve check let run
2300 // to completion. Whether the solve returned a completed step or
2301 // bailed mid-escalation, if the deadline tripped, stop here with
2302 // the time-limit status *before* the `!ok` branch below would
2303 // otherwise route a deadline-aborted solve into restoration.
2304 // `data.curr` is untouched by the step computation, so it still
2305 // holds the last accepted iterate.
2306 if let Some(ret) = self.deadline_status() {
2307 return IterateOutcome::Terminate(ret);
2308 }
2309 if !ok {
2310 // Mirror upstream `IpIpoptAlg.cpp:417-430`: a failed
2311 // step computation puts the algorithm in emergency
2312 // mode, which calls `BacktrackingLineSearch::
2313 // ActivateFallbackMechanism` (cpp:1312-1328). When a
2314 // restoration phase is configured, the next pass of
2315 // `ComputeAcceptableTrialPoint` sees `goto_resto` at
2316 // cpp:299-306 and hands control to restoration. Only
2317 // when neither restoration nor an acceptor-level
2318 // fallback is available does upstream throw
2319 // `STEP_COMPUTATION_FAILED`.
2320 if self.restoration.is_some() {
2321 return self.invoke_restoration_debugged();
2322 }
2323 return IterateOutcome::Terminate(SolverReturn::ErrorInStepComputation);
2324 }
2325 if std::env::var_os("POUNCE_DBG_DELTA").is_some() {
2326 let d = self.data.borrow();
2327 let it = d.iter_count;
2328 if let Some(delta) = d.delta.as_ref() {
2329 use crate::iterates_vector::IteratesVector;
2330 use pounce_linalg::{Vector, compound_vector::CompoundVector};
2331 let dv: &IteratesVector = delta;
2332 tracing::debug!(target: "pounce::algorithm",
2333 "[PN_DELTA] iter={} mu={:.6e} dx_amax={:.6e} ds_amax={:.6e} dyc_amax={:.6e} dyd_amax={:.6e} dzL_amax={:.6e} dzU_amax={:.6e} dvL_amax={:.6e} dvU_amax={:.6e}",
2334 it, d.curr_mu,
2335 dv.x.amax(), dv.s.amax(), dv.y_c.amax(), dv.y_d.amax(),
2336 dv.z_l.amax(), dv.z_u.amax(), dv.v_l.amax(), dv.v_u.amax()
2337 );
2338 if let Some(cdx) = dv.x.as_any().downcast_ref::<CompoundVector>() {
2339 tracing::debug!(target: "pounce::algorithm",
2340 "[PN_DELTA] iter={} dx_blocks_amax: orig={:.6e} nc={:.6e} pc={:.6e} nd={:.6e} pd={:.6e}",
2341 it,
2342 cdx.comp(0).amax(),
2343 cdx.comp(1).amax(),
2344 cdx.comp(2).amax(),
2345 cdx.comp(3).amax(),
2346 cdx.comp(4).amax(),
2347 );
2348 tracing::debug!(target: "pounce::algorithm",
2349 "[PN_DELTA] iter={} dx_blocks_nrm2: orig={:.6e} nc={:.6e} pc={:.6e} nd={:.6e} pd={:.6e}",
2350 it,
2351 cdx.comp(0).nrm2(),
2352 cdx.comp(1).nrm2(),
2353 cdx.comp(2).nrm2(),
2354 cdx.comp(3).nrm2(),
2355 cdx.comp(4).nrm2(),
2356 );
2357 tracing::debug!(target: "pounce::algorithm",
2358 "[PN_DELTA] iter={} dx_blocks_asum: orig={:.6e} nc={:.6e} pc={:.6e} nd={:.6e} pd={:.6e}",
2359 it,
2360 cdx.comp(0).asum(),
2361 cdx.comp(1).asum(),
2362 cdx.comp(2).asum(),
2363 cdx.comp(3).asum(),
2364 cdx.comp(4).asum(),
2365 );
2366 // Argmax of orig block via dot with sign — print first few values.
2367 if let Some(dv_orig) =
2368 cdx.comp(0)
2369 .as_any()
2370 .downcast_ref::<pounce_linalg::dense_vector::DenseVector>()
2371 {
2372 let v = dv_orig.values();
2373 let mut imax = 0usize;
2374 let mut amax = 0.0f64;
2375 for (i, &x) in v.iter().enumerate() {
2376 if x.abs() > amax {
2377 amax = x.abs();
2378 imax = i;
2379 }
2380 }
2381 tracing::debug!(target: "pounce::algorithm",
2382 "[PN_DELTA] iter={} dx_orig argmax: i={} v={:.17e} (n={})",
2383 it,
2384 imax,
2385 v[imax],
2386 v.len()
2387 );
2388 }
2389 }
2390 let p = &d.perturbations;
2391 tracing::debug!(target: "pounce::algorithm",
2392 "[PN_DELTA] iter={} pert: dx={:.6e} ds={:.6e} dc={:.6e} dd={:.6e}",
2393 it, p.delta_x, p.delta_s, p.delta_c, p.delta_d
2394 );
2395 drop(d);
2396 let cq = self.cq.borrow();
2397 let gf = cq.curr_grad_f();
2398 let gl = cq.curr_grad_lag_x();
2399 let cc = cq.curr_c();
2400 let cd = cq.curr_d_minus_s();
2401 let sx = cq.curr_sigma_x();
2402 let ss = cq.curr_sigma_s();
2403 tracing::debug!(target: "pounce::algorithm",
2404 "[PN_DELTA] iter={} cq: gradf_amax={:.6e} gradf_nrm2={:.6e} gradlag_amax={:.6e} gradlag_nrm2={:.6e} c_amax={:.6e} c_nrm2={:.6e} d_amax={:.6e} d_nrm2={:.6e} sigx_amax={:.6e} sigx_nrm2={:.6e} sigs_amax={:.6e} sigs_nrm2={:.6e}",
2405 it,
2406 gf.amax(), gf.nrm2(),
2407 gl.amax(), gl.nrm2(),
2408 cc.amax(), cc.nrm2(),
2409 cd.amax(), cd.nrm2(),
2410 sx.amax(), sx.nrm2(),
2411 ss.amax(), ss.nrm2(),
2412 );
2413 if let Some(cgf) = gf.as_any().downcast_ref::<CompoundVector>() {
2414 tracing::debug!(target: "pounce::algorithm",
2415 "[PN_DELTA] iter={} gradf_blocks_amax: orig={:.6e} nc={:.6e} pc={:.6e} nd={:.6e} pd={:.6e}",
2416 it,
2417 cgf.comp(0).amax(),
2418 cgf.comp(1).amax(),
2419 cgf.comp(2).amax(),
2420 cgf.comp(3).amax(),
2421 cgf.comp(4).amax(),
2422 );
2423 }
2424 if let Some(curr) = self.data.borrow().curr.clone() {
2425 tracing::debug!(target: "pounce::algorithm",
2426 "[PN_DELTA] iter={} bound_mults: zL_amax={:.6e} zU_amax={:.6e} vL_amax={:.6e} vU_amax={:.6e} s_amax={:.6e} s_nrm2={:.6e} x_amax={:.6e} x_nrm2={:.6e}",
2427 it,
2428 curr.z_l.amax(), curr.z_u.amax(),
2429 curr.v_l.amax(), curr.v_u.amax(),
2430 curr.s.amax(), curr.s.nrm2(),
2431 curr.x.amax(), curr.x.nrm2(),
2432 );
2433 if let Some(czl) = curr.z_l.as_any().downcast_ref::<CompoundVector>() {
2434 tracing::debug!(target: "pounce::algorithm",
2435 "[PN_DELTA] iter={} zL_blocks_amax: orig={:.6e} nc={:.6e} pc={:.6e} nd={:.6e} pd={:.6e}",
2436 it,
2437 czl.comp(0).amax(),
2438 czl.comp(1).amax(),
2439 czl.comp(2).amax(),
2440 czl.comp(3).amax(),
2441 czl.comp(4).amax(),
2442 );
2443 }
2444 if let Some(czu) = curr.z_u.as_any().downcast_ref::<CompoundVector>() {
2445 tracing::debug!(target: "pounce::algorithm", "[PN_DELTA] iter={} zU_ncomps={}", it, czu.n_comps());
2446 for ic in 0..czu.n_comps() {
2447 tracing::debug!(target: "pounce::algorithm",
2448 "[PN_DELTA] iter={} zU_block[{}]_amax={:.6e} dim={}",
2449 it,
2450 ic,
2451 czu.comp(ic).amax(),
2452 czu.comp(ic).dim()
2453 );
2454 }
2455 }
2456 }
2457 if let Some(csx) = sx.as_any().downcast_ref::<CompoundVector>() {
2458 tracing::debug!(target: "pounce::algorithm",
2459 "[PN_DELTA] iter={} sigx_blocks_amax: orig={:.6e} nc={:.6e} pc={:.6e} nd={:.6e} pd={:.6e}",
2460 it,
2461 csx.comp(0).amax(),
2462 csx.comp(1).amax(),
2463 csx.comp(2).amax(),
2464 csx.comp(3).amax(),
2465 csx.comp(4).amax(),
2466 );
2467 }
2468 drop(cq);
2469 let d = self.data.borrow();
2470 // Also dump curr.x_orig argmax
2471 if let Some(curr) = d.curr.as_ref() {
2472 if let Some(cx) = curr.x.as_any().downcast_ref::<CompoundVector>() {
2473 if let Some(xo) =
2474 cx.comp(0)
2475 .as_any()
2476 .downcast_ref::<pounce_linalg::dense_vector::DenseVector>()
2477 {
2478 let v = xo.values();
2479 let mut imax = 0usize;
2480 let mut amax = 0.0f64;
2481 for (i, &x) in v.iter().enumerate() {
2482 if x.abs() > amax {
2483 amax = x.abs();
2484 imax = i;
2485 }
2486 }
2487 tracing::debug!(target: "pounce::algorithm", "[PN_DELTA] iter={} curr_x_orig argmax: i={} v={:.17e} amax={:.17e} nrm2={:.17e}",
2488 it, imax, v[imax], xo.amax(), xo.nrm2());
2489 }
2490 }
2491 }
2492 }
2493 }
2494 }
2495
2496 // Capture KKT-factorization diagnostics for the debugger before
2497 // the line search runs. Only when a debugger is installed. The
2498 // inertia/status fields are cheap and always captured; the matrix
2499 // triplets and `LDLᵀ` factor are O(nnz) assemblies, so they're
2500 // captured only while the debugger is stepping (`wants_kkt_capture`)
2501 // — a detached/free-running debugger drops them to keep the run
2502 // cheap. `kkt_debug` is overwritten every iteration and never
2503 // cleared at `iter_start`, so a stepping session always has the
2504 // previous iteration's system to look back at via `viz kkt`/`viz L`.
2505 if let Some(hook) = self.debug.as_ref() {
2506 let capture_heavy = hook.borrow().wants_kkt_capture();
2507 let captured_iter = self.data.borrow().iter_count;
2508 let info = self.search_dir.as_ref().map(|sd| {
2509 let pd = sd.pd_solver_mut();
2510 let aug = pd.aug_solver();
2511 let provides = aug.provides_inertia();
2512 crate::ipopt_data::KktDebug {
2513 iter: captured_iter,
2514 dim: aug.system_dim(),
2515 n_neg: if provides {
2516 aug.number_of_neg_evals()
2517 } else {
2518 -1
2519 },
2520 provides_inertia: provides,
2521 status: format!("{:?}", aug.last_solve_status()),
2522 matrix: if capture_heavy {
2523 aug.kkt_triplets()
2524 } else {
2525 None
2526 },
2527 l_factor: if capture_heavy {
2528 aug.l_factor(true)
2529 } else {
2530 None
2531 },
2532 }
2533 });
2534 self.data.borrow_mut().kkt_debug = info;
2535 }
2536
2537 // Sub-iteration checkpoint: the Newton step `δ` (data.delta) and
2538 // the applied regularization are now available, before the line
2539 // search consumes them.
2540 if let Some(o) = self.debug_stop(crate::debug::Checkpoint::AfterSearchDirection) {
2541 return o;
2542 }
2543
2544 // Fine-grained time-budget gate (pounce#242). The KKT
2545 // factorization (and any inertia-correction / quality-escalation
2546 // refactorizations) is the single most expensive step of a large
2547 // solve, and it has just finished. Check the deadline here so an
2548 // over-budget solve returns its current best iterate *before*
2549 // spending a line search and another whole iteration — bounding
2550 // the overshoot to roughly one search-direction computation
2551 // instead of a full outer iteration. `data.curr` is untouched by
2552 // the step computation (the trial lives in `data.trial`), so it
2553 // still holds the last accepted iterate.
2554 if let Some(ret) = self.deadline_status() {
2555 return IterateOutcome::Terminate(ret);
2556 }
2557
2558 // 6. Acceptable trial point — run the line search if we have a
2559 // primal/dual step on `data.delta`. Wrap in a guard so all
2560 // early-return paths (ErrorInStepComputation, InternalError,
2561 // restoration entry) still stop the timer.
2562 let _ls_guard = timing.compute_acceptable_trial_point.guard();
2563 let have_delta = self.data.borrow().delta.is_some();
2564 if have_delta {
2565 let delta = match self.data.borrow().delta.as_ref().cloned() {
2566 Some(d) => d,
2567 None => {
2568 return IterateOutcome::Terminate(SolverReturn::ErrorInStepComputation);
2569 }
2570 };
2571 // Cap alpha by the primal fraction-to-the-boundary so the
2572 // first trial cannot push slacks past their bounds, and by
2573 // the dual FTB so bound multipliers stay positive. Mirrors
2574 // upstream `IpBacktrackingLineSearch::FindAcceptableTrialPoint`'s
2575 // calls to `IpCq.primal_frac_to_the_bound` /
2576 // `IpCq.dual_frac_to_the_bound` with τ = `curr_tau`.
2577 let tau = self.data.borrow().curr_tau;
2578 let alpha_p_max = self.cq.borrow().aff_step_alpha_primal_max(&delta, tau);
2579 let alpha_d_max = self.cq.borrow().aff_step_alpha_dual_max(&delta, tau);
2580
2581 // Tiny-step gate — port of `IpBacktrackingLineSearch.cpp:363`
2582 // and the handling block at lines 382-435. When the search
2583 // direction is so small that any nonzero α would just
2584 // bounce inside floating-point noise, take the FTB step
2585 // unchecked and skip the line search; that's the only way
2586 // to hit `STOP_AT_TINY_STEP` cleanly when the iterate is
2587 // already at a converged point but `nlp_error > tol` due to
2588 // scaling or unbounded duals.
2589 if self.detect_tiny_step(&delta) {
2590 let alpha_p = alpha_p_max;
2591 let alpha_d = alpha_d_max;
2592 let curr = match self.data.borrow().curr.clone() {
2593 Some(c) => c,
2594 None => return IterateOutcome::Terminate(SolverReturn::InternalError),
2595 };
2596 let trial_iv = scaled_step_unchecked(&curr, &delta, alpha_p, alpha_d);
2597 {
2598 let mut d = self.data.borrow_mut();
2599 d.set_trial(trial_iv);
2600 d.info_alpha_primal = alpha_p;
2601 d.info_alpha_dual = alpha_d;
2602 d.info_ls_count = 0;
2603 if self.tiny_step_last_iteration {
2604 d.info_alpha_primal_char = 'T';
2605 d.tiny_step_flag = true;
2606 } else {
2607 d.info_alpha_primal_char = 't';
2608 }
2609 }
2610 let dy_amax = delta.y_c.amax().max(delta.y_d.amax());
2611 self.tiny_step_last_iteration = dy_amax < self.tiny_step_y_tol;
2612 } else {
2613 self.tiny_step_last_iteration = false;
2614 let alpha_init = self.alpha_init.min(alpha_p_max);
2615 let alpha_dual = self.alpha_init.min(alpha_d_max);
2616 let outcome = self.bundle.line_search.find_acceptable_trial_point(
2617 &self.data,
2618 &self.cq,
2619 &delta,
2620 alpha_init,
2621 alpha_dual,
2622 self.nlp.as_ref(),
2623 self.search_dir.as_mut(),
2624 );
2625 match outcome {
2626 Outcome::Accepted => {
2627 // A normal LS-accepted step breaks any in-flight
2628 // restoration cycle — clear the cycle detector
2629 // so the next resto entry starts fresh.
2630 self.last_resto_entry_x = None;
2631 self.last_resto_entry_s = None;
2632 self.last_resto_recovery_x = None;
2633 self.last_resto_recovery_s = None;
2634 self.resto_no_outer_progress_count = 0;
2635 // Intentionally *not* clearing
2636 // `resto_near_feasible_count` here: DECONVBNE's
2637 // cycle interleaves R-recoveries with 2-3
2638 // LS-accepted 'f'/'h' steps (which return
2639 // `Outcome::Accepted` but accomplish no real
2640 // outer progress — alpha drops to 1e-6 and
2641 // inf_du remains pinned at 1.9e7), so resetting
2642 // on every accept would zero the counter every
2643 // cycle and never fire. The counter persists
2644 // for the duration of the run and trips after
2645 // 3 cumulative near-feasible entries; legitimate
2646 // solves enter resto at most once at near-
2647 // feasibility (POLAK6, HAIFAM) and stay under
2648 // the limit.
2649 }
2650 Outcome::TinyStep | Outcome::Failed => {
2651 // Debugger stop: the line search rejected the step
2652 // (tiny-step floor or all backtracks failed), before
2653 // we fall into restoration. Lets a "why did the line
2654 // search give up?" inspection happen at the failing
2655 // point distinctly from the restoration entry.
2656 if let Some(o) = self.debug_stop(crate::debug::Checkpoint::StepRejected) {
2657 return o;
2658 }
2659 // Upstream `IpBacktrackingLineSearch.cpp` raises
2660 // `LINE_SEARCH_FAILED` when α drops below
2661 // `alpha_min` or all retries reject, which in
2662 // turn triggers `ActivateLineSearch` →
2663 // restoration.
2664 return self.invoke_restoration_debugged();
2665 }
2666 Outcome::Deadline => {
2667 // The time budget was crossed inside the line
2668 // search (pounce#242). No trial was promoted, so
2669 // `data.curr` still holds the best iterate; stop
2670 // with the matching time-limit status. Re-derive
2671 // wall vs CPU from the deadline (it can only still
2672 // be exceeded — time is monotonic).
2673 return IterateOutcome::Terminate(
2674 self.deadline_status()
2675 .unwrap_or(SolverReturn::WallTimeExceeded),
2676 );
2677 }
2678 }
2679 }
2680 }
2681
2682 // End the line-search/trial timer here so the bookkeeping in
2683 // steps 7-8 below is attributed to `accept_trial_point` (which
2684 // mirrors upstream's split: filter update and FTB reset are
2685 // accept-side, not line-search-side).
2686 _ls_guard.stop();
2687
2688 // 7. Accept trial point (promotes `trial` to `curr` if set).
2689 // The acceptor's filter has already been augmented (when
2690 // appropriate) inside `find_acceptable_trial_point` via
2691 // `update_for_next_iteration`, mirroring upstream's call
2692 // chain in `IpBacktrackingLineSearch.cpp:839`.
2693 let _accept_guard = timing.accept_trial_point.guard();
2694
2695 // 7a. Safe-slack bound adjustment. Before promoting `trial`, move
2696 // any `x_L/x_U/d_L/d_U` whose trial slack fell below
2697 // `eps*min(1,mu)` so the slack becomes representable (port of
2698 // the bound-adjustment block in
2699 // `IpoptAlgorithm::AcceptTrialPoint`, `IpIpoptAlg.cpp:664-706`).
2700 self.adjust_variable_bounds_for_small_slacks();
2701
2702 self.data.borrow_mut().accept_trial_point();
2703
2704 // 8. Bound multiplier kappa_sigma reset.
2705 self.correct_bound_multiplier();
2706
2707 // Sub-iteration checkpoint: the trial point was accepted; α and
2708 // the new iterate are in place (before the loop's iter bookkeeping
2709 // and the next `IterStart`).
2710 drop(_accept_guard);
2711 if let Some(o) = self.debug_stop(crate::debug::Checkpoint::AfterStep) {
2712 return o;
2713 }
2714
2715 IterateOutcome::Continue
2716 }
2717
2718 /// Port of `IpBacktrackingLineSearch::DetectTinyStep`
2719 /// (`IpBacktrackingLineSearch.cpp:1219-1278`). Returns true iff
2720 /// `max_i |δx_i|/(1+|x_i|) ≤ tiny_step_tol`,
2721 /// `max_i |δs_i|/(1+|s_i|) ≤ tiny_step_tol`, AND
2722 /// `curr_constraint_violation ≤ 1e-4`. Disabled when
2723 /// `tiny_step_tol == 0`.
2724 fn detect_tiny_step(&self, delta: &crate::iterates_vector::IteratesVector) -> bool {
2725 if self.tiny_step_tol == 0.0 {
2726 return false;
2727 }
2728 let curr = match self.data.borrow().curr.clone() {
2729 Some(c) => c,
2730 None => return false,
2731 };
2732
2733 // |x_i|+1
2734 let mut tmp = curr.x.make_new_copy();
2735 tmp.element_wise_abs();
2736 tmp.add_scalar(1.0);
2737 // |δx_i|/(|x_i|+1) ; checked via Amax of (δx ./ (|x|+1)).
2738 let mut tmp2 = delta.x.make_new_copy();
2739 tmp2.element_wise_divide(&*tmp);
2740 if tmp2.amax() > self.tiny_step_tol {
2741 return false;
2742 }
2743
2744 if curr.s.dim() > 0 {
2745 let mut tmp = curr.s.make_new_copy();
2746 tmp.element_wise_abs();
2747 tmp.add_scalar(1.0);
2748 let mut tmp2 = delta.s.make_new_copy();
2749 tmp2.element_wise_divide(&*tmp);
2750 if tmp2.amax() > self.tiny_step_tol {
2751 return false;
2752 }
2753 }
2754
2755 let cviol = self.cq.borrow().curr_constraint_violation();
2756 if cviol > 1e-4 {
2757 return false;
2758 }
2759 true
2760 }
2761
2762 /// Drive the restoration phase after a line-search failure.
2763 /// Returns `IterateOutcome::Continue` if the restoration driver
2764 /// recovered (the algorithm carries on from the recovered iterate);
2765 /// otherwise terminates with [`SolverReturn::RestorationFailure`].
2766 /// Mirrors upstream's
2767 /// `IpBacktrackingLineSearch::ActivateLineSearch` → `PerformRestoration`
2768 /// chain.
2769 fn invoke_restoration(&mut self) -> IterateOutcome {
2770 // Snapshot the outer reference iterate's `(theta, barr)` and
2771 // build the orig-progress callback the inner IPM will consult
2772 // at every iteration (mirrors upstream
2773 // `IpRestoFilterConvCheck::SetOrigLSAcceptor` plus
2774 // `IpFilterLSAcceptor::Reset`'s `reference_*_` snapshot).
2775 let reference_theta = self.cq.borrow().curr_constraint_violation();
2776 let reference_barr = self.cq.borrow().curr_barrier_obj();
2777
2778 if std::env::var("POUNCE_DBG_RESTO").is_ok() {
2779 let iter = self.data.borrow().iter_count;
2780 tracing::debug!(target: "pounce::algorithm",
2781 "RESTO_ENTRY iter={} theta={:.6e} barr={:.6e} near_feas_ct={}",
2782 iter, reference_theta, reference_barr, self.resto_near_feasible_count,
2783 );
2784 }
2785
2786 // Port gap: upstream refuses to enter restoration from an acceptable
2787 // point, and this was missing. `IpBacktrackingLineSearch.cpp:557-570`,
2788 // in the `if (!accept)` arm that hands off to restoration:
2789 //
2790 // if( CurrentIsAcceptable() )
2791 // {
2792 // THROW_EXCEPTION(ACCEPTABLE_POINT_REACHED,
2793 // "Restoration phase called at acceptable point.");
2794 // }
2795 //
2796 // The rationale is the obvious one: restoration reduces the constraint
2797 // violation, so from a point that already passes the acceptable-level
2798 // tolerances it has nothing to reduce, and entering can only risk a
2799 // reportable solution.
2800 //
2801 // What the gap cost, measured on mittelmann `qcqp1000-1nc` (n=1000):
2802 // the line search fails at iteration 187 on a point carrying the
2803 // published optimum (`-2.6628866e+07`, matching ipopt-ma57 to 9
2804 // significant figures) with overall NLP error `6.0e-8` — two orders
2805 // inside `acceptable_tol`. Restoration walked it to `theta 5e-3` and
2806 // ground out 2780 further iterations without recovering, so a solved
2807 // problem reported a failure.
2808 //
2809 // The predicate is upstream's, unmodified: acceptability alone. A
2810 // strict `constr_viol_tol` gate was tried on top and is both a
2811 // deviation and useless — at their restoration entries `qcqp1000-1nc`
2812 // sits at `theta = 6.0e-8`, `csfi2` at `1.5e-7`, `eigena2` at
2813 // `2.1e-10`, all strictly feasible, one by six orders. Nothing
2814 // observable at the doorway separates a restoration that recovers from
2815 // one that does not, which is why upstream does not try to.
2816 //
2817 // Placed ahead of the cycle detectors below rather than beside
2818 // upstream's `PrepareRestoPhaseStart()`: those detectors are a
2819 // pounce-side addition, and an acceptable point should be reported
2820 // regardless of cycle state. Filter augmentation is skipped on this
2821 // path, which is immaterial — the run stops here.
2822 //
2823 // `current_is_acceptable_with_state` is the full triplet, never
2824 // `theta` alone: gh #274, a perfectly feasible point can be
2825 // arbitrarily far from stationary (`min -exp(x) s.t. x >= 0` reaches
2826 // here with `inf_pr = 1.7e-10` and `inf_du = 8.8e+47`), and the
2827 // triplet carries `acceptable_dual_inf_tol` to reject it. The
2828 // finiteness check mirrors the one below (CUTE `himmelbj` reaches a
2829 // near-feasible point where `f` evaluates to NaN) and matches
2830 // upstream's own `curr_f` precondition for acceptability.
2831 let (entry_f_finite, entry_nlp_err) = {
2832 let cq = self.cq.borrow();
2833 (cq.curr_f().is_finite(), cq.curr_nlp_error())
2834 };
2835 //
2836 // What the guard still did not ask is whether the solve was *converging*
2837 // (gh #534). It reads the entry point and nothing about the trajectory
2838 // that reached it, so it stops a contracting endgame and a dead stall
2839 // with equal confidence. On `eigena2` it fires while the dual
2840 // infeasibility is quartering every iteration on unit steps
2841 // (`1.19e-5 → 2.96e-6 → 7.38e-7 → 1.84e-7`), three iterations short of a
2842 // strict certificate that costs nothing but those three iterations.
2843 // [`Self::may_defer_acceptable_decline`] adds that missing question,
2844 // and only that: when the answer is no — `eigenb2`'s tail rises, and
2845 // `csfi2`'s last two iterations are flat to three digits — the guard
2846 // fires exactly as before.
2847 if entry_f_finite
2848 && self.bundle.conv_check.current_is_acceptable_with_state(
2849 entry_nlp_err,
2850 &self.data,
2851 &self.cq,
2852 )
2853 {
2854 if self.may_defer_acceptable_decline() {
2855 tracing::debug!(target: "pounce::algorithm",
2856 "[POUNCE] deferring the restoration decline at theta {:.3e}: the entry \
2857 point passes the acceptable-level tolerances (nlp_err {:.3e}) but the \
2858 NLP error has contracted every iteration over the last {} \
2859 ({:.3e} -> {:.3e}); continuing for up to {} iterations, with that point \
2860 held as the floor (gh #534).",
2861 reference_theta, entry_nlp_err, DECLINE_PROGRESS_SAMPLES - 1,
2862 self.nlp_err_recent[0], entry_nlp_err, DECLINE_CONTINUATION_BUDGET,
2863 );
2864 } else {
2865 // The window is on the line because "why did the guard not
2866 // defer?" is the first question anyone reading this trace has
2867 // (gh #534), and reconstructing it from the iteration table
2868 // means recomputing the scaled aggregate by hand.
2869 tracing::debug!(target: "pounce::algorithm",
2870 "[POUNCE] declining restoration at theta {:.3e}: the entry point already \
2871 passes the acceptable-level tolerances (nlp_err {:.3e}); reporting it \
2872 rather than risking it in restoration. Recent NLP errors {} \
2873 (contracting: {}).",
2874 reference_theta, entry_nlp_err,
2875 self.nlp_err_window_str(), self.nlp_err_contracting(),
2876 );
2877 return IterateOutcome::Terminate(SolverReturn::StopAtAcceptablePoint);
2878 }
2879 }
2880
2881 // No-progress restoration cycle detector. Two layered checks
2882 // surface as `ErrorInStepComputation` instead of cycling to
2883 // `max_iter` exhaustion (mirrors the *intent* of upstream
2884 // `IpBacktrackingLineSearch.cpp:580-600`'s almost-feasible
2885 // resto guard):
2886 //
2887 // 1. *Static cycle*: entry-to-entry — when the curr `(x, s)`
2888 // at this entry is essentially identical to the snapshot
2889 // from the previous entry, the inner resto-IPM is
2890 // returning recovered iterates indistinguishable from
2891 // entry, AND the outer didn't move either. Fires
2892 // immediately. Catches QCNEW, EQC, MESH, POLAK6, S365,
2893 // S365MOD, SIPOW2M, PFIT4.
2894 //
2895 // 2. *Slow-progress cycle*: recovery-to-entry — when curr at
2896 // this entry is essentially identical to the *recovery*
2897 // iterate from the previous resto, the outer made no
2898 // progress between resto invocations even though resto's
2899 // inner moved substantively. Counted, fires after 5
2900 // consecutive entries. Catches ACOPR14, ACOPR30, TRO3X3
2901 // while letting MAKELA3, HAIFAM, HALDMADS, ROBOT,
2902 // TENBARS2 — which need 2-3 productive resto entries
2903 // before LS accepts — pass through.
2904 //
2905 // A productive single-restoration sequence (BT8, HIMMELBJ,
2906 // LINSPANH, LSNNODOC, ODFITS, OET3) clears both snapshots via
2907 // `Outcome::Accepted` between entries and is unaffected.
2908 let curr = self
2909 .data
2910 .borrow()
2911 .curr
2912 .as_ref()
2913 .expect("curr set before invoke_restoration")
2914 .clone();
2915 // Helper: when the cycle detector fires and the orig cv is a
2916 // violation the *user* calls a violation (e.g. PFIT1's 2.73e-2),
2917 // the outer is stuck at a feasibility-stationary point and the
2918 // honest exit is `LocalInfeasibility`. Below that threshold the
2919 // iterate is primal-feasible by the user's own declaration, so there
2920 // is no infeasibility to certify — the failure is numerical, not
2921 // algorithmic, and `ErrorInStepComputation` is retained.
2922 //
2923 // The threshold is `constr_viol_tol`, and *only* `constr_viol_tol`
2924 // (gh #508). The question this ternary asks — "is this violation
2925 // real?" — is a question about the constraint violation, so it has to
2926 // be asked with the option that declares what a violated constraint
2927 // is. The previous form, `max(100·tol, 1e-4)`, was built from `tol`, a
2928 // tolerance on the **KKT error**: different quantity, different units,
2929 // and it never consulted `constr_viol_tol` at all. Two consequences,
2930 // both measured on `min (x-5)² s.t. x²+δ = 0` (infeasible for every
2931 // δ>0, reported violation exactly δ):
2932 //
2933 // * sweeping `constr_viol_tol` over four orders moved the boundary
2934 // not at all — at `constr_viol_tol = 1e-3` a violation of `1e-4`,
2935 // comfortably inside the user's declared feasibility tolerance,
2936 // still exited 500;
2937 // * sweeping `tol` moved it a great deal, and in the wrong
2938 // direction: at `tol = 1e-4` the `1e-2` threshold swallowed every
2939 // δ from `3e-4` to `1e-2` — a model infeasible by a full percent
2940 // answered "your solver broke". Loosening `tol` is the standard
2941 // user reaction to a struggling solve, so the failure widened
2942 // exactly when the user tried to help.
2943 //
2944 // No `infeas_viol_kappa` margin on top, unlike the rapid-infeasibility
2945 // pre-filter in `conv_check`. That detector fires *during* the solve
2946 // off a streak heuristic and needs the margin to avoid convicting an
2947 // iterate that is still converging; here restoration has already
2948 // demonstrably cycled, so the certainty comes from the cycle evidence
2949 // rather than from extra violation headroom. Widening to
2950 // `kappa·constr_viol_tol` would move the default threshold from `1e-4`
2951 // to `1e-2` and hand back 500 on the whole band in between.
2952 //
2953 // The comparison is `>=`, not `>`. A violation landing exactly on the
2954 // threshold is a violation at the user's declared tolerance, and the
2955 // reproducer above hits the boundary to the digit (`δ = 1e-4` at the
2956 // default `constr_viol_tol`), where `>` returned 500 for a model
2957 // infeasible by precisely the amount the user said was too much.
2958 //
2959 // The violation is measured **unscaled**. `reference_theta` is the
2960 // row-scaled residual, but the floor below is an absolute, user-facing
2961 // magnitude, so comparing the two mixes unit systems — and on a problem
2962 // whose rows are scaled down the scaled residual can never clear it.
2963 // `infeasible_equalities.nl` is the worked example: a square 2x2 system
2964 // with a true violation of 2.0 that NLP scaling reports as 6.67e-7, so
2965 // this test read `6.67e-7 > 1e-4` = false and a blatantly infeasible
2966 // model exited `Error_In_Step_Computation` (AMPL 500, Pyomo
2967 // `internalSolverError`). Square problems have no restoration-side
2968 // locally-infeasible gate — `strict` carves them out so the outer gets
2969 // another shot — so this cycle exit *is* their safety net, and it was
2970 // disabled by the unit mismatch. Same user-visible family as gh #372.
2971 //
2972 // Note this also moves from a 1-norm (`curr_constraint_violation`) to a
2973 // max-norm. Max-norm <= 1-norm, so the test is marginally stricter
2974 // about declaring infeasibility on an unscaled problem — the safe
2975 // direction for a verdict this consequential.
2976 //
2977 // `theta > 0` in front of the `>=` is not redundant: the options layer
2978 // registers `constr_viol_tol` with a *strict* lower bound of zero, but
2979 // a library embedder setting `ConvCheckOptions` directly is not bound
2980 // by that, and `0 >= 0` would turn an exactly-feasible iterate into an
2981 // infeasibility certificate. A zero violation never proves anything.
2982 let cycle_viol_tol = self.bundle.conv_check.constr_viol_tol_or_default();
2983 let reference_theta_unscaled = self.cq.borrow().curr_unscaled_primal_infeasibility_max();
2984 let cycle_exit =
2985 if reference_theta_unscaled > 0.0 && reference_theta_unscaled >= cycle_viol_tol {
2986 SolverReturn::LocalInfeasibility
2987 } else {
2988 SolverReturn::ErrorInStepComputation
2989 };
2990 let static_cycle = if let (Some(prev_x), Some(prev_s)) = (
2991 self.last_resto_entry_x.as_ref(),
2992 self.last_resto_entry_s.as_ref(),
2993 ) {
2994 let dx_rel = relative_distance(&*curr.x, &**prev_x);
2995 let ds_rel = relative_distance(&*curr.s, &**prev_s);
2996 if std::env::var_os("POUNCE_DBG_RESTO_CYCLE").is_some() {
2997 tracing::debug!(target: "pounce::algorithm",
2998 "[PN_RESTO_CYCLE] entry-vs-entry dx_rel={:.6e} ds_rel={:.6e}",
2999 dx_rel, ds_rel
3000 );
3001 }
3002 dx_rel <= 1e-10 && ds_rel <= 1e-10
3003 } else {
3004 false
3005 };
3006 if static_cycle {
3007 // Prefer the last acceptable point over the cycle error —
3008 // the borrows above are released, so the `&mut self` helper
3009 // is free to roll back.
3010 return self.terminate_acceptable_or(cycle_exit);
3011 }
3012 let recovery_cycle = if let (Some(prev_x), Some(prev_s)) = (
3013 self.last_resto_recovery_x.as_ref(),
3014 self.last_resto_recovery_s.as_ref(),
3015 ) {
3016 let dx_rel = relative_distance(&*curr.x, &**prev_x);
3017 let ds_rel = relative_distance(&*curr.s, &**prev_s);
3018 if std::env::var_os("POUNCE_DBG_RESTO_CYCLE").is_some() {
3019 tracing::debug!(target: "pounce::algorithm",
3020 "[PN_RESTO_CYCLE] entry-vs-recovery dx_rel={:.6e} ds_rel={:.6e} count={}",
3021 dx_rel, ds_rel, self.resto_no_outer_progress_count
3022 );
3023 }
3024 dx_rel <= 1e-10 && ds_rel <= 1e-10
3025 } else {
3026 false
3027 };
3028 if recovery_cycle {
3029 self.resto_no_outer_progress_count =
3030 self.resto_no_outer_progress_count.saturating_add(1);
3031 // 10-strike limit: tuned to give OET7-style traces room
3032 // to break through (inner inf_pr still decreasing across
3033 // strikes) while still bounding DECONVBNE-style cycles
3034 // (which need a guard but tolerate a wider window —
3035 // ~3 outer steps per cycle, so 10 strikes ≈ 30 outer
3036 // iters, well below the 2987-iter pathological run).
3037 if self.resto_no_outer_progress_count >= 10 {
3038 // Prefer the last acceptable point over the cycle error;
3039 // borrows are released, so the `&mut self` helper is free.
3040 return self.terminate_acceptable_or(cycle_exit);
3041 }
3042 } else {
3043 self.resto_no_outer_progress_count = 0;
3044 }
3045 // Near-feasible resto re-entry detector — matches the *intent*
3046 // of upstream `IpBacktrackingLineSearch.cpp:580-600`'s almost-
3047 // feasible-resto guard with a looser cv threshold. When the
3048 // outer enters restoration with the constraint violation
3049 // already at or below `tol`, the resto sub-IPM will produce a
3050 // recovered iterate that's at most marginally more feasible,
3051 // and any post-recovery σ-blowup from the next outer KKT solve
3052 // will re-trigger resto on the next iteration. Counting these
3053 // entries surfaces the cycle as `StopAtAcceptablePoint` —
3054 // primal feasibility is already met, only the dual residual
3055 // remains. Catches DECONVBNE: pounce ran 2987 iters before
3056 // this guard (cycle of ~30-inner-resto + 3 outer per cycle);
3057 // upstream solves in 505 iters via a different x trajectory.
3058 // Single-entry productive restos (BT8, HIMMELBJ, ODFITS) and
3059 // sub-tol-but-recoverable starts pass through under the 3-
3060 // strike limit.
3061 let outer_tol = self.bundle.conv_check.tol_or_default();
3062 if reference_theta <= outer_tol {
3063 self.resto_near_feasible_count = self.resto_near_feasible_count.saturating_add(1);
3064 if self.resto_near_feasible_count >= 3 {
3065 // Constraint feasibility is met, but a near-feasible iterate is
3066 // only "acceptable" if its objective is finite. CUTE `himmelbj`
3067 // reaches a point with cv ≈ 2e-9 where f evaluates to NaN; that
3068 // must surface as Invalid_Number_Detected rather than be
3069 // reported as Solved_To_Acceptable_Level with a `nan` objective.
3070 if !self.cq.borrow().curr_f().is_finite() {
3071 return IterateOutcome::Terminate(SolverReturn::InvalidNumberDetected);
3072 }
3073 // Constraint feasibility alone does not make a point
3074 // acceptable. `reference_theta` measures only the *primal*
3075 // residual, so a perfectly feasible iterate can still be
3076 // arbitrarily far from stationary — which is exactly what an
3077 // unbounded objective looks like from here: the constraints
3078 // stay satisfied while the iterates run off toward -inf.
3079 //
3080 // `min -exp(x) s.t. x >= 0` re-enters restoration with
3081 // `inf_pr = 1.7e-10` and `inf_du = 8.8e+47`; before gh #274
3082 // the finiteness check was the only gate, `-8.8e47` is
3083 // finite, and the solve was reported as
3084 // `Solved_To_Acceptable_Level` with `solve_result_num = 100`.
3085 // Pyomo maps that into the *solved* family and loads the
3086 // diverging iterate as an optimal solution.
3087 //
3088 // So require the point to pass the full acceptable-level
3089 // triplet (which includes `acceptable_dual_inf_tol`) before
3090 // claiming acceptability. When it does not, surface
3091 // `cycle_exit` — the same honest status the other two
3092 // restoration-cycle exits in this function use.
3093 let nlp_err = self.cq.borrow().curr_nlp_error();
3094 if !self
3095 .bundle
3096 .conv_check
3097 .current_is_acceptable_with_state(nlp_err, &self.data, &self.cq)
3098 {
3099 tracing::debug!(target: "pounce::algorithm",
3100 "[POUNCE] near-feasible restoration re-entry at theta {:.3e} \
3101 but the point fails the acceptable-level tolerances \
3102 (nlp_err {:.3e}); reporting {:?} rather than \
3103 Solved_To_Acceptable_Level (gh#274).",
3104 reference_theta, nlp_err, cycle_exit,
3105 );
3106 return IterateOutcome::Terminate(cycle_exit);
3107 }
3108 return IterateOutcome::Terminate(SolverReturn::StopAtAcceptablePoint);
3109 }
3110 } else {
3111 self.resto_near_feasible_count = 0;
3112 }
3113 self.last_resto_entry_x = Some(curr.x.make_new_copy());
3114 self.last_resto_entry_s = Some(curr.s.make_new_copy());
3115
3116 // Augment the outer's filter with the resto-entry envelope —
3117 // mirrors upstream `IpBacktrackingLineSearch.cpp:566`:
3118 // `acceptor_->PrepareRestoPhaseStart()`. Adds
3119 // `((1-γ_θ)·θ_entry, φ_entry - γ_φ·θ_entry)` to the filter so
3120 // that after restoration recovers, the outer's Newton step is
3121 // forced by the filter to make real progress vs the entry
3122 // point. Without this, the outer accepts null-progress 'h'
3123 // steps and re-enters restoration on the next iteration (root
3124 // cause of DECONVBNE's 323 R-accepts vs ipopt's 21).
3125 self.bundle
3126 .line_search
3127 .acceptor_mut()
3128 .prepare_resto_phase_start(reference_theta, reference_barr);
3129
3130 let orig_progress_cb = self.bundle.line_search.acceptor().make_orig_progress_check(
3131 reference_theta,
3132 reference_barr,
3133 5.0,
3134 );
3135
3136 let (Some(nlp), Some(sd), Some(resto)) = (
3137 self.nlp.as_ref(),
3138 self.search_dir.as_mut(),
3139 self.restoration.as_mut(),
3140 ) else {
3141 return IterateOutcome::Terminate(SolverReturn::RestorationFailure);
3142 };
3143 resto.set_orig_progress_check(orig_progress_cb);
3144 // Forward the shared debugger so it can step the inner solve.
3145 resto.set_debug_hook(self.debug.as_ref().map(Rc::clone));
3146 let mut pd_guard = sd.pd_solver_mut();
3147 let aug = pd_guard.aug_solver_mut();
3148 // Audit counters (pounce#12). Increment call count + outer-iter
3149 // count (one outer iter is consumed per restoration call) and
3150 // wall-time around the inner call. Inner iter count is read
3151 // after via the trait accessor.
3152 self.resto_calls = self.resto_calls.saturating_add(1);
3153 self.resto_outer_iters = self.resto_outer_iters.saturating_add(1);
3154 let resto_t0 = std::time::Instant::now();
3155 let outcome = resto.perform_restoration(&self.data, &self.cq, nlp, aug);
3156 drop(pd_guard);
3157 self.resto_wall_secs += resto_t0.elapsed().as_secs_f64();
3158 self.resto_inner_iters = self
3159 .resto_inner_iters
3160 .saturating_add(resto.last_inner_iter_count());
3161 // pounce#244: the restoration inner IPM shares the outer solve's
3162 // `Deadline` (both its convergence check and — post-#244 — its KKT
3163 // solves consult it), so a budget crossing inside restoration
3164 // terminates the inner solve with a time-limit status. Surface that
3165 // as the time limit directly instead of letting the `Failed` arm map
3166 // it onto `RestorationFailure` / `StopAtAcceptablePoint`. `data.curr`
3167 // is the last accepted outer iterate — restoration stages its
3168 // recovered point onto `trial`, not `curr`, and we return before
3169 // promoting it — so this hands back a valid iterate.
3170 if let Some(ret) = self.deadline_status() {
3171 return IterateOutcome::Terminate(ret);
3172 }
3173 match outcome {
3174 RestorationOutcome::Recovered => {
3175 // Mirror upstream `IpBacktrackingLineSearch.cpp:624-631`:
3176 // a successful restoration clears the line search's
3177 // cross-iteration globalization counters. Upstream runs
3178 // restoration inside `FindAcceptableTrialPoint` so those
3179 // assignments are inline; pounce runs it here, so the
3180 // reset has to be driven from here. Without it
3181 // `watchdog_shortened_iter` survives a restoration
3182 // episode and runs of shortened steps on either side of
3183 // one accumulate as if consecutive, arming the watchdog
3184 // where upstream would not. See
3185 // `BacktrackingLineSearch::reset_after_restoration`.
3186 self.bundle.line_search.reset_after_restoration();
3187 // The driver has staged the recovered point on
3188 // `data.trial`; apply the safe-slack bound adjustment
3189 // (as the main accept path does), then promote it and
3190 // continue iterating.
3191 self.adjust_variable_bounds_for_small_slacks();
3192 self.data.borrow_mut().accept_trial_point();
3193 // Snapshot the recovery iterate for the slow-cycle
3194 // detector at the top of the next `invoke_restoration`.
3195 // Compared against next-entry curr, dx_rel ≈ ‖α·d‖ —
3196 // measures purely the outer step. See header comment
3197 // on the cycle detector above.
3198 let recovered = self
3199 .data
3200 .borrow()
3201 .curr
3202 .as_ref()
3203 .expect("accept_trial_point sets curr")
3204 .clone();
3205 self.last_resto_recovery_x = Some(recovered.x.make_new_copy());
3206 self.last_resto_recovery_s = Some(recovered.s.make_new_copy());
3207 // Mirror upstream `IpoptAlgorithm::AcceptTrialPoint`
3208 // (`IpIpoptAlg.cpp:917-963`): kappa_sigma clamp on the
3209 // four bound-multiplier vectors. Upstream applies this
3210 // unconditionally inside AcceptTrialPoint, so the
3211 // post-restoration path inherits it; pounce factored
3212 // the clamp out of the data swap so we must call it
3213 // explicitly here. Without it the all-1 multiplier
3214 // reset (`bound_mult_reset_threshold`) leaves z*s far
3215 // from mu at the recovered iterate, blowing up the
3216 // next KKT solve's σ = z/s diagonal.
3217 self.correct_bound_multiplier();
3218 IterateOutcome::Continue
3219 }
3220 RestorationOutcome::Failed => {
3221 // Mirrors upstream `IpBacktrackingLineSearch.cpp:611-623`:
3222 // when `PerformRestoration` returns false, attempt to
3223 // roll back to the most recent acceptable iterate before
3224 // surfacing failure. If a snapshot is available we exit
3225 // cleanly with `StopAtAcceptablePoint` (mapped by the
3226 // application layer to `Solved_To_Acceptable_Level`),
3227 // matching the upstream `ACCEPTABLE_POINT_REACHED`
3228 // throw. Without a snapshot we surface
3229 // `RestorationFailure` — unless the restoration left the
3230 // iterate diverging (`|x|_∞ > diverging_iterates_tol`), in
3231 // which case we surface `DivergingIterates` to mirror the
3232 // outcome upstream produces on pathological problems like
3233 // MESH (where ipopt reports `Diverging_Iterates` and
3234 // pounce previously reported `Restoration_Failed` with an
3235 // obj of −3.6e+33). As in the running guard above, a large
3236 // `|x|` is only reported as unbounded when it is
3237 // structurally consistent with an unbounded feasible region
3238 // and the divergence is genuine — either it has persisted
3239 // (the running guard's growth-and-descent streak, which only
3240 // accumulates on a real recession ray; issues #248 / #252) or
3241 // blown past the absolute runaway backstop; otherwise the
3242 // failure is a plain `RestorationFailure`, never a spurious
3243 // `Unbounded`.
3244 if self.restore_acceptable_point() {
3245 IterateOutcome::Terminate(SolverReturn::StopAtAcceptablePoint)
3246 } else {
3247 let diverging = {
3248 let data = self.data.borrow();
3249 match data.curr.as_ref() {
3250 Some(curr) => {
3251 let amax = curr.x.amax();
3252 amax > self.diverging_iterates_tol
3253 && self.divergence_is_true_unboundedness(&*curr.x)
3254 && (amax >= Self::DIVERGENCE_ABS_RUNAWAY
3255 || self.divergence_streak >= Self::DIVERGENCE_PERSIST_ITERS)
3256 }
3257 None => false,
3258 }
3259 };
3260 if diverging {
3261 IterateOutcome::Terminate(SolverReturn::DivergingIterates)
3262 } else {
3263 IterateOutcome::Terminate(SolverReturn::RestorationFailure)
3264 }
3265 }
3266 }
3267 RestorationOutcome::LocallyInfeasible => {
3268 // Mirrors upstream's catch of `LOCALLY_INFEASIBLE` thrown
3269 // from `IpRestoConvCheck.cpp:240` — the resto sub-IPM
3270 // settled at a stationary point of `||c(x)||_1` whose
3271 // residual is still well above `tol`. Without this
3272 // detection the outer would re-enter restoration on the
3273 // unchanged iterate forever.
3274 //
3275 // gh #505: consult the acceptable-point stash, for the same
3276 // reason the conv-check arm above does and the cycle exits
3277 // already did. This is the *third* site that produced
3278 // `LocalInfeasibility`, and the only one not gated on
3279 // `infeas_max_streak` — which matters, because on the reported
3280 // instance raising that knob to 15 did not move the run by a
3281 // single iteration, so the verdict there is not the outer
3282 // detector's. Whichever route reaches it, a solve that passed
3283 // through an acceptable iterate must not discard it.
3284 //
3285 // Inert on genuinely infeasible models by the same argument as
3286 // the other two: nothing is stashed unless the whole acceptable
3287 // triplet passed, so `terminate_acceptable_or` falls through to
3288 // the verdict unchanged.
3289 self.terminate_local_infeasibility()
3290 }
3291 }
3292 }
3293
3294 /// Safe-slack bound adjustment, applied to the staged `trial`
3295 /// iterate before it is promoted to `curr`. When one or more trial
3296 /// slacks fell below `eps*min(1,mu)`, [`IpoptCalculatedQuantities::
3297 /// adjusted_trial_bounds`] returns the moved `x_L/x_U/d_L/d_U`; we
3298 /// install them on the NLP so the slack becomes representable. Port
3299 /// of the bound-adjustment block in `IpoptAlgorithm::AcceptTrialPoint`
3300 /// (`IpIpoptAlg.cpp:664-706`).
3301 fn adjust_variable_bounds_for_small_slacks(&mut self) {
3302 // Compute the moved bounds (releases the CQ/NLP borrows on return).
3303 let adjusted = {
3304 let trial_set = self.data.borrow().trial.is_some();
3305 if !trial_set {
3306 return;
3307 }
3308 self.cq.borrow().adjusted_trial_bounds()
3309 };
3310 let Some(bounds) = adjusted else {
3311 return;
3312 };
3313 tracing::debug!(
3314 target: "pounce::algorithm",
3315 "slack_move: {} slack(s) too small, adjusting variable bound(s) at iter {}",
3316 bounds.adjusted,
3317 self.data.borrow().iter_count,
3318 );
3319 let nlp = Rc::clone(self.cq.borrow().nlp());
3320 nlp.borrow_mut().adjust_variable_bounds(
3321 &*bounds.x_l,
3322 &*bounds.x_u,
3323 &*bounds.d_l,
3324 &*bounds.d_u,
3325 );
3326 }
3327
3328 /// Port of `IpIpoptAlg::correct_bound_multiplier`
3329 /// (`IpIpoptAlg.cpp:1055-1134`). Clamp each bound multiplier
3330 /// component into `[mu/(kappa_sigma * s_i), kappa_sigma * mu / s_i]`
3331 /// for all four bound-multiplier vectors.
3332 fn correct_bound_multiplier(&mut self) {
3333 if self.kappa_sigma < 1.0 {
3334 return;
3335 }
3336 let mu = self.data.borrow().curr_mu;
3337 let curr = match self.data.borrow().curr.clone() {
3338 Some(c) => c,
3339 None => return,
3340 };
3341
3342 let cq = self.cq.borrow();
3343
3344 let z_l_new = clamp_against_slack(&*curr.z_l, &*cq.curr_slack_x_l(), mu, self.kappa_sigma);
3345 let z_u_new = clamp_against_slack(&*curr.z_u, &*cq.curr_slack_x_u(), mu, self.kappa_sigma);
3346 let v_l_new = clamp_against_slack(&*curr.v_l, &*cq.curr_slack_s_l(), mu, self.kappa_sigma);
3347 let v_u_new = clamp_against_slack(&*curr.v_u, &*cq.curr_slack_s_u(), mu, self.kappa_sigma);
3348 drop(cq);
3349
3350 let new_iv = crate::iterates_vector::IteratesVector::new(
3351 curr.x.clone(),
3352 curr.s.clone(),
3353 curr.y_c.clone(),
3354 curr.y_d.clone(),
3355 z_l_new,
3356 z_u_new,
3357 v_l_new,
3358 v_u_new,
3359 );
3360 self.data.borrow_mut().set_curr(new_iv);
3361 }
3362
3363 /// Outer entry point — port of `IpoptAlgorithm::Optimize()`. Calls
3364 /// the iterate-initializer once, then loops `iterate()` until a
3365 /// terminal status. The exception → SolverReturn mapping
3366 /// (TINY_STEP_DETECTED → STEP_BECOMES_TINY,
3367 /// RESTORATION_FAILED → RESTORATION_FAILURE, etc.) lands in
3368 /// Phase 9 alongside the restoration phase.
3369 /// Run the solve and finalize its result.
3370 ///
3371 /// A thin wrapper on purpose. The gh #200 fallback must see **every** exit
3372 /// of the driver loop, and wiring it into individual termination sites was
3373 /// tried and failed — there are sixteen, and the ones easiest to overlook
3374 /// are the ones most likely to matter. Keeping the loop in a separate
3375 /// function means every `return` inside it, present or future, flows through
3376 /// [`Self::honour_refused_certificate`] by construction rather than by the
3377 /// author remembering to.
3378 ///
3379 /// This got more important once the fallback started changing the status in
3380 /// *both* directions: it can now hand back `StopAtAcceptablePoint` for a
3381 /// `Success` it was given. Anything reading `result` before the hook is
3382 /// reading a status that is not the one reported.
3383 pub fn optimize(&mut self) -> SolverReturn {
3384 let result = self.optimize_inner();
3385
3386 // gh #200: a refused certificate outranks any non-success verdict the
3387 // continued run reached, and an earlier refusal can outrank the
3388 // continued run's own certificate. Applied here, once.
3389 let result = self.honour_refused_certificate(result);
3390
3391 // pounce#250 follow-up: the dual-divergence guard's diversion to
3392 // restoration is a bet, and a lost bet must not return a worse point
3393 // than the solve already had in hand. Applied here, once, for the same
3394 // reason the #200 hook is — every `return` in the loop flows through
3395 // this point by construction.
3396 let result = self.honour_best_acceptable_after_dual_guard(result);
3397
3398 // gh #534: deferring the acceptable-point restoration decline is also a
3399 // bet, and this is the net under it — a continuation that did not beat
3400 // the point the guard would have returned hands that point back. Last of
3401 // the three, so it compares against whatever the hooks above settled on.
3402 let result = self.honour_decline_floor(result);
3403
3404 // Terminal post-mortem checkpoint. Skipped when the user already
3405 // asked to stop (they were just at a prompt); otherwise the
3406 // debugger gets a last look at the final/failing iterate.
3407 if !matches!(result, SolverReturn::UserRequestedStop) {
3408 self.fire_debug_terminal(result);
3409 }
3410 result
3411 }
3412
3413 fn optimize_inner(&mut self) -> SolverReturn {
3414 // Top-level span for the whole solve; every iteration / linear
3415 // solve / restoration event nests under it (pounce#71).
3416 let _solve_span = tracing::info_span!("solve").entered();
3417
3418 // Shared timing accumulator — every phase below records into it.
3419 let timing = self.data.borrow().timing.clone();
3420
3421 // Install the shared accumulator on the augmented-system solver
3422 // so its factor / back-solve calls are attributed to
3423 // `linear_system_factorization` / `linear_system_back_solve`.
3424 // Same pattern for the diagnostics state when present, so KKT
3425 // dump sites can consult per-iter gating.
3426 if let Some(sd) = self.search_dir.as_mut() {
3427 sd.pd_solver_mut()
3428 .aug_solver_mut()
3429 .set_timing_stats(std::rc::Rc::clone(&timing));
3430 if let Some(diag) = self.diagnostics.as_ref() {
3431 sd.pd_solver_mut()
3432 .aug_solver_mut()
3433 .set_diagnostics(Rc::clone(diag));
3434 }
3435 }
3436
3437 // 0a. Strategy initialization — port of upstream's
3438 // `IpoptAlgorithm::InitializeImpl` calls. The mu update needs
3439 // `data.curr_mu`/`curr_tau` seeded before the iterate
3440 // initializer runs (`CalculateSafeSlack` reads them).
3441 self.bundle.mu_update.initialize(&self.data);
3442
3443 // 0b. Iterate initializer. Requires NLP; without one the caller
3444 // must have populated `data.curr` themselves.
3445 if let Some(nlp) = self.nlp.as_ref() {
3446 // The initializer needs an aug-system solver for the
3447 // least-square multiplier branch; until that's wired we
3448 // route through whatever the search-direction calculator
3449 // owns when present. For the stub flow we skip the LSM
3450 // path by giving the initializer a dummy solver only if
3451 // the search_dir is present (otherwise the init function
3452 // is responsible for not consulting it).
3453 if let Some(sd) = self.search_dir.as_mut() {
3454 timing.initialize_iterates.start();
3455 let mut pd_guard = sd.pd_solver_mut();
3456 let aug_solver = pd_guard.aug_solver_mut();
3457 let ok = self
3458 .bundle
3459 .init
3460 .set_initial_iterates(&self.data, &self.cq, nlp, aug_solver);
3461 drop(pd_guard);
3462 timing.initialize_iterates.end();
3463 if !ok {
3464 return SolverReturn::InvalidProblemDefinition;
3465 }
3466 }
3467 }
3468
3469 // 0c. Seed `IpoptData::w` with the initial-iterate Hessian.
3470 // Redundant with the iter-body `update_hessian` call (which
3471 // now runs BEFORE `update_barrier_parameter`) but kept to
3472 // cover any code path that consults `data.w` between
3473 // `set_initial_iterates` and the first `iterate()` call
3474 // (e.g. the iter-0 trace dump below).
3475 if self.data.borrow().curr.is_some() {
3476 timing.update_hessian.start();
3477 let _ = self.bundle.hess.update_hessian(&self.data, &self.cq);
3478 timing.update_hessian.end();
3479 }
3480
3481 // Track-A iterate-trace dumper. Activated by
3482 // `IPOPT_ITER_DUMP_PATH`; otherwise no-op. See `iter_dump.rs`.
3483 let mut dumper = IterDumper::from_env();
3484 // Iter 0 record — captures the initialised iterate before any
3485 // step. Mirrors upstream's "after InitializeIterates(), before
3486 // the loop" emission point.
3487 if let Some(d) = dumper.as_mut() {
3488 d.write_record(&self.data, &self.cq);
3489 }
3490
3491 // Advance the diagnostics iter counter so the first `iterate()`
3492 // body reports as iter 0 (matches `data.iter_count`). Subsequent
3493 // bumps live at the bottom of the loop alongside the iter_count
3494 // bookkeeping.
3495 if let Some(diag) = self.diagnostics.as_ref() {
3496 diag.bump_iter();
3497 // Iter-0 iterate row (issue #68). Same hook point as
3498 // the binary IterDumper above; emits only when
3499 // `--dump iterates:*` is configured.
3500 emit_iterate_record(diag.as_ref(), &self.data, &self.cq);
3501 }
3502
3503 // Iter 0 intermediate callback — upstream fires once after
3504 // `InitializeIterates` before the loop body starts so users
3505 // observe the initial point.
3506 if !self.fire_intermediate() {
3507 return SolverReturn::UserRequestedStop;
3508 }
3509 if self.fire_debug(crate::debug::Checkpoint::IterStart) == crate::debug::DebugAction::Stop {
3510 return SolverReturn::UserRequestedStop;
3511 }
3512
3513 // pounce#246: bound the initialization / restoration-entry window.
3514 // Everything above — `mu_update.initialize`, `set_initial_iterates`
3515 // (which for a bad warm start can grind in its least-square /
3516 // feasibility setup), and the initial `update_hessian` — runs
3517 // *before* the first `iterate()`, whose convergence check and
3518 // post-`compute_search_direction` gate are the earliest deadline
3519 // checks (#242/#244/#245). A solve handed a poor warm start could
3520 // therefore spend the whole budget here and only consult the
3521 // deadline once it reached the first outer-iteration /
3522 // KKT-factorization boundary. Consult it now, before the loop, so a
3523 // bad-start init stall returns promptly with the time-limit status
3524 // (best-so-far being the initialised iterate) instead of running to
3525 // a multiple of the budget. This also bounds the *restoration
3526 // entry*: the nested restoration IPM shares this `Deadline` and runs
3527 // the same `optimize_inner`, so a budget already crossed by the time
3528 // the inner solve starts up terminates it here rather than after its
3529 // own first iterate.
3530 if let Some(ret) = self.deadline_status() {
3531 return ret;
3532 }
3533
3534 let result = loop {
3535 match self.iterate() {
3536 IterateOutcome::Terminate(ret) => break ret,
3537 IterateOutcome::Continue => {
3538 // Source the local counter from `data.iter_count`
3539 // each pass so a pre-seeded counter (e.g. the inner
3540 // restoration IPM at `outer.iter + 1`, matching
3541 // upstream `IpRestoMinC_1Nrm.cpp:181`) and any
3542 // restoration step that set
3543 // `data.iter_count = inner.iter_count - 1`
3544 // (mirroring `IpRestoMinC_1Nrm.cpp:Set_iter_count`)
3545 // are honored — without this the local counter
3546 // would advance from its pre-restoration value,
3547 // ignoring the inner-IPM iterations.
3548 let mut iter_count: Index = self.data.borrow().iter_count;
3549 iter_count += 1;
3550 // Do NOT short-circuit to `MaxiterExceeded` here: bump the
3551 // counter and loop, letting the next `iterate()` run its
3552 // convergence check (`OptimalityErrorConvergenceCheck`,
3553 // which tests the component tolerances *before* its own
3554 // `iter_count >= max_iter` gate at
3555 // `conv_check/opt_error.rs`). Breaking before that call
3556 // skipped the convergence test on the iterate produced by
3557 // the final permitted step, so a solve converging on
3558 // exactly the `max_iter`-th iterate reported
3559 // `Maximum_Iterations_Exceeded` where upstream Ipopt —
3560 // which runs `CheckConvergence` at the top of its loop,
3561 // convergence-first — reports success. The check is
3562 // guaranteed to terminate the loop: once `iter_count`
3563 // reaches `max_iter`, `check_convergence_with_state`
3564 // returns either `Converged`/`ConvergedToAcceptable` or
3565 // `MaxIterExceeded`, never `Continue` (L1).
3566 self.data.borrow_mut().iter_count = iter_count;
3567 // Keep the diagnostics counter in lock-step with
3568 // `data.iter_count` so KKT-dump gating reflects the
3569 // about-to-execute iteration.
3570 if let Some(diag) = self.diagnostics.as_ref() {
3571 diag.bump_iter();
3572 // Per-iter iterate row (issue #68). Mirrors
3573 // the binary IterDumper hook below.
3574 emit_iterate_record(diag.as_ref(), &self.data, &self.cq);
3575 }
3576 // Per-iteration record — emitted after the
3577 // iter_count bump so the recorded `iter` field
3578 // matches `IpData().iter_count()` at the moment of
3579 // emission, identical to upstream's writer.
3580 if let Some(d) = dumper.as_mut() {
3581 d.write_record(&self.data, &self.cq);
3582 }
3583 // Per-iteration intermediate callback — fired with
3584 // an `IntermediateContext` guard so downstream
3585 // inspector entry points (the C API
3586 // `GetIpoptCurrent*` family) see live state for the
3587 // duration of the user callback.
3588 if !self.fire_intermediate() {
3589 break SolverReturn::UserRequestedStop;
3590 }
3591 if self.fire_debug(crate::debug::Checkpoint::IterStart)
3592 == crate::debug::DebugAction::Stop
3593 {
3594 break SolverReturn::UserRequestedStop;
3595 }
3596 }
3597 }
3598 };
3599
3600 result
3601 }
3602}
3603
3604/// A termination certificate the masked-scale veto refused (gh #200), with
3605/// everything the fallback needs to undo the refusal verbatim.
3606///
3607/// One struct rather than a field per component. The fallback is only correct if
3608/// these all describe *the same iterate*, and parallel `Option`s make
3609/// "objective recorded, iterate missing" representable — which was reachable:
3610/// the iterate is cloned out of `data.curr` and can come back `None`, while the
3611/// objective and barrier parameter were written unconditionally. Capture is now
3612/// all-or-nothing, so the disagreement cannot be constructed.
3613#[derive(Clone)]
3614struct VetoSnapshot {
3615 /// The refused iterate itself.
3616 iterate: crate::iterates_vector::IteratesVector,
3617 /// Iteration at which the refusal happened.
3618 ///
3619 /// Needed to identify which refusal is the *baseline-equivalent* one. The
3620 /// baseline stops at the first iterate where it would terminate, so when
3621 /// both a strict and an acceptable-level refusal are on record it is the
3622 /// chronologically earlier one that says what the baseline returned — not
3623 /// the stricter one. The later refusal sits on the continued trajectory,
3624 /// which the baseline never walked, so comparing against it compares
3625 /// against a point that was never on offer.
3626 iter: Index,
3627 /// Scaled objective there, so the refused point can be compared against
3628 /// whatever the continued run reached without re-evaluating it.
3629 obj: Number,
3630 /// Barrier parameter there.
3631 ///
3632 /// `curr_mu` lives on `IpoptData` rather than in the `IteratesVector`, so
3633 /// restoring the iterate does not rewind it, and `stats.final_mu` is read
3634 /// after the restore — leaving the continued run's barrier parameter
3635 /// reported next to the refused run's `x`. That pair feeds a warm-started
3636 /// corrector's `mu_init` and reaches callers as `info["mu"]`, so it must
3637 /// describe the point actually returned. (Not currently observable: `mu` has
3638 /// bottomed out at its floor in every fallback case reachable so far, making
3639 /// the two values coincide. Kept correct rather than left to depend on that.)
3640 mu: Number,
3641 /// Max-norm unscaled KKT error there, so the tiebreak can see what
3642 /// `apply_kkt_fidelity_gate` will see. Recorded at refusal time because the
3643 /// gate runs post-solve, long after this iterate is gone.
3644 unscaled_kkt: Number,
3645 /// Unscaled max-norm constraint violation there — the same quantity the
3646 /// `acceptable_constr_viol_tol` gate is defined against
3647 /// (`curr_unscaled_primal_infeasibility_max`, cf. gh #261). The
3648 /// best-acceptable fallback ranks candidates by `(feasible_enough,
3649 /// objective)` rather than objective alone, so a point outside a capped
3650 /// feasibility band can never displace one inside it on objective grounds;
3651 /// without this field the ranking has no feasibility term and, under a
3652 /// user-widened `acceptable_constr_viol_tol`, will trade feasibility for
3653 /// objective and hand back a verifiably infeasible point under a success
3654 /// status (gh #267). Unused by the gh #200 masked-scale paths, which key on
3655 /// objective only.
3656 constr_viol: Number,
3657 /// Objective scaling factor in force when `obj` was recorded.
3658 ///
3659 /// `obj` is a *scaled* objective, so comparing it against the continued
3660 /// run's is only meaningful under the same factor, sign included. Held so
3661 /// that assumption is asserted rather than trusted: periodic or adaptive
3662 /// rescaling is a natural thing to add for exactly the ill-scaled problems
3663 /// this mechanism targets, and it would silently turn the comparison into
3664 /// noise.
3665 obj_scale: Number,
3666}
3667
3668/// Internal result of one [`IpoptAlgorithm::iterate`] call. Mirrors the
3669/// upstream try/catch around `IpoptAlg::Optimize` — anything that's not
3670/// `Continue` carries the [`SolverReturn`] that the outer loop will
3671/// surface to `IpoptApplication`.
3672enum IterateOutcome {
3673 Continue,
3674 Terminate(SolverReturn),
3675}
3676
3677/// Feasibility-aware ranking core for the best-acceptable fallback (gh #267).
3678///
3679/// `true` iff candidate `(a_obj, a_viol)` ranks **strictly** better than
3680/// `(b_obj, b_viol)` under the `(band_clamped_viol, objective)` key, where each
3681/// violation is clamped up to `band` before it is compared. Lexicographic: the
3682/// point with the smaller clamped violation wins outright; only when the clamped
3683/// violations tie does the lower objective decide. A non-finite objective ranks
3684/// worst (never wins, always loses to a finite one); a non-finite violation is
3685/// treated as infinitely infeasible.
3686///
3687/// The clamp is what makes this a **total order** rather than a two-class
3688/// partition, and it is the whole gh #280 fix. Clamping the violation *up* to
3689/// `band` collapses every point inside the feasibility band to the single value
3690/// `band`, so within the band those points tie on feasibility and objective
3691/// decides — the intended, gh #267-preserving behaviour. Outside the band the
3692/// clamp is the identity, so the *actual* violation decides and the
3693/// less-infeasible point always wins. The earlier `(feasible_enough, objective)`
3694/// key was a two-class partition: once **both** points sat outside the band it
3695/// fell through to a bare `a_obj < b_obj`, reading neither violation — the exact
3696/// pre-#267 objective-only rule, which lets a strictly-more-infeasible point win
3697/// on objective (gh #280). Under the clamped key a strictly-more-infeasible
3698/// point can never rank better, at any band.
3699///
3700/// Pure and total, so the fallback's "never worse off" guarantee is a theorem
3701/// this function's unit tests prove by cases — host-independent by construction,
3702/// unlike an end-to-end objective comparison across two live nonconvex solves
3703/// (the trap gh #267 caught). [`IpoptAlgorithm::ranks_better`] supplies `band`
3704/// as `min(acceptable_constr_viol_tol, FEASIBLE_ENOUGH_CAP)`; both the record
3705/// and the read side route through here, so they cannot disagree.
3706/// Append `value` to a fixed-capacity oldest-first window, dropping the oldest
3707/// sample once the window is full (gh #534).
3708fn push_sample(buf: &mut [Number; DECLINE_PROGRESS_SAMPLES], len: &mut usize, value: Number) {
3709 if *len < DECLINE_PROGRESS_SAMPLES {
3710 buf[*len] = value;
3711 *len += 1;
3712 } else {
3713 buf.rotate_left(1);
3714 buf[DECLINE_PROGRESS_SAMPLES - 1] = value;
3715 }
3716}
3717
3718/// Whether every consecutive pair in `samples` (oldest first) contracted by at
3719/// least `ratio` — the gh #534 progress test, as a pure function.
3720///
3721/// A sample that is not finite, or a predecessor that is not strictly positive,
3722/// fails the window: neither is evidence of progress, and a zero predecessor
3723/// makes the ratio meaningless. A `ratio` of `1` admits any non-increasing
3724/// window, and a large one admits every finite window — which is how
3725/// `resto_decline_progress_ratio` doubles as the "drop the progress
3726/// requirement" switch.
3727///
3728/// Pure and total for the same reason [`ranks_better_within_band`] is: the two
3729/// traces the issue records — `eigena2` quartering and `eigenb2` rising — decide
3730/// what this must do, and a unit test can hold it to them exactly.
3731fn window_is_contracting(samples: &[Number], ratio: Number) -> bool {
3732 samples.windows(2).all(|w| {
3733 let (prev, next) = (w[0], w[1]);
3734 prev.is_finite() && prev > 0.0 && next.is_finite() && next <= ratio * prev
3735 })
3736}
3737
3738fn ranks_better_within_band(
3739 a_obj: Number,
3740 a_viol: Number,
3741 b_obj: Number,
3742 b_viol: Number,
3743 band: Number,
3744) -> bool {
3745 if !a_obj.is_finite() {
3746 return false;
3747 }
3748 if !b_obj.is_finite() {
3749 return true;
3750 }
3751 // Clamp each violation up to `band`: everything inside the feasibility band
3752 // maps to the single value `band` (so objective decides there), while outside
3753 // it the actual violation is kept (so the less-infeasible point strictly
3754 // wins). A non-finite violation is infinitely infeasible. This is a total
3755 // order — a strictly-more-infeasible point can never rank better (gh #280).
3756 let clamped = |v: Number| {
3757 if v.is_finite() {
3758 v.max(band)
3759 } else {
3760 Number::INFINITY
3761 }
3762 };
3763 let (a_key, b_key) = (clamped(a_viol), clamped(b_viol));
3764 if a_key != b_key {
3765 // The less-infeasible point wins outright, whatever the objectives.
3766 return a_key < b_key;
3767 }
3768 // Same clamped feasibility (both inside the band, or an exact tie outside it):
3769 // lower objective wins (the original within-band behaviour).
3770 a_obj < b_obj
3771}
3772
3773/// `||a - b||_2 / (1 + ||b||_2)`. Used by the restoration cycle
3774/// detector in [`IpoptAlgorithm::invoke_restoration`] to test whether
3775/// the outer iterate has moved between two consecutive restoration
3776/// entries.
3777fn relative_distance(a: &dyn Vector, b: &dyn Vector) -> Number {
3778 if a.dim() == 0 {
3779 return 0.0;
3780 }
3781 let mut diff = a.make_new_copy();
3782 diff.axpy(-1.0, b);
3783 diff.nrm2() / (1.0 + b.nrm2())
3784}
3785
3786/// `out = curr + α_p · δ` for the primal/equality blocks and
3787/// `out = curr + α_d · δ` for the bound multipliers, returned as a
3788/// fresh frozen `IteratesVector`. Mirrors `scaled_step` in the line
3789/// search; duplicated here for the tiny-step branch which bypasses
3790/// the line-search driver.
3791fn scaled_step_unchecked(
3792 curr: &crate::iterates_vector::IteratesVector,
3793 delta: &crate::iterates_vector::IteratesVector,
3794 alpha_primal: Number,
3795 alpha_dual: Number,
3796) -> crate::iterates_vector::IteratesVector {
3797 let mut out = curr.make_new_zeroed();
3798 out.add_one_vector(1.0, curr, 0.0);
3799 out.x.axpy(alpha_primal, &*delta.x);
3800 out.s.axpy(alpha_primal, &*delta.s);
3801 out.y_c.axpy(alpha_primal, &*delta.y_c);
3802 out.y_d.axpy(alpha_primal, &*delta.y_d);
3803 out.z_l.axpy(alpha_dual, &*delta.z_l);
3804 out.z_u.axpy(alpha_dual, &*delta.z_u);
3805 out.v_l.axpy(alpha_dual, &*delta.v_l);
3806 out.v_u.axpy(alpha_dual, &*delta.v_u);
3807 out.freeze()
3808}
3809
3810/// Allocate a fresh `Rc<dyn Vector>` with `kappa_sigma_clamp`
3811/// applied component-wise against the supplied `slack`. Inputs are
3812/// borrowed; the original `z` is never mutated. Ports the per-vector
3813/// piece of `IpIpoptAlg.cpp:1080-1133`.
3814fn clamp_against_slack(
3815 z: &dyn Vector,
3816 slack: &dyn Vector,
3817 mu: Number,
3818 kappa_sigma: Number,
3819) -> Rc<dyn Vector> {
3820 debug_assert_eq!(z.dim(), slack.dim());
3821 let n = z.dim() as usize;
3822 // Flatten both z and slack into contiguous slices so the
3823 // elementwise clamp doesn't care whether the inputs are
3824 // [`DenseVector`] (regular IPM path) or [`CompoundVector`]
3825 // (resto IPM path). The result is reconstructed into a
3826 // same-shape Vector via `Vector::make_new` + a flat-write
3827 // helper so the caller sees a vector with the same blocking as
3828 // its input.
3829 let mut buf = vec![0.0_f64; n];
3830 flat_read_into(z, &mut buf);
3831 let s_vals = flat_read_owned(slack);
3832 let _ = kappa_sigma_clamp(&mut buf, &s_vals, mu, kappa_sigma);
3833 let mut out: Box<dyn Vector> = z.make_new();
3834 flat_write_into(&mut *out, &buf);
3835 Rc::from(out)
3836}
3837
3838pub(crate) fn flat_read_into(v: &dyn Vector, dst: &mut [Number]) {
3839 if let Some(dv) = v
3840 .as_any()
3841 .downcast_ref::<pounce_linalg::dense_vector::DenseVector>()
3842 {
3843 let vs = dv.expanded_values();
3844 dst.copy_from_slice(&vs);
3845 return;
3846 }
3847 if let Some(cv) = v.as_any().downcast_ref::<pounce_linalg::CompoundVector>() {
3848 let mut off = 0usize;
3849 for k in 0..cv.n_comps() {
3850 let blk = cv.comp(k);
3851 let dim = blk.dim() as usize;
3852 let dblk = blk
3853 .as_any()
3854 .downcast_ref::<pounce_linalg::dense_vector::DenseVector>()
3855 .expect("clamp_against_slack: CompoundVector blocks must be DenseVectors");
3856 let vs = dblk.expanded_values();
3857 dst[off..off + dim].copy_from_slice(&vs);
3858 off += dim;
3859 }
3860 return;
3861 }
3862 panic!("clamp_against_slack: unsupported Vector kind");
3863}
3864
3865pub(crate) fn flat_read_owned(v: &dyn Vector) -> Vec<Number> {
3866 let mut out = vec![0.0; v.dim() as usize];
3867 flat_read_into(v, &mut out);
3868 out
3869}
3870
3871pub(crate) fn flat_write_into(v: &mut dyn Vector, src: &[Number]) {
3872 if let Some(dv) = v
3873 .as_any_mut()
3874 .downcast_mut::<pounce_linalg::dense_vector::DenseVector>()
3875 {
3876 dv.set_values(src);
3877 return;
3878 }
3879 if let Some(cv) = v
3880 .as_any_mut()
3881 .downcast_mut::<pounce_linalg::CompoundVector>()
3882 {
3883 let mut off = 0usize;
3884 for k in 0..cv.n_comps() {
3885 let blk = cv.comp_mut(k);
3886 let dim = blk.dim() as usize;
3887 let dblk = blk
3888 .as_any_mut()
3889 .downcast_mut::<pounce_linalg::dense_vector::DenseVector>()
3890 .expect("clamp_against_slack: CompoundVector blocks must be DenseVectors");
3891 dblk.set_values(&src[off..off + dim]);
3892 off += dim;
3893 }
3894 return;
3895 }
3896 panic!("clamp_against_slack: unsupported Vector kind");
3897}
3898
3899/// Per-element kappa-sigma clamp — the elementwise arithmetic at the
3900/// heart of `IpIpoptAlg.cpp:correct_bound_multiplier` (lines
3901/// 1090-1133). For each index `i`:
3902///
3903/// ```text
3904/// slack_i = max(slack_i, tiny_double) // avoid /0
3905/// z_lo_i = mu / (kappa_sigma * slack_i)
3906/// z_hi_i = kappa_sigma * mu / slack_i
3907/// z_i ← clamp(z_i, z_lo_i, z_hi_i)
3908/// ```
3909///
3910/// Returns the maximum elementwise correction magnitude (matching
3911/// upstream's `Max(max_correction_up, max_correction_low)`).
3912///
3913/// `kappa_sigma < 1` short-circuits to the identity per upstream's
3914/// guard at line 1065.
3915pub fn kappa_sigma_clamp(
3916 z: &mut [Number],
3917 slack: &[Number],
3918 mu: Number,
3919 kappa_sigma: Number,
3920) -> Number {
3921 debug_assert_eq!(z.len(), slack.len());
3922 if kappa_sigma < 1.0 {
3923 return 0.0;
3924 }
3925 let mut max_correction = 0.0_f64;
3926 for (zi, &si) in z.iter_mut().zip(slack.iter()) {
3927 let s_safe = si.max(Number::MIN_POSITIVE);
3928 let lo = mu / (kappa_sigma * s_safe);
3929 let hi = kappa_sigma * mu / s_safe;
3930 let clamped = zi.clamp(lo, hi);
3931 let delta = (clamped - *zi).abs();
3932 if delta > max_correction {
3933 max_correction = delta;
3934 }
3935 *zi = clamped;
3936 }
3937 max_correction
3938}
3939
3940#[cfg(test)]
3941mod tests {
3942 use super::*;
3943
3944 #[test]
3945 fn kappa_sigma_below_one_is_identity() {
3946 let mut z = vec![1.0, 2.0, 3.0];
3947 let slack = [1.0, 1.0, 1.0];
3948 let m = kappa_sigma_clamp(&mut z, &slack, 1.0, 0.5);
3949 assert_eq!(m, 0.0);
3950 assert_eq!(z, [1.0, 2.0, 3.0]);
3951 }
3952
3953 #[test]
3954 fn within_band_is_unchanged() {
3955 // mu=1, kappa=10, slack=1 → band [0.1, 10]. z=1 → unchanged.
3956 let mut z = vec![1.0];
3957 let slack = [1.0];
3958 let m = kappa_sigma_clamp(&mut z, &slack, 1.0, 10.0);
3959 assert_eq!(m, 0.0);
3960 assert_eq!(z, [1.0]);
3961 }
3962
3963 #[test]
3964 fn above_upper_clamped_down() {
3965 // mu=1, kappa=10, slack=1 → upper = 10. z=100 → 10.
3966 let mut z = vec![100.0];
3967 let slack = [1.0];
3968 let m = kappa_sigma_clamp(&mut z, &slack, 1.0, 10.0);
3969 assert!((m - 90.0).abs() < 1e-13);
3970 assert_eq!(z, [10.0]);
3971 }
3972
3973 #[test]
3974 fn below_lower_clamped_up() {
3975 // mu=1, kappa=10, slack=1 → lower = 0.1. z=0.001 → 0.1.
3976 let mut z = vec![0.001];
3977 let slack = [1.0];
3978 let m = kappa_sigma_clamp(&mut z, &slack, 1.0, 10.0);
3979 assert!((m - 0.099).abs() < 1e-13);
3980 assert!((z[0] - 0.1).abs() < 1e-15);
3981 }
3982
3983 #[test]
3984 fn returns_max_over_components() {
3985 let mut z = vec![100.0, 0.001];
3986 let slack = [1.0, 1.0];
3987 let m = kappa_sigma_clamp(&mut z, &slack, 1.0, 10.0);
3988 assert!((m - 90.0).abs() < 1e-13);
3989 assert_eq!(z[0], 10.0);
3990 assert!((z[1] - 0.1).abs() < 1e-15);
3991 }
3992
3993 #[test]
3994 fn slack_clamped_to_min_positive_avoids_division_by_zero() {
3995 let mut z = vec![1e100];
3996 let slack = [0.0];
3997 let _ = kappa_sigma_clamp(&mut z, &slack, 1.0, 10.0);
3998 assert!(z[0].is_finite() || z[0] == 1e100);
3999 }
4000
4001 /// The restoration slot is exercised structurally:
4002 /// `IpoptAlgorithm::with_restoration` accepts a
4003 /// `Box<dyn RestorationPhase>` and the trait's default
4004 /// `perform_restoration` returns `Failed`. End-to-end coverage
4005 /// (iterate() → line-search-Failed → restoration → recovered)
4006 /// lands in the Phase 9 integration suite alongside the nested
4007 /// IPM driver.
4008 struct _DummyResto;
4009 impl RestorationPhase for _DummyResto {}
4010
4011 // --------------------------------------------------------------
4012 // Best-acceptable fallback ranking (gh #267).
4013 //
4014 // These prove the "never worse off" guarantee host-independently, by
4015 // cases, on the pure ranking core — the property the earlier end-to-end
4016 // `hair_trigger_*` objective comparison could only *approximate* on one
4017 // host's basin luck (gh #267's secondary finding). `band` here stands in
4018 // for the resolved `min(acceptable_constr_viol_tol, FEASIBLE_ENOUGH_CAP)`.
4019 // --------------------------------------------------------------
4020
4021 #[test]
4022 fn ranks_better_is_a_strict_order_within_a_feasibility_class() {
4023 let band = 1e-2;
4024 // Both feasible: lower objective wins, strictly.
4025 assert!(ranks_better_within_band(-2.0, 1e-4, -1.0, 1e-4, band));
4026 assert!(!ranks_better_within_band(-1.0, 1e-4, -2.0, 1e-4, band));
4027 // Ties are not "strictly better" in either direction — so an equal
4028 // returned point is never displaced, matching the read side's
4029 // keep-current-on-tie contract.
4030 assert!(!ranks_better_within_band(-1.0, 1e-4, -1.0, 5e-3, band));
4031 assert!(!ranks_better_within_band(-1.0, 5e-3, -1.0, 1e-4, band));
4032 // Both infeasible: the less-infeasible point wins. Here it also has the
4033 // lower objective, so this held under the old objective-only fall-through
4034 // too — `ranks_better_puts_feasibility_first_among_two_infeasibles` is the
4035 // case that separates the two rules (gh #280).
4036 assert!(ranks_better_within_band(-2.0, 5.0, -1.0, 9.0, band));
4037 }
4038
4039 #[test]
4040 fn ranks_better_puts_feasibility_first_among_two_infeasibles() {
4041 // The gh #280 hole: once BOTH points sit outside the (capped) band the
4042 // old `(feasible_enough, objective)` partition read a_ok == b_ok == false
4043 // and fell through to `a_obj < b_obj` — objective alone, the exact
4044 // pre-#267 rule — so a strictly-MORE-infeasible point could win by having
4045 // a better objective. This is the deb7 swap the fallback made: incumbent
4046 // at viol 5.292e-1, recorded point at viol 9.951e-1 with a 36%-better
4047 // objective. The less-infeasible point must win regardless of objective.
4048 let band = 1e-2;
4049 // Less-infeasible incumbent, WORSE objective — must still win.
4050 assert!(ranks_better_within_band(
4051 89.0, 5.292e-1, 56.9, 9.951e-1, band
4052 ));
4053 // The more-infeasible, better-objective point must NOT win — the swap
4054 // gh #280 forbids.
4055 assert!(!ranks_better_within_band(
4056 56.9, 9.951e-1, 89.0, 5.292e-1, band
4057 ));
4058 // A strictly-more-infeasible point never replaces the incumbent even with
4059 // an arbitrarily better objective.
4060 assert!(!ranks_better_within_band(-1e12, 9.0, 0.0, 5.0, band));
4061 assert!(ranks_better_within_band(0.0, 5.0, -1e12, 9.0, band));
4062 }
4063
4064 #[test]
4065 fn ranks_better_puts_feasibility_first_regardless_of_objective() {
4066 let band = 1e-2;
4067 // The gh #267 case in miniature: a feasible point outranks an
4068 // arbitrarily-lower-objective infeasible one, and vice-versa.
4069 assert!(ranks_better_within_band(0.0, 1e-4, -1e9, 9.94, band));
4070 assert!(!ranks_better_within_band(-1e9, 9.94, 0.0, 1e-4, band));
4071 // Exactly at the band is still feasible_enough; just past it is not.
4072 assert!(ranks_better_within_band(
4073 1.0,
4074 band,
4075 -1.0,
4076 band * 1.000_001,
4077 band
4078 ));
4079 }
4080
4081 #[test]
4082 fn ranks_better_never_lets_a_widened_band_matter_past_the_cap() {
4083 // The band the method feeds is capped at FEASIBLE_ENOUGH_CAP, so a
4084 // point beyond the cap is never feasible_enough however loose the
4085 // user's `acceptable_constr_viol_tol`. Model that by passing the capped
4086 // band: the near-optimal-but-mildly-infeasible endpoint (viol 1.13e-4,
4087 // within the cap) must beat the grossly-infeasible lower-objective
4088 // point (viol 9.94, past it) — the exact swap the fix forbids.
4089 let band = IpoptAlgorithm::FEASIBLE_ENOUGH_CAP;
4090 assert!(ranks_better_within_band(
4091 -2303.99, 1.13e-4, -2307.32, 9.94, band
4092 ));
4093 assert!(!ranks_better_within_band(
4094 -2307.32, 9.94, -2303.99, 1.13e-4, band
4095 ));
4096 }
4097
4098 #[test]
4099 fn ranks_better_ranks_a_nonfinite_objective_worst() {
4100 let band = 1e-2;
4101 for bad in [f64::NAN, f64::INFINITY, f64::NEG_INFINITY] {
4102 // A non-finite objective never ranks better than a finite one...
4103 assert!(!ranks_better_within_band(bad, 0.0, 0.0, 9.9, band));
4104 // ...and always loses to one, even on worse feasibility — so a
4105 // finite recorded point is restored over a NaN-objective return,
4106 // preserving the pre-fix `!(curr <= best)` NaN behaviour.
4107 assert!(ranks_better_within_band(0.0, 9.9, bad, 0.0, band));
4108 }
4109 }
4110
4111 #[test]
4112 fn ranks_better_treats_a_nonfinite_violation_as_infeasible() {
4113 let band = 1e-2;
4114 // A NaN violation is never feasible_enough, so a genuinely feasible
4115 // point outranks it regardless of objective.
4116 assert!(ranks_better_within_band(0.0, 0.0, -100.0, f64::NAN, band));
4117 assert!(!ranks_better_within_band(-100.0, f64::NAN, 0.0, 0.0, band));
4118 }
4119
4120 /// gh #534, the case the guard was stopping: `eigena2`'s dual infeasibility
4121 /// quarters on unit steps for four straight iterations, three short of a
4122 /// strict certificate. Quoted from the issue's own iteration table.
4123 #[test]
4124 fn eigena2_endgame_reads_as_contracting() {
4125 let eigena2 = [1.19e-05, 2.96e-06, 7.38e-07, 1.84e-07];
4126 assert!(window_is_contracting(
4127 &eigena2,
4128 DEFAULT_DECLINE_PROGRESS_RATIO
4129 ));
4130 }
4131
4132 /// gh #534, the case the guard was right about: `eigenb2`'s tail *rises*
4133 /// on heavily backtracked steps. The issue calls it a plausible genuine
4134 /// stall, so the progress test must refuse it and leave the guard alone.
4135 #[test]
4136 fn eigenb2_stall_does_not_read_as_contracting() {
4137 let eigenb2 = [1.88e-07, 2.69e-07, 2.89e-07, 2.93e-07];
4138 assert!(!window_is_contracting(
4139 &eigenb2,
4140 DEFAULT_DECLINE_PROGRESS_RATIO
4141 ));
4142 }
4143
4144 /// gh #534: `csfi2`'s window, measured on this build at the guard. Three
4145 /// healthy contractions and then a flat step — the solve has stopped
4146 /// moving, so the deferral must not fire however good the earlier ratios
4147 /// look. This is the shape every live guard firing reachable from the
4148 /// in-repo corpus has, which is why the whole window is tested and not
4149 /// just its first ratios.
4150 #[test]
4151 fn csfi2_flat_final_step_does_not_read_as_contracting() {
4152 let csfi2 = [3.267e0, 1.845e-6, 8.468e-8, 8.524e-8];
4153 assert!(!window_is_contracting(
4154 &csfi2,
4155 DEFAULT_DECLINE_PROGRESS_RATIO
4156 ));
4157 // ... and it is the *last* step that decides: drop it and the same
4158 // trace passes, which is exactly the distinction the test exists for.
4159 assert!(window_is_contracting(
4160 &csfi2[..3],
4161 DEFAULT_DECLINE_PROGRESS_RATIO
4162 ));
4163 }
4164
4165 /// gh #534: a large ratio drops the progress requirement, so the decline is
4166 /// deferred on any window. That is the "bypass the guard and see how far the
4167 /// solve gets" switch the issue asks for. A ratio of exactly `1` is the
4168 /// weaker "no backsliding" reading and still refuses `csfi2`, whose last
4169 /// step rises.
4170 #[test]
4171 fn a_large_ratio_accepts_a_stalled_window() {
4172 let csfi2 = [3.267e0, 1.845e-6, 8.468e-8, 8.524e-8];
4173 assert!(!window_is_contracting(&csfi2, 1.0));
4174 assert!(window_is_contracting(&[1e-8, 1e-8, 1e-8, 1e-8], 1.0));
4175 assert!(window_is_contracting(&csfi2, 1e20));
4176 // Still not a licence to read garbage as progress.
4177 assert!(!window_is_contracting(
4178 &[1.0, Number::NAN, 1e-9, 1e-12],
4179 1e20
4180 ));
4181 }
4182
4183 /// gh #534 edge cases: the ratio must never be evaluated against a
4184 /// non-positive or non-finite predecessor.
4185 #[test]
4186 fn degenerate_windows_never_read_as_contracting() {
4187 let r = DEFAULT_DECLINE_PROGRESS_RATIO;
4188 // A zero predecessor makes the ratio meaningless (0 <= 0.5*0 would
4189 // otherwise read as "contracting" forever).
4190 assert!(!window_is_contracting(&[0.0, 0.0, 0.0, 0.0], r));
4191 assert!(!window_is_contracting(&[1e-9, 0.0, 0.0, 0.0], r));
4192 assert!(!window_is_contracting(
4193 &[Number::INFINITY, 1e-3, 1e-6, 1e-9],
4194 r
4195 ));
4196 assert!(!window_is_contracting(&[1e-3, 1e-6, 1e-9, Number::NAN], r));
4197 // A genuine run down to exactly zero is progress, not a degenerate
4198 // window — the predecessor is positive at every step.
4199 assert!(window_is_contracting(&[1e-3, 1e-6, 1e-9, 0.0], r));
4200 }
4201
4202 /// gh #534: the window slides one sample per outer iteration and holds the
4203 /// most recent [`DECLINE_PROGRESS_SAMPLES`]. A short history is never a full
4204 /// window, which is what stops the first restoration entry of a solve from
4205 /// being deferred on no evidence at all — `nlp_err_contracting` requires
4206 /// `len == DECLINE_PROGRESS_SAMPLES` before it consults the samples.
4207 #[test]
4208 fn progress_window_slides_oldest_out() {
4209 let mut buf = [Number::NAN; DECLINE_PROGRESS_SAMPLES];
4210 let mut len = 0usize;
4211 for e in [1e-1, 1e-2, 1e-3] {
4212 push_sample(&mut buf, &mut len, e);
4213 }
4214 assert_eq!(len, 3);
4215 push_sample(&mut buf, &mut len, 1e-4);
4216 assert_eq!(len, DECLINE_PROGRESS_SAMPLES);
4217 assert_eq!(buf, [1e-1, 1e-2, 1e-3, 1e-4]);
4218 assert!(window_is_contracting(&buf, DEFAULT_DECLINE_PROGRESS_RATIO));
4219 // One flat iteration slides the oldest sample out and withdraws the
4220 // verdict.
4221 push_sample(&mut buf, &mut len, 1e-4);
4222 assert_eq!(len, DECLINE_PROGRESS_SAMPLES);
4223 assert_eq!(buf, [1e-2, 1e-3, 1e-4, 1e-4]);
4224 assert!(!window_is_contracting(&buf, DEFAULT_DECLINE_PROGRESS_RATIO));
4225 }
4226
4227 /// gh #505: no route may conclude `LocalInfeasibility` on its own.
4228 ///
4229 /// Three routes reach that verdict, and two of them independently shipped
4230 /// the same defect — building the terminate outcome directly, so the
4231 /// acceptable-point stash was never consulted and a good point the solve
4232 /// already had in hand was discarded. They were found one at a time,
4233 /// because nothing tied them together.
4234 ///
4235 /// The route a solve takes is an internal detail; the user sees one status
4236 /// either way. So what that status means is decided in one place — every
4237 /// route goes through [`IpoptAlgorithm::terminate_local_infeasibility`],
4238 /// or, for the cycle exits whose fallback is chosen between two statuses
4239 /// at the call site, through `terminate_acceptable_or`. Both consult the
4240 /// stash.
4241 ///
4242 /// **This is a tripwire, not a proof.** It is a substring scan of this
4243 /// file's source for the bare `IterateOutcome::Terminate(SolverReturn::
4244 /// LocalInfeasibility)` construction. A rustfmt line break through that
4245 /// expression, a `let` binding for the status, or a construction in
4246 /// another module all evade it — `application.rs` names the same variant
4247 /// on the SQP and ℓ₁ elastic paths and is deliberately out of scope. What
4248 /// it does buy is that the *obvious* way to add a fourth bare exit here
4249 /// fails loudly and points at the helper, which is the mistake that was
4250 /// actually made twice.
4251 ///
4252 /// The needle is assembled at runtime so this test's own source cannot
4253 /// satisfy the pattern it is checking for; an earlier version counted its
4254 /// own lines and failed against clean code.
4255 #[test]
4256 fn no_route_concludes_local_infeasibility_alone() {
4257 let needle = format!(
4258 "IterateOutcome::Terminate(SolverReturn::{})",
4259 "LocalInfeasibility"
4260 );
4261 let offenders: Vec<usize> = include_str!("ipopt_alg.rs")
4262 .lines()
4263 .enumerate()
4264 .filter(|(_, l)| {
4265 let t = l.trim_start();
4266 !t.starts_with("//") && !t.starts_with("///")
4267 })
4268 .filter(|(_, l)| l.contains(&needle))
4269 .map(|(i, _)| i + 1)
4270 .collect();
4271 assert!(
4272 offenders.is_empty(),
4273 "line(s) {offenders:?} build the local-infeasibility verdict directly. \
4274 Call `terminate_local_infeasibility()` instead — it consults the \
4275 acceptable-point stash first, so a solve that already passed through an \
4276 acceptable iterate returns that point rather than a hard failure. Two \
4277 routes shipped this bug before the helper existed (gh #505)."
4278 );
4279 }
4280}