pounce_algorithm/alg_builder.rs
1//! Algorithm builder — port of `Algorithm/IpAlgBuilder.{hpp,cpp}`.
2//!
3//! Reads `OptionsList`, walks the dependency order documented in
4//! `ref/Ipopt/AGENT_REFERENCE/ARCHITECTURE.md` §"BuildBasicAlgorithm",
5//! and assembles the strategy objects needed by `IpoptAlgorithm`:
6//!
7//! * `SymLinearSolver` (MA57 / MUMPS / FERAL) → `AugSystemSolver`
8//! (`StdAugSystemSolver`) → `PdSystemSolver` (`PdFullSpaceSolver`)
9//! → `SearchDirCalculator` (`PdSearchDirCalc`).
10//! * `BacktrackingLsAcceptor` (filter / penalty / cg-penalty) →
11//! `BacktrackingLineSearch`.
12//! * `MuUpdate` (monotone / adaptive[+oracle]).
13//! * `ConvCheck` (`OptErrorConvCheck`).
14//! * `IterateInitializer` (default / warm-start) and
15//! `EqMultCalculator` (`LeastSquareMults`).
16//! * `HessianUpdater` (exact / limited-memory).
17//! * `IterationOutput` (`OrigIterationOutput`).
18//! * `NLPScalingObject` (none / user / gradient-based / equilibration-based).
19//!
20//! Phase 7 ships the option-driven dispatch surface; the assembled
21//! `IpoptAlgorithm` lands once each strategy's arithmetic does.
22
23use crate::conv_check::opt_error::OptErrorConvCheck;
24use crate::eq_mult::least_square::LeastSquareMults;
25use crate::hess::exact::ExactHessianUpdater;
26use crate::hess::lim_mem_quasi_newton::{LimMemQuasiNewtonUpdater, UpdateType};
27use crate::init::default::DefaultIterateInitializer;
28use crate::init::warm_start::WarmStartIterateInitializer;
29use crate::kkt::aug_system_solver::AugSystemSolver;
30use crate::kkt::low_rank_aug_system_solver::LowRankAugSystemSolver;
31use crate::kkt::pd_full_space_solver::PdFullSpaceSolver;
32use crate::kkt::pd_search_dir_calc::PdSearchDirCalc;
33use crate::kkt::perturbation_handler::PdPerturbationHandler;
34use crate::kkt::std_aug_system_solver::StdAugSystemSolver;
35use crate::line_search::backtracking::BacktrackingLineSearch;
36use crate::line_search::filter_acceptor::FilterLsAcceptor;
37use crate::line_search::ls_acceptor::BacktrackingLsAcceptor;
38use crate::line_search::penalty_acceptor::PenaltyLsAcceptor;
39use crate::mu::adaptive::{AdaptiveMuUpdate, MuOracleKind};
40use crate::mu::monotone::MonotoneMuUpdate;
41use crate::output::orig::OrigIterationOutput;
42use pounce_common::types::{Index, Number};
43use pounce_linsol::{SparseSymLinearSolverInterface, TSymLinearSolver};
44use std::cell::RefCell;
45use std::rc::Rc;
46
47/// Backend factory — the application supplies one before calling
48/// [`AlgorithmBuilder::build`]. Mirrors upstream's
49/// `SymLinearSolverFactory` knob in `IpAlgBuilder.cpp`. The default
50/// factory wires in FERAL; MA57 is selectable when the `ma57` cargo
51/// feature is enabled.
52pub type LinearBackendFactory =
53 Box<dyn FnMut(LinearSolverChoice) -> Box<dyn SparseSymLinearSolverInterface>>;
54
55/// Top-level algorithm choice. `InteriorPoint` is pounce's default
56/// (the existing `IpoptAlgorithm`); `ActiveSetSqp` is the
57/// Phase 5b SQP driver in `crate::sqp::SqpAlgorithm`, which uses
58/// `pounce-qp` for QP subproblem solves and reuses
59/// `FilterLsAcceptor` for globalization.
60#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
61pub enum AlgorithmChoice {
62 #[default]
63 InteriorPoint,
64 ActiveSetSqp,
65}
66
67#[derive(Debug, Clone, Copy, PartialEq, Eq)]
68pub enum LinearSolverChoice {
69 Ma57,
70 Feral,
71}
72
73/// Symmetric scaling method applied to the augmented KKT system by
74/// [`TSymLinearSolver`]. Mirrors the `linear_system_scaling` option
75/// in `IpAlgBuilder.cpp:302-318` and the `RuizTSymScalingMethod` /
76/// `Mc19TSymScalingMethod` strategies in upstream Ipopt.
77///
78/// * `None` (default) — no scaling; `TSymLinearSolver` runs with a
79/// null scaling method. Matches upstream's default.
80/// * `Ruiz` — iterative symmetric ∞-norm equilibration (Ruiz, 2001).
81/// Implemented in `pounce_linsol::RuizTSymScalingMethod`.
82/// * `Mc19` — Curtis-Reid (HSL MC19) scaling. Not yet implemented;
83/// falls back to `None` with a warning.
84#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
85pub enum LinearSystemScalingChoice {
86 #[default]
87 None,
88 Ruiz,
89 Mc19,
90}
91
92#[derive(Debug, Clone, Copy, PartialEq, Eq)]
93pub enum MuStrategyChoice {
94 Monotone,
95 Adaptive,
96}
97
98#[derive(Debug, Clone, Copy, PartialEq, Eq)]
99pub enum HessianApproxChoice {
100 Exact,
101 LimitedMemory,
102}
103
104#[derive(Debug, Clone, Copy, PartialEq, Eq)]
105pub enum LineSearchChoice {
106 Filter,
107 CgPenalty,
108 Penalty,
109}
110
111/// Assembled strategy bundle. Phase 7 ships the structural bundle;
112/// `IpoptAlgorithm::new` reads from this when it lands.
113pub struct AlgorithmBundle {
114 pub mu_update: Box<dyn crate::mu::r#trait::MuUpdate>,
115 pub conv_check: Box<dyn crate::conv_check::r#trait::ConvCheck>,
116 pub init: Box<dyn crate::init::r#trait::IterateInitializer>,
117 pub eq_mult: Box<dyn crate::eq_mult::r#trait::EqMultCalculator>,
118 pub hess: Box<dyn crate::hess::r#trait::HessianUpdater>,
119 pub line_search: BacktrackingLineSearch,
120 pub iter_output: Box<dyn crate::output::r#trait::IterationOutput>,
121 /// `Some` when the builder was given a [`LinearBackendFactory`];
122 /// `None` for the bare structural bundle that pre-Phase-6 unit
123 /// tests still rely on.
124 pub search_dir: Option<PdSearchDirCalc>,
125}
126
127/// Knobs read off `OptionsList` and baked into the assembled
128/// `OptErrorConvCheck`. Defaults mirror
129/// `IpOptErrorConvCheck.cpp:RegisterOptions`.
130#[derive(Debug, Clone)]
131pub struct ConvCheckOptions {
132 pub tol: Number,
133 pub dual_inf_tol: Number,
134 pub constr_viol_tol: Number,
135 pub compl_inf_tol: Number,
136 pub acceptable_tol: Number,
137 pub acceptable_dual_inf_tol: Number,
138 pub acceptable_constr_viol_tol: Number,
139 pub acceptable_compl_inf_tol: Number,
140 pub acceptable_obj_change_tol: Number,
141 pub acceptable_iter: Index,
142 pub max_iter: Index,
143 pub max_cpu_time: Number,
144 pub max_wall_time: Number,
145 pub infeas_stationarity_tol: Number,
146 pub infeas_viol_kappa: Number,
147 pub infeas_max_streak: Index,
148 /// Objective-scale floor below which a strict termination certificate is
149 /// refused while the unscaled KKT error is still above `acceptable_tol`
150 /// (gh #200). `0` disables the mechanism.
151 pub obj_scale_certificate_threshold: Number,
152 /// Safety factor on the per-row floor the **strict** gate uses to decide
153 /// when a constraint residual is finer than the row can represent
154 /// (gh #528). `0` disables the floor, restoring upstream Ipopt's
155 /// bare-absolute primal term.
156 pub primal_noise_floor_kappa: Number,
157 /// Fraction of `acceptable_tol` the KKT error and the objective may drift
158 /// across the acceptable-level streak's window while the streak still
159 /// counts as settled (gh #533). `0` disables the progress test, leaving
160 /// acceptable-level termination the bare consecutive-count criterion.
161 pub acceptable_progress_kappa: Number,
162 /// Safety factor on the scale-relative floor under `dual_inf_tol` the
163 /// **strict** gate judges the dual infeasibility against (gh #532). `0`
164 /// disables the floor, restoring upstream Ipopt's bare-absolute bound.
165 pub dual_inf_scale_kappa: Number,
166}
167
168impl Default for ConvCheckOptions {
169 fn default() -> Self {
170 Self {
171 tol: 1e-8,
172 dual_inf_tol: 1.0,
173 constr_viol_tol: 1e-4,
174 compl_inf_tol: 1e-4,
175 acceptable_tol: 1e-6,
176 acceptable_dual_inf_tol: 1e10,
177 acceptable_constr_viol_tol: 1e-2,
178 acceptable_compl_inf_tol: 1e-2,
179 acceptable_obj_change_tol: 1e20,
180 acceptable_iter: 15,
181 max_iter: 3000,
182 max_cpu_time: 1e6,
183 max_wall_time: 1e6,
184 infeas_stationarity_tol: 1e-8,
185 infeas_viol_kappa: 1e2,
186 infeas_max_streak: 5,
187 obj_scale_certificate_threshold: 1e-4,
188 primal_noise_floor_kappa: 64.0,
189 acceptable_progress_kappa: 1e-1,
190 dual_inf_scale_kappa: 1.0,
191 }
192 }
193}
194
195#[derive(Debug, Clone)]
196pub struct AlgorithmBuilder {
197 /// Top-level algorithm dispatch. Default `InteriorPoint` ⇒
198 /// `build_with_backend` returns the existing `AlgorithmBundle`
199 /// (consumed by `IpoptAlgorithm`). `ActiveSetSqp` ⇒ caller
200 /// must use `build_sqp_with_backend` to assemble the Phase 5b
201 /// `SqpAlgorithm`. The two builder methods sit side by side
202 /// because the assembled algorithm shape differs (IPM bundle
203 /// vs SQP struct).
204 pub algorithm: AlgorithmChoice,
205 pub linear_solver: LinearSolverChoice,
206 /// Symmetric scaling method for the augmented KKT system. Wired
207 /// into [`TSymLinearSolver`] by [`Self::build_with_backend`].
208 /// Mirrors upstream `linear_system_scaling` (`IpAlgBuilder.cpp:538-560`).
209 pub linear_system_scaling: LinearSystemScalingChoice,
210 /// Lazy-vs-eager scaling toggle (`linear_scaling_on_demand`,
211 /// `IpTSymLinearSolver.cpp:50-58`). Only consulted when
212 /// `linear_system_scaling != None`. Upstream default is `true`
213 /// (compute scaling only on the first solve that fails / shows
214 /// poor conditioning); pounce mirrors that. Set to `false` to
215 /// scale every factorization.
216 pub linear_scaling_on_demand: bool,
217 pub mu_strategy: MuStrategyChoice,
218 /// Selector forwarded to [`AdaptiveMuUpdate`] when
219 /// `mu_strategy = Adaptive`. Ignored for `Monotone`. Defaults to
220 /// `QualityFunction` per upstream's `RegisterOptions` default.
221 pub mu_oracle: MuOracleKind,
222 pub hessian_approximation: HessianApproxChoice,
223 pub limited_memory_update_type: UpdateType,
224 /// History length for the limited-memory quasi-Newton approximation
225 /// (`limited_memory_max_history`). Defaults to upstream's 6.
226 pub limited_memory_max_history: i32,
227 /// `limited_memory_init_val_max` / `_min` — the clamp on the initial
228 /// Hessian scalar σ before the rank-2 updates. Upstream defaults 1e8
229 /// / 1e-8, which `LimMemQuasiNewtonUpdater` has carried as hard-coded
230 /// fields and consumed in `initial_hessian_scalar` all along; only
231 /// the read sites were missing (gh#483, #191 round 2).
232 pub limited_memory_init_val_max: Number,
233 pub limited_memory_init_val_min: Number,
234 pub line_search_method: LineSearchChoice,
235 pub warm_start_init_point: bool,
236 /// `mehrotra_algorithm` — when true, [`PdSearchDirCalc`] folds
237 /// the Mehrotra second-order complementarity term into the
238 /// search-direction RHS. Mirrors upstream's
239 /// `IpAlgBuilder.cpp:Mehrotra` flag. Requires `mu_strategy =
240 /// Adaptive` so that an affine step is computed each iteration;
241 /// [`Self::build_with_backend`] does not enforce this — the
242 /// option-parser in `application.rs` is responsible for the
243 /// cascading defaults (`mu_oracle = probing` etc.).
244 pub mehrotra_algorithm: bool,
245 /// `fast_step_computation` — when true, [`PdSearchDirCalc`] accepts
246 /// the search direction without the residual check and allows an
247 /// inexact linear solve. Mirrors upstream's flag of the same name,
248 /// default `no`. The field existed and was consumed from the day the
249 /// search-direction calculator landed, hard-coded to `false`; only
250 /// the option's read site was missing, so setting it did nothing
251 /// (gh#483 follow-up, #191 round 2).
252 pub fast_step_computation: bool,
253 /// `kappa_sigma` — factor bounding how far the bound multipliers may
254 /// deviate from their primal estimates. The clamp
255 /// (`kappa_sigma_clamp`) runs after every accepted step; `< 1`
256 /// disables the correction. Mirrors `IpIpoptAlg.cpp` (Eqn. (16)),
257 /// default `1e10`. Baked onto [`crate::ipopt_alg::IpoptAlgorithm`] by
258 /// the solve path.
259 pub kappa_sigma: Number,
260 /// `kappa_d` — weight of the linear damping term added to the barrier
261 /// objective/gradient (and dual-infeasibility) to handle one-sided
262 /// bounds. Mirrors `IpIpoptCalculatedQuantities.cpp`, default `1e-5`.
263 /// Baked onto [`crate::ipopt_cq::IpoptCalculatedQuantities`] by the
264 /// solve path.
265 pub kappa_d: Number,
266 /// `tiny_step_tol` — relative primal step size below which the full
267 /// step is accepted without line search; repeated tiny steps
268 /// terminate the solve. Mirrors `IpBacktrackingLineSearch.cpp`,
269 /// default `10·EPSILON`. Baked onto
270 /// [`crate::ipopt_alg::IpoptAlgorithm`] by the solve path.
271 pub tiny_step_tol: Number,
272 /// `tiny_step_y_tol` — dual-step threshold; when both primal and dual
273 /// steps are tiny in consecutive iterations the algorithm stops at the
274 /// best attainable accuracy. Default `1e-2`.
275 pub tiny_step_y_tol: Number,
276 /// `diverging_iterates_tol` — if `max_i |x_i|` exceeds this the solve
277 /// aborts as diverging. Default `1e20`.
278 pub diverging_iterates_tol: Number,
279 /// `dual_diverging_streak` (pounce#246) — consecutive growing-dual-
280 /// infeasibility iterations before the dual-divergence guard routes to
281 /// restoration. **Default `0` (off).**
282 ///
283 /// It defaulted to `15` when introduced, on the strength of a reported
284 /// emfl050 bad-warm-start grind. That justification did not survive being
285 /// reproduced: the measurement was caller-side JAX compilation, and the
286 /// build predating the guard solves both emfl050 instances to the same
287 /// optimum in the same time (pounce#246 / pounce#250). What remained was a
288 /// knife-edge, non-monotone effect on four of 1284 MINLPLib models — so it
289 /// is opt-in rather than imposed. See `upstream_options.rs` for the full
290 /// account.
291 pub dual_diverging_streak: Index,
292 /// `resto_decline_deferrals` (gh #534) — how many times the
293 /// acceptable-point restoration decline may be deferred on a solve whose
294 /// NLP error is still contracting. Default `1`; `0` restores the pre-#534
295 /// behaviour (decline immediately, always). See `upstream_options.rs`.
296 pub resto_decline_deferrals: Index,
297 /// `resto_decline_progress_ratio` (gh #534) — required per-iteration
298 /// contraction of the NLP error before a decline is deferred. Default
299 /// `0.5`; at or above `1` the progress requirement is dropped entirely.
300 pub resto_decline_progress_ratio: Number,
301 /// `kkt_fidelity_tol` (pounce#173). Read by the algorithm as well as by the
302 /// post-solve gate, because the #200 fallback's tiebreak has to rank the two
303 /// candidate points by the status each will be *reported* under. Default
304 /// `0.0` (gate disabled).
305 pub kkt_fidelity_tol: Number,
306 pub conv_check: ConvCheckOptions,
307 pub mu: MuOptions,
308 pub line_search: LineSearchOptions,
309 pub refinement: RefinementOptions,
310 pub perturbation: PerturbationOptions,
311 pub resto: RestoOptions,
312 pub output: OutputOptions,
313 pub warm: WarmStartOptions,
314 /// SQP-specific options (consulted only when
315 /// `algorithm = ActiveSetSqp`).
316 pub sqp: crate::sqp::SqpOptions,
317 /// QP-subproblem-solver options for the active-set SQP path
318 /// (`pounce_qp::QpOptions`), threaded into the `SqpAlgorithm` via
319 /// `with_qp_options`. Consulted only when `algorithm = ActiveSetSqp`.
320 /// Populated from the `sqp_qp_*` CLI options by
321 /// `application::apply_qp_subproblem_options`.
322 pub sqp_qp: pounce_qp::QpOptions,
323 pub init: InitOptions,
324 /// Optional block-triangular / Schur KKT partition (pounce#180 item 2):
325 /// `(schur_indices, feral_cfg)`. When `Some` and the IPM path is selected
326 /// with the feral linear solver and an exact Hessian, `build_with_backend`
327 /// wraps the standard aug-system solver in a
328 /// [`crate::kkt::SchurAugSystemSolver`] over the given KKT-space indices.
329 /// The Schur solver falls back to the standard solver transparently when
330 /// the partition is unsuitable. Set via [`Self::set_kkt_schur`].
331 pub kkt_schur: Option<(Vec<usize>, pounce_feral::FeralConfig)>,
332}
333
334/// Knobs read off `OptionsList` and baked into
335/// [`DefaultIterateInitializer`]. Defaults mirror
336/// `IpDefaultIterateInitializer.cpp:RegisterOptions`. The Mehrotra
337/// cascade in `application.rs` overrides `bound_push`, `bound_frac`,
338/// and `bound_mult_init_val` to upstream's more-aggressive values
339/// (`10`, `0.2`, `1.0`).
340#[derive(Debug, Clone)]
341pub struct InitOptions {
342 pub bound_push: Number,
343 pub bound_frac: Number,
344 pub slack_bound_push: Number,
345 pub slack_bound_frac: Number,
346 pub constr_mult_init_max: Number,
347 pub bound_mult_init_val: Number,
348 /// `bound_mult_init_method`: `"constant"` (default) or `"mu-based"`
349 /// (matches upstream's `IpDefaultIterateInitializer.cpp`).
350 pub bound_mult_init_method: String,
351 /// `least_square_init_primal` — replace the user's starting `x`
352 /// with the min-norm primal that satisfies the linearized
353 /// constraints. Used by the Mehrotra cascade in `application.rs`
354 /// to drop iter-0 primal infeasibility on LP-shaped problems.
355 /// Mirrors upstream `IpDefaultIterateInitializer.cpp:200-222`.
356 pub least_square_init_primal: bool,
357}
358
359impl Default for InitOptions {
360 fn default() -> Self {
361 Self {
362 bound_push: 1e-2,
363 bound_frac: 1e-2,
364 slack_bound_push: 1e-2,
365 slack_bound_frac: 1e-2,
366 constr_mult_init_max: 1e3,
367 bound_mult_init_val: 1.0,
368 bound_mult_init_method: "constant".into(),
369 least_square_init_primal: false,
370 }
371 }
372}
373
374/// Knobs read off `OptionsList` and baked into
375/// [`WarmStartIterateInitializer`]. Defaults mirror
376/// `IpWarmStartIterateInitializer.cpp:RegisterOptions`.
377///
378/// Wired today: `mult_init_max` (clamps |y_c|, |y_d| and caps z/v
379/// blocks) and `target_mu` (overrides `data.curr_mu` at iter 0).
380/// The remaining knobs (`bound_push`, `bound_frac`, `slack_bound_push`,
381/// `slack_bound_frac`, `mult_bound_push`, `entire_iterate`,
382/// `same_structure`) are stored on the initializer but not yet
383/// consumed — `WarmStartIterateInitializer::set_initial_iterates`
384/// currently trusts the caller-populated `data.curr` rather than
385/// re-running the upstream `push_variables` machinery.
386#[derive(Debug, Clone)]
387pub struct WarmStartOptions {
388 pub bound_push: Number,
389 pub bound_frac: Number,
390 pub slack_bound_push: Number,
391 pub slack_bound_frac: Number,
392 pub mult_bound_push: Number,
393 pub mult_init_max: Number,
394 pub target_mu: Number,
395 pub entire_iterate: bool,
396 pub same_structure: bool,
397 /// The value a NaN-seeded bound multiplier takes: NaN in a
398 /// user-supplied `z`/`v` seed means "unseeded, use the default".
399 /// Threaded from `builder.init.bound_mult_init_val` at build time
400 /// so the Mehrotra override and any user setting stay the single
401 /// source of truth.
402 pub bound_mult_init_val: Number,
403}
404
405impl Default for WarmStartOptions {
406 fn default() -> Self {
407 Self {
408 bound_push: 1e-3,
409 bound_frac: 1e-3,
410 slack_bound_push: 1e-3,
411 slack_bound_frac: 1e-3,
412 mult_bound_push: 1e-3,
413 mult_init_max: 1e6,
414 target_mu: 0.0,
415 entire_iterate: false,
416 same_structure: false,
417 // seeded from the init options so the default has one
418 // home; build() re-resolves it from the live init options
419 // anyway (see `resolved_warm_options`)
420 bound_mult_init_val: InitOptions::default().bound_mult_init_val,
421 }
422 }
423}
424
425/// The warm-start options as the initializer actually receives them:
426/// `bound_mult_init_val` comes from the (option-read,
427/// Mehrotra-resolved) init options, never from `WarmStartOptions`'s
428/// own copy. Split out of `build()` so the threading is testable.
429pub(crate) fn resolved_warm_options(
430 warm: &WarmStartOptions,
431 init: &InitOptions,
432) -> WarmStartOptions {
433 let mut w = warm.clone();
434 w.bound_mult_init_val = init.bound_mult_init_val;
435 w
436}
437
438/// Knobs read off `OptionsList` and baked into the assembled
439/// `MonotoneMuUpdate` or `AdaptiveMuUpdate`. Defaults mirror
440/// `IpMonotoneMuUpdate.cpp` / `IpAdaptiveMuUpdate.cpp:RegisterOptions`.
441/// `mu_max` defaults to the sentinel `-1`; positive values are baked
442/// into both updaters at build time (adaptive interprets `-1` as
443/// "lazy-init from `mu_max_fact * avrg_compl`").
444#[derive(Debug, Clone)]
445pub struct MuOptions {
446 pub mu_init: Number,
447 pub mu_max: Number,
448 pub mu_max_fact: Number,
449 pub mu_min: Number,
450 pub mu_target: Number,
451 pub mu_linear_decrease_factor: Number,
452 pub mu_superlinear_decrease_power: Number,
453 pub mu_allow_fast_monotone_decrease: bool,
454 pub barrier_tol_factor: Number,
455 /// `sigma_max` / `sigma_min` — clamp on the centering parameter σ
456 /// chosen by `QualityFunctionMuOracle`. Only consumed when
457 /// `mu_strategy=adaptive` and `mu_oracle=quality-function`.
458 /// Defaults from `IpQualityFunctionMuOracle.cpp:RegisterOptions`.
459 pub sigma_max: Number,
460 pub sigma_min: Number,
461 /// `adaptive_mu_globalization` — globalization strategy for the
462 /// adaptive μ-selection mode. Mirrors
463 /// `IpAdaptiveMuUpdate.cpp:RegisterOptions`. Default is
464 /// `ObjConstrFilter`; the Mehrotra cascade switches to
465 /// `NeverMonotoneMode` to disable globalization entirely.
466 pub adaptive_mu_globalization: crate::mu::adaptive::AdaptiveMuGlobalization,
467 /// `quality_function_norm_type` — norm used inside the quality
468 /// function to aggregate the three KKT components. Forwarded to
469 /// `QualityFunctionMuOracle` when `mu_oracle=quality-function`.
470 pub quality_function_norm_type: crate::mu::oracle::quality_function::NormType,
471 /// `quality_function_centrality` — centrality penalty term added
472 /// to the quality function.
473 pub quality_function_centrality: crate::mu::oracle::quality_function::CentralityType,
474 /// `quality_function_balancing_term` — balancing penalty term in
475 /// the quality function (kicks in when complementarity is far
476 /// below infeasibilities).
477 pub quality_function_balancing_term: crate::mu::oracle::quality_function::BalancingTermType,
478 /// `quality_function_max_section_steps` — cap on golden-section
479 /// iterations when picking σ. Default 8.
480 pub quality_function_max_section_steps: i32,
481 /// `quality_function_section_sigma_tol` — width tolerance in
482 /// σ-space for golden section. Default 1e-2.
483 pub quality_function_section_sigma_tol: Number,
484 /// `quality_function_section_qf_tol` — relative flatness
485 /// tolerance for golden section. Default 0.0.
486 pub quality_function_section_qf_tol: Number,
487 /// `adaptive_mu_safeguard_factor` — guard for the LOQO fallback
488 /// in adaptive mode. Default 0.0.
489 pub adaptive_mu_safeguard_factor: Number,
490 /// `adaptive_mu_monotone_init_factor` — multiplier on the
491 /// average complementarity when seeding monotone mode after a
492 /// free-mode bailout. Default 0.8.
493 pub adaptive_mu_monotone_init_factor: Number,
494 /// `adaptive_mu_restore_previous_iterate` — restore the most
495 /// recent free-mode iterate when switching to fixed mode.
496 /// Default `false`.
497 pub adaptive_mu_restore_previous_iterate: bool,
498 /// `adaptive_mu_kkterror_red_iters` — window length for the
499 /// `KKT_ERROR` globalization history. Default 4.
500 pub adaptive_mu_kkterror_red_iters: usize,
501 /// `adaptive_mu_kkterror_red_fact` — required relative reduction
502 /// of the KKT error over the window. Default 0.9999.
503 pub adaptive_mu_kkterror_red_fact: Number,
504 /// `adaptive_mu_kkt_norm_type` — norm used to score the iterate
505 /// in adaptive globalization decisions.
506 pub adaptive_mu_kkt_norm_type: crate::mu::adaptive::AdaptiveMuKktNorm,
507 /// `probing_iterate_quality_factor` (default 1e4, pounce-specific
508 /// — see pounce#58). When the probing (Mehrotra) μ-oracle is
509 /// about to read `curr_avrg_compl()` for its `mu_curr` input, a
510 /// single imbalanced `(s_i, z_i)` pair can inflate the average
511 /// 5+ orders above the stored `data.curr_mu`. The oracle then
512 /// returns `σ · mu_curr` ≫ previous μ, throwing the iterate out
513 /// of the convergence neighborhood. This guard short-circuits
514 /// that case by signalling restoration when the ratio
515 /// `curr_avrg_compl / curr_mu` exceeds the factor. Set to 0 or
516 /// any non-positive value to disable.
517 pub probing_iterate_quality_factor: Number,
518}
519
520impl Default for MuOptions {
521 fn default() -> Self {
522 Self {
523 mu_init: 0.1,
524 mu_max: -1.0,
525 mu_max_fact: 1e3,
526 mu_min: 1e-11,
527 mu_target: 0.0,
528 mu_linear_decrease_factor: 0.2,
529 mu_superlinear_decrease_power: 1.5,
530 mu_allow_fast_monotone_decrease: true,
531 barrier_tol_factor: 10.0,
532 sigma_max: 1e2,
533 sigma_min: 1e-6,
534 adaptive_mu_globalization:
535 crate::mu::adaptive::AdaptiveMuGlobalization::ObjConstrFilter,
536 quality_function_norm_type:
537 crate::mu::oracle::quality_function::NormType::TwoNormSquared,
538 quality_function_centrality: crate::mu::oracle::quality_function::CentralityType::None,
539 quality_function_balancing_term:
540 crate::mu::oracle::quality_function::BalancingTermType::None,
541 quality_function_max_section_steps: 8,
542 quality_function_section_sigma_tol: 1e-2,
543 quality_function_section_qf_tol: 0.0,
544 adaptive_mu_safeguard_factor: 0.0,
545 adaptive_mu_monotone_init_factor: 0.8,
546 adaptive_mu_restore_previous_iterate: false,
547 adaptive_mu_kkterror_red_iters: 4,
548 adaptive_mu_kkterror_red_fact: 0.9999,
549 adaptive_mu_kkt_norm_type: crate::mu::adaptive::AdaptiveMuKktNorm::TwoNormSquared,
550 probing_iterate_quality_factor: 1e4,
551 }
552 }
553}
554
555/// Knobs baked into the assembled [`BacktrackingLineSearch`]. Defaults
556/// mirror `IpBacktrackingLineSearch.cpp:RegisterOptions`.
557#[derive(Debug, Clone)]
558pub struct LineSearchOptions {
559 pub watchdog_shortened_iter_trigger: Index,
560 pub watchdog_trial_iter_max: Index,
561 /// `soft_resto_pderror_reduction_factor` — required relative
562 /// reduction in the primal-dual error for a soft-resto step.
563 /// `0` disables the soft restoration phase.
564 pub soft_resto_pderror_reduction_factor: Number,
565 /// `max_soft_resto_iters` — cap on consecutive soft-resto
566 /// iterations before full restoration is forced.
567 pub max_soft_resto_iters: Index,
568 /// `accept_every_trial_step` — short-circuits the filter / alpha
569 /// loop and accepts the full fraction-to-the-boundary step every
570 /// outer iteration. Mirrors upstream's
571 /// `IpBacktrackingLineSearch::accept_every_trial_step_`. Drops
572 /// global convergence guarantees; only safe for problems where the
573 /// Newton step is already a descent step (LPs, convex QPs). The
574 /// Mehrotra cascade in `application.rs` flips this on.
575 pub accept_every_trial_step: bool,
576 /// `alpha_for_y` — policy for the equality-multiplier (y_c / y_d)
577 /// step length. Upstream default is `Primal`; the Mehrotra cascade
578 /// switches to `BoundMult`.
579 pub alpha_for_y: crate::line_search::backtracking::AlphaForY,
580
581 // Filter switching / Armijo / margin constants baked onto the
582 // assembled [`crate::line_search::filter_acceptor::FilterLsAcceptor`]
583 // (only when `line_search_method = Filter`). All were registered but
584 // never read (#191); defaults mirror `IpFilterLSAcceptor.cpp`.
585 /// `eta_phi` — relaxation factor in the Armijo condition (Eqn. (20)).
586 pub eta_phi: Number,
587 /// `theta_min_fact` — constraint-violation threshold factor in the
588 /// switching rule.
589 pub theta_min_fact: Number,
590 /// `theta_max_fact` — upper-bound factor for constraint violation in
591 /// the filter (Eqn. (21)).
592 pub theta_max_fact: Number,
593 /// `theta_max_row_scale_kappa` — multiplier on the constraint-row
594 /// count used as the floor of the `theta_max` reference.
595 /// **Opt-in**: default `0`, which is upstream's bare
596 /// `max(1, theta_0)` floor bit-for-bit. Set to `1` on a large model
597 /// that stalls from a feasible start. See
598 /// [`FilterLsAcceptor::theta_max_row_scale_kappa`].
599 pub theta_max_row_scale_kappa: Number,
600 /// `theta_max_adaptive_trigger` — consecutive line searches whose
601 /// every trial was refused at the `theta_max` gate before the
602 /// ceiling is raised. `0` disables the rule. See
603 /// [`FilterLsAcceptor::theta_max_adaptive_trigger`] (pounce#546).
604 pub theta_max_adaptive_trigger: u32,
605 /// Geometric factor applied to `theta_max` on each adaptive raise.
606 /// See [`FilterLsAcceptor::theta_max_adaptive_factor`].
607 pub theta_max_adaptive_factor: Number,
608 /// Cap on adaptive raises per solve, which is what keeps `theta_max`
609 /// finite. See [`FilterLsAcceptor::theta_max_adaptive_max_raises`].
610 pub theta_max_adaptive_max_raises: u32,
611 /// `gamma_phi` — filter margin factor for the barrier function
612 /// (Eqn. (18a)).
613 pub gamma_phi: Number,
614 /// `gamma_theta` — filter margin factor for the constraint violation
615 /// (Eqn. (18b)).
616 pub gamma_theta: Number,
617 /// `s_phi` — exponent for the linear barrier model in the switching
618 /// rule (Eqn. (19)).
619 pub s_phi: Number,
620 /// `s_theta` — exponent for the current constraint violation in the
621 /// switching rule (Eqn. (19)).
622 pub s_theta: Number,
623 /// `alpha_min_frac` — safety factor for the minimal step size before
624 /// switching to restoration (gamma_alpha, Eqn. (23)).
625 pub alpha_min_frac: Number,
626 /// `obj_max_inc` — max acceptable increase (orders of magnitude) of
627 /// the barrier objective for a trial point.
628 pub obj_max_inc: Number,
629 /// `max_filter_resets` — maximum number of filter resets allowed
630 /// (`0` disables the reset heuristic).
631 pub max_filter_resets: Index,
632 /// `filter_reset_trigger` — successive filter-rejected iterations that
633 /// trigger a filter reset.
634 pub filter_reset_trigger: Index,
635
636 // Second-order-correction constants baked onto the assembled
637 // [`BacktrackingLineSearch`]. Registered but never read (#191);
638 // defaults mirror `IpBacktrackingLineSearch.cpp`.
639 /// `max_soc` — max second-order-correction trial steps per iteration;
640 /// `0` disables SOC.
641 pub max_soc: Index,
642 /// `kappa_soc` — sufficient-reduction factor for a SOC step to be
643 /// continued.
644 pub kappa_soc: Number,
645 /// `soc_method` — `0` (paper method) or `1` (alpha-on-rhs variant).
646 pub soc_method: Index,
647}
648
649impl Default for LineSearchOptions {
650 fn default() -> Self {
651 Self {
652 watchdog_shortened_iter_trigger: 10,
653 watchdog_trial_iter_max: 3,
654 soft_resto_pderror_reduction_factor: 1.0 - 1e-4,
655 max_soft_resto_iters: 10,
656 accept_every_trial_step: false,
657 alpha_for_y: crate::line_search::backtracking::AlphaForY::Primal,
658 eta_phi: 1e-8,
659 theta_min_fact: 1e-4,
660 theta_max_fact: 1e4,
661 theta_max_row_scale_kappa: 0.0,
662 theta_max_adaptive_trigger: 3,
663 theta_max_adaptive_factor: 100.0,
664 theta_max_adaptive_max_raises: 4,
665 gamma_phi: 1e-8,
666 gamma_theta: 1e-5,
667 s_phi: 2.3,
668 s_theta: 1.1,
669 alpha_min_frac: 0.05,
670 obj_max_inc: 5.0,
671 max_filter_resets: 5,
672 filter_reset_trigger: 5,
673 max_soc: 4,
674 kappa_soc: 0.99,
675 soc_method: 0,
676 }
677 }
678}
679
680/// Inertia-correction / regularization knobs baked onto the assembled
681/// [`crate::kkt::perturbation_handler::PdPerturbationHandler`]. Field
682/// names use the option names; they map to the handler's `delta_xs_*` /
683/// `delta_cd_*` fields. Defaults mirror
684/// `IpPDPerturbationHandler.cpp:RegisterOptions`. All were registered but
685/// never read (#191).
686#[derive(Debug, Clone)]
687pub struct PerturbationOptions {
688 /// `max_hessian_perturbation` → `delta_xs_max`.
689 pub max_hessian_perturbation: Number,
690 /// `min_hessian_perturbation` → `delta_xs_min`.
691 pub min_hessian_perturbation: Number,
692 /// `perturb_inc_fact_first` → `delta_xs_first_inc_fact`.
693 pub perturb_inc_fact_first: Number,
694 /// `perturb_inc_fact` → `delta_xs_inc_fact`.
695 pub perturb_inc_fact: Number,
696 /// `perturb_dec_fact` → `delta_xs_dec_fact`.
697 pub perturb_dec_fact: Number,
698 /// `first_hessian_perturbation` → `delta_xs_init`.
699 pub first_hessian_perturbation: Number,
700 /// `jacobian_regularization_value` → `delta_cd_val`.
701 pub jacobian_regularization_value: Number,
702 /// `jacobian_regularization_exponent` → `delta_cd_exp`.
703 pub jacobian_regularization_exponent: Number,
704 /// `perturb_always_cd` — always regularize the c/d (Jacobian) block.
705 pub perturb_always_cd: bool,
706}
707
708impl Default for PerturbationOptions {
709 fn default() -> Self {
710 Self {
711 max_hessian_perturbation: 1e20,
712 min_hessian_perturbation: 1e-20,
713 perturb_inc_fact_first: 100.0,
714 perturb_inc_fact: 8.0,
715 perturb_dec_fact: 1.0 / 3.0,
716 first_hessian_perturbation: 1e-4,
717 jacobian_regularization_value: 1e-8,
718 jacobian_regularization_exponent: 0.25,
719 perturb_always_cd: false,
720 }
721 }
722}
723
724/// Restoration-phase knobs carried on the outer builder and copied into
725/// the `RestoAlgorithmBuilder` when the restoration factory is minted
726/// (`pounce-restoration`). The restoration builder is constructed with
727/// defaults by each frontend and never options-configured, so these were
728/// registered but never read (#191). Defaults mirror upstream's
729/// restoration `RegisterOptions`.
730#[derive(Debug, Clone)]
731pub struct RestoOptions {
732 /// `bound_mult_reset_threshold` — reset bound multipliers to 1 after
733 /// restoration if the largest exceeds this.
734 pub bound_mult_reset_threshold: Number,
735 /// `constr_mult_reset_threshold` — ignore the least-square constraint
736 /// multiplier estimate after restoration if its norm exceeds this
737 /// (`0` keeps the estimate).
738 pub constr_mult_reset_threshold: Number,
739 /// `resto_penalty_parameter` — penalty on the slack 1-norm in the
740 /// restoration objective (`rho`).
741 pub resto_penalty_parameter: Number,
742 /// `resto_proximity_weight` — proximity-term weight (`eta_factor`;
743 /// `η = eta_factor · sqrt(μ)`).
744 pub resto_proximity_weight: Number,
745 /// `required_infeasibility_reduction` — the restoration sub-solve
746 /// keeps iterating until the *original* NLP's infeasibility has been
747 /// reduced to at most this fraction of its value at restoration entry
748 /// (`κ_resto` in `IpRestoConvCheck.cpp:58`). `0` disables the guard,
749 /// i.e. restoration runs until the sub-NLP itself converges.
750 pub required_infeasibility_reduction: Number,
751 /// `evaluate_orig_obj_at_resto_trial` — evaluate the *original*
752 /// objective at every restoration trial point, so an iterate the
753 /// restoration problem likes but the original cannot evaluate is
754 /// rejected there rather than after the phase exits. Upstream default
755 /// `yes`. `RestoAlgorithmBuilder` has consumed this since it landed;
756 /// only the read site was missing (gh#483, #191 round 2).
757 pub evaluate_orig_obj_at_resto_trial: bool,
758 /// `expect_infeasible_problem` — enter restoration sooner and demand
759 /// more infeasibility reduction before leaving it. Upstream default
760 /// `no`. Same story: consumed, never read.
761 pub expect_infeasible_problem: bool,
762 /// `start_with_resto` — switch to restoration in the first iteration.
763 /// Upstream default `no`. Same story.
764 pub start_with_resto: bool,
765}
766
767impl Default for RestoOptions {
768 fn default() -> Self {
769 Self {
770 bound_mult_reset_threshold: 1e3,
771 constr_mult_reset_threshold: 0.0,
772 resto_penalty_parameter: 1e3,
773 resto_proximity_weight: 1.0,
774 required_infeasibility_reduction: 0.9,
775 evaluate_orig_obj_at_resto_trial: true,
776 expect_infeasible_problem: false,
777 start_with_resto: false,
778 }
779 }
780}
781
782/// Iterative-refinement knobs baked onto the assembled
783/// [`crate::kkt::pd_full_space_solver::PdFullSpaceSolver`]. Defaults
784/// mirror `IpPDFullSpaceSolver.cpp:RegisterOptions`. All were registered
785/// but never read (#191).
786#[derive(Debug, Clone)]
787pub struct RefinementOptions {
788 /// `min_refinement_steps` — minimum iterative-refinement steps per
789 /// linear solve.
790 pub min_refinement_steps: Index,
791 /// `max_refinement_steps` — maximum iterative-refinement steps.
792 pub max_refinement_steps: Index,
793 /// `residual_ratio_max` — refine until the residual test ratio drops
794 /// below this (or `max_refinement_steps` is reached).
795 pub residual_ratio_max: Number,
796 /// `residual_ratio_singular` — above this ratio after failed
797 /// refinement, the system is declared singular.
798 pub residual_ratio_singular: Number,
799 /// `residual_improvement_factor` — minimum per-step reduction of the
800 /// residual test ratio before refinement is aborted.
801 pub residual_improvement_factor: Number,
802}
803
804impl Default for RefinementOptions {
805 fn default() -> Self {
806 Self {
807 min_refinement_steps: 1,
808 max_refinement_steps: 10,
809 residual_ratio_max: 1e-10,
810 residual_ratio_singular: 1e-5,
811 residual_improvement_factor: 0.999_999_999,
812 }
813 }
814}
815
816/// Knobs baked into the assembled [`OrigIterationOutput`]. Defaults
817/// mirror `IpOrigIterationOutput.cpp:RegisterOptions` /
818/// `IpAlgorithmRegOp.cpp`.
819#[derive(Debug, Clone)]
820pub struct OutputOptions {
821 pub print_frequency_iter: Index,
822 pub print_frequency_time: Number,
823 /// `print_info_string` (default `false`). When on, the iter row
824 /// ends with the contents of `IpoptData::info_string` so users
825 /// can read the per-iteration diagnostic tags.
826 pub print_info_string: bool,
827 /// `inf_pr_output` — `"original"` (default) prints the unscaled
828 /// NLP primal infeasibility; `"internal"` prints the internal
829 /// reformulated violation. Only meaningful once NLP-side scaling
830 /// is in play; until then both modes produce the same number.
831 pub inf_pr_output_internal: bool,
832}
833
834impl Default for OutputOptions {
835 fn default() -> Self {
836 Self {
837 print_frequency_iter: 1,
838 print_frequency_time: 0.0,
839 print_info_string: false,
840 inf_pr_output_internal: false,
841 }
842 }
843}
844
845impl Default for AlgorithmBuilder {
846 fn default() -> Self {
847 Self {
848 algorithm: AlgorithmChoice::default(),
849 linear_solver: LinearSolverChoice::Feral,
850 linear_system_scaling: LinearSystemScalingChoice::None,
851 linear_scaling_on_demand: true,
852 mu_strategy: MuStrategyChoice::Monotone,
853 mu_oracle: MuOracleKind::QualityFunction,
854 hessian_approximation: HessianApproxChoice::Exact,
855 limited_memory_update_type: UpdateType::Bfgs,
856 limited_memory_max_history: 6,
857 limited_memory_init_val_max: 1e8,
858 limited_memory_init_val_min: 1e-8,
859 line_search_method: LineSearchChoice::Filter,
860 warm_start_init_point: false,
861 mehrotra_algorithm: false,
862 fast_step_computation: false,
863 kappa_sigma: 1e10,
864 kappa_d: 1e-5,
865 tiny_step_tol: 10.0 * Number::EPSILON,
866 tiny_step_y_tol: 1e-2,
867 diverging_iterates_tol: 1e20,
868 dual_diverging_streak: 0,
869 resto_decline_deferrals: 1,
870 resto_decline_progress_ratio: 0.5,
871 kkt_fidelity_tol: 0.0,
872 conv_check: ConvCheckOptions::default(),
873 mu: MuOptions::default(),
874 line_search: LineSearchOptions::default(),
875 refinement: RefinementOptions::default(),
876 perturbation: PerturbationOptions::default(),
877 resto: RestoOptions::default(),
878 output: OutputOptions::default(),
879 warm: WarmStartOptions::default(),
880 sqp: crate::sqp::SqpOptions::default(),
881 sqp_qp: pounce_qp::QpOptions::default(),
882 init: InitOptions::default(),
883 kkt_schur: None,
884 }
885 }
886}
887
888impl AlgorithmBuilder {
889 pub fn new() -> Self {
890 Self::default()
891 }
892
893 /// Install a Schur KKT partition (pounce#180 item 2). `schur_indices` are
894 /// KKT-space indices (`0..dim`, the `x,s,c,d` block order the aug-system
895 /// solver assembles); `cfg` configures the per-block feral solvers. Only
896 /// honored on the IPM + feral + exact-Hessian path by
897 /// [`Self::build_with_backend`]; ignored otherwise.
898 pub fn set_kkt_schur(&mut self, schur_indices: Vec<usize>, cfg: pounce_feral::FeralConfig) {
899 self.kkt_schur = Some((schur_indices, cfg));
900 }
901
902 /// Assemble the strategy bundle without a search-direction
903 /// calculator. Used by structural unit tests that don't want to
904 /// pull in a linear-solver backend.
905 pub fn build(&self) -> AlgorithmBundle {
906 self.build_inner(None)
907 }
908
909 /// Same as [`Self::build`] but also constructs the
910 /// `SymLinearSolver → AugSystemSolver → PdFullSpaceSolver →
911 /// PdSearchDirCalc` chain via the supplied `factory`.
912 pub fn build_with_backend(&self, mut factory: LinearBackendFactory) -> AlgorithmBundle {
913 let backend = factory(self.linear_solver);
914 let scaling: Option<Box<dyn pounce_linsol::TSymScalingMethod>> =
915 match self.linear_system_scaling {
916 LinearSystemScalingChoice::None => None,
917 LinearSystemScalingChoice::Ruiz => {
918 Some(Box::new(pounce_linsol::RuizTSymScalingMethod::new()))
919 }
920 LinearSystemScalingChoice::Mc19 => {
921 tracing::warn!(target: "pounce::algorithm",
922 "pounce: linear_system_scaling=mc19 not yet implemented; using no scaling"
923 );
924 None
925 }
926 };
927 let linsol = TSymLinearSolver::new(backend, scaling, self.linear_scaling_on_demand);
928 let inner_aug = StdAugSystemSolver::new(linsol);
929 // Limited-memory mode publishes the Hessian as a
930 // `LowRankUpdateSymMatrix`; wrap the standard solver in the
931 // Sherman-Morrison-Woodbury low-rank solver so the augmented
932 // system factorizes only the diagonal `B0` and the quasi-Newton
933 // update is applied as a rank-`m` correction (`O(n·m)` memory).
934 let is_lbfgs = matches!(
935 self.hessian_approximation,
936 HessianApproxChoice::LimitedMemory
937 );
938 let aug_solver: Box<dyn AugSystemSolver> = if is_lbfgs {
939 Box::new(LowRankAugSystemSolver::new(Box::new(inner_aug)))
940 } else if let Some((indices, cfg)) = self.kkt_schur.clone() {
941 // Block-triangular / Schur KKT path (pounce#180 item 2). Only on the
942 // exact-Hessian feral path — the Schur backend is feral-specific,
943 // and the L-BFGS low-rank Woodbury wrapper owns the (2,2) block.
944 // The Schur solver falls back to `StdAugSystemSolver` transparently
945 // when the partition is unsuitable, so a stray hook never breaks a
946 // solve; we gate on `linear_solver == Feral` here to avoid silently
947 // ignoring a user's explicit MA57 selection.
948 if matches!(self.linear_solver, LinearSolverChoice::Feral) {
949 Box::new(crate::kkt::SchurAugSystemSolver::new(
950 inner_aug, indices, cfg,
951 ))
952 } else {
953 Box::new(inner_aug)
954 }
955 } else {
956 Box::new(inner_aug)
957 };
958 // Inertia-correction / Jacobian-regularization constants (#191):
959 // registered but previously never read. Defaults equal the
960 // registered defaults. `perturb_always_cd` goes through the setter
961 // because it also rebuilds the initial jac-degeneracy state.
962 let mut ph = PdPerturbationHandler::new();
963 ph.delta_xs_max = self.perturbation.max_hessian_perturbation;
964 ph.delta_xs_min = self.perturbation.min_hessian_perturbation;
965 ph.delta_xs_first_inc_fact = self.perturbation.perturb_inc_fact_first;
966 ph.delta_xs_inc_fact = self.perturbation.perturb_inc_fact;
967 ph.delta_xs_dec_fact = self.perturbation.perturb_dec_fact;
968 ph.delta_xs_init = self.perturbation.first_hessian_perturbation;
969 ph.delta_cd_val = self.perturbation.jacobian_regularization_value;
970 ph.delta_cd_exp = self.perturbation.jacobian_regularization_exponent;
971 ph.set_perturb_always_cd(self.perturbation.perturb_always_cd);
972 let perturb = Rc::new(RefCell::new(ph));
973 let mut pd_solver = PdFullSpaceSolver::new(aug_solver, perturb);
974 // Iterative-refinement constants (#191): registered but previously
975 // never read, so overrides were silently dropped. Defaults equal
976 // the registered defaults.
977 pd_solver.min_refinement_steps = self.refinement.min_refinement_steps;
978 pd_solver.max_refinement_steps = self.refinement.max_refinement_steps;
979 pd_solver.residual_ratio_max = self.refinement.residual_ratio_max;
980 pd_solver.residual_ratio_singular = self.refinement.residual_ratio_singular;
981 pd_solver.residual_improvement_factor = self.refinement.residual_improvement_factor;
982 let mut search_dir = PdSearchDirCalc::new(pd_solver);
983 search_dir.mehrotra_algorithm = self.mehrotra_algorithm;
984 search_dir.fast_step_computation = self.fast_step_computation;
985 self.build_inner(Some(search_dir))
986 }
987
988 /// Phase 5b assembly path for the SQP algorithm. Consults
989 /// `self.algorithm`: when `ActiveSetSqp`, constructs an
990 /// `SqpAlgorithm` using the supplied backend factory for the
991 /// QP subproblem solver; otherwise returns `None` so the
992 /// caller can fall back to the IPM `build_with_backend`.
993 ///
994 /// Sister to `build_with_backend`: the SQP algorithm doesn't
995 /// share `AlgorithmBundle`'s shape (no mu_update / no IPM
996 /// line search), so the two paths return different types.
997 pub fn build_sqp_with_backend(
998 &self,
999 mut factory: LinearBackendFactory,
1000 ) -> Option<crate::sqp::SqpAlgorithm> {
1001 if !matches!(self.algorithm, AlgorithmChoice::ActiveSetSqp) {
1002 return None;
1003 }
1004 let backend = factory(self.linear_solver);
1005 let qp_solver = pounce_qp::ParametricActiveSetSolver::new(backend);
1006 Some(
1007 crate::sqp::SqpAlgorithm::new(qp_solver, self.sqp.clone())
1008 .with_qp_options(self.sqp_qp.clone()),
1009 )
1010 }
1011
1012 fn build_inner(&self, search_dir: Option<PdSearchDirCalc>) -> AlgorithmBundle {
1013 let mu_update: Box<dyn crate::mu::r#trait::MuUpdate> = match self.mu_strategy {
1014 MuStrategyChoice::Monotone => {
1015 let mut m = MonotoneMuUpdate::new();
1016 m.mu_init = self.mu.mu_init;
1017 // `mu_max` sentinel `-1` keeps the monotone default
1018 // (1e5); only override on a user-supplied positive.
1019 if self.mu.mu_max > 0.0 {
1020 m.mu_max = self.mu.mu_max;
1021 }
1022 m.mu_min = self.mu.mu_min;
1023 m.mu_target = self.mu.mu_target;
1024 m.mu_linear_decrease_factor = self.mu.mu_linear_decrease_factor;
1025 m.mu_superlinear_decrease_power = self.mu.mu_superlinear_decrease_power;
1026 m.mu_allow_fast_monotone_decrease = self.mu.mu_allow_fast_monotone_decrease;
1027 m.barrier_tol_factor = self.mu.barrier_tol_factor;
1028 m.compl_inf_tol = self.conv_check.compl_inf_tol;
1029 Box::new(m)
1030 }
1031 MuStrategyChoice::Adaptive => {
1032 let mut adaptive = AdaptiveMuUpdate::new();
1033 adaptive.mu_oracle = self.mu_oracle;
1034 adaptive.mu_init = self.mu.mu_init;
1035 // Adaptive treats `mu_max == -1` as "lazy init from
1036 // `mu_max_fact * curr_avrg_compl`" — forward the
1037 // sentinel as-is.
1038 adaptive.mu_max = self.mu.mu_max;
1039 adaptive.mu_max_fact = self.mu.mu_max_fact;
1040 adaptive.mu_min = self.mu.mu_min;
1041 adaptive.compl_inf_tol = self.conv_check.compl_inf_tol;
1042 adaptive.mu_linear_decrease_factor = self.mu.mu_linear_decrease_factor;
1043 adaptive.mu_superlinear_decrease_power = self.mu.mu_superlinear_decrease_power;
1044 adaptive.barrier_tol_factor = self.mu.barrier_tol_factor;
1045 adaptive.sigma_min = self.mu.sigma_min;
1046 adaptive.sigma_max = self.mu.sigma_max;
1047 adaptive.adaptive_mu_globalization = self.mu.adaptive_mu_globalization;
1048 adaptive.qf_norm_type = self.mu.quality_function_norm_type;
1049 adaptive.qf_centrality_type = self.mu.quality_function_centrality;
1050 adaptive.qf_balancing_term = self.mu.quality_function_balancing_term;
1051 adaptive.qf_max_section_steps = self.mu.quality_function_max_section_steps;
1052 adaptive.qf_section_sigma_tol = self.mu.quality_function_section_sigma_tol;
1053 adaptive.qf_section_qf_tol = self.mu.quality_function_section_qf_tol;
1054 adaptive.probing_iterate_quality_factor = self.mu.probing_iterate_quality_factor;
1055 adaptive.adaptive_mu_safeguard_factor = self.mu.adaptive_mu_safeguard_factor;
1056 adaptive.adaptive_mu_monotone_init_factor =
1057 self.mu.adaptive_mu_monotone_init_factor;
1058 adaptive.restore_accepted_iterate = self.mu.adaptive_mu_restore_previous_iterate;
1059 adaptive.adaptive_mu_kkterror_red_iters = self.mu.adaptive_mu_kkterror_red_iters;
1060 adaptive.adaptive_mu_kkterror_red_fact = self.mu.adaptive_mu_kkterror_red_fact;
1061 adaptive.adaptive_mu_kkt_norm = self.mu.adaptive_mu_kkt_norm_type;
1062 Box::new(adaptive)
1063 }
1064 };
1065
1066 let acceptor: Box<dyn BacktrackingLsAcceptor> = match self.line_search_method {
1067 LineSearchChoice::Filter => {
1068 // Filter switching / Armijo / margin constants (#191):
1069 // registered but previously never read. Set them on the
1070 // concrete acceptor before boxing; defaults equal the
1071 // registered defaults, so a run that doesn't set them is
1072 // unchanged.
1073 let mut f = FilterLsAcceptor::default();
1074 f.eta_phi = self.line_search.eta_phi;
1075 f.theta_min_fact = self.line_search.theta_min_fact;
1076 f.theta_max_fact = self.line_search.theta_max_fact;
1077 f.theta_max_row_scale_kappa = self.line_search.theta_max_row_scale_kappa;
1078 f.theta_max_adaptive_trigger = self.line_search.theta_max_adaptive_trigger;
1079 f.theta_max_adaptive_factor = self.line_search.theta_max_adaptive_factor;
1080 f.theta_max_adaptive_max_raises = self.line_search.theta_max_adaptive_max_raises;
1081 f.gamma_phi = self.line_search.gamma_phi;
1082 f.gamma_theta = self.line_search.gamma_theta;
1083 f.s_phi = self.line_search.s_phi;
1084 f.s_theta = self.line_search.s_theta;
1085 f.alpha_min_frac = self.line_search.alpha_min_frac;
1086 f.obj_max_inc = self.line_search.obj_max_inc;
1087 f.max_filter_resets = self.line_search.max_filter_resets;
1088 f.filter_reset_trigger = self.line_search.filter_reset_trigger;
1089 Box::new(f)
1090 }
1091 LineSearchChoice::Penalty => Box::new(PenaltyLsAcceptor::default()),
1092 // CG-penalty acceptor lands with the rest of the
1093 // CG-penalty path; fall back to the penalty acceptor's
1094 // surface for now.
1095 LineSearchChoice::CgPenalty => Box::new(PenaltyLsAcceptor::default()),
1096 };
1097 let mut line_search = BacktrackingLineSearch::new(acceptor);
1098 line_search.watchdog_shortened_iter_trigger =
1099 self.line_search.watchdog_shortened_iter_trigger;
1100 line_search.watchdog_trial_iter_max = self.line_search.watchdog_trial_iter_max;
1101 line_search.soft_resto_pderror_reduction_factor =
1102 self.line_search.soft_resto_pderror_reduction_factor;
1103 line_search.max_soft_resto_iters = self.line_search.max_soft_resto_iters;
1104 line_search.accept_every_trial_step = self.line_search.accept_every_trial_step;
1105 line_search.alpha_for_y = self.line_search.alpha_for_y;
1106 // Second-order-correction constants (#191): registered but
1107 // previously never read. Same direct-field pattern as the
1108 // watchdog knobs above.
1109 line_search.max_soc = self.line_search.max_soc;
1110 line_search.kappa_soc = self.line_search.kappa_soc;
1111 line_search.soc_method = self.line_search.soc_method;
1112
1113 let conv_check: Box<dyn crate::conv_check::r#trait::ConvCheck> =
1114 Box::new(OptErrorConvCheck {
1115 tol: self.conv_check.tol,
1116 dual_inf_tol: self.conv_check.dual_inf_tol,
1117 constr_viol_tol: self.conv_check.constr_viol_tol,
1118 compl_inf_tol: self.conv_check.compl_inf_tol,
1119 acceptable_tol: self.conv_check.acceptable_tol,
1120 acceptable_dual_inf_tol: self.conv_check.acceptable_dual_inf_tol,
1121 acceptable_constr_viol_tol: self.conv_check.acceptable_constr_viol_tol,
1122 acceptable_compl_inf_tol: self.conv_check.acceptable_compl_inf_tol,
1123 acceptable_obj_change_tol: self.conv_check.acceptable_obj_change_tol,
1124 acceptable_iter: self.conv_check.acceptable_iter,
1125 max_iter: self.conv_check.max_iter,
1126 max_cpu_time: self.conv_check.max_cpu_time,
1127 max_wall_time: self.conv_check.max_wall_time,
1128 acceptable_count: 0,
1129 last_acceptable_obj: None,
1130 infeas_stationarity_tol: self.conv_check.infeas_stationarity_tol,
1131 infeas_viol_kappa: self.conv_check.infeas_viol_kappa,
1132 infeas_max_streak: self.conv_check.infeas_max_streak,
1133 infeas_streak: 0,
1134 obj_scale_certificate_threshold: self.conv_check.obj_scale_certificate_threshold,
1135 primal_noise_floor_kappa: self.conv_check.primal_noise_floor_kappa,
1136 acceptable_progress_kappa: self.conv_check.acceptable_progress_kappa,
1137 acceptable_window: std::collections::VecDeque::new(),
1138 acceptable_progress_refusals: 0,
1139 dual_inf_scale_kappa: self.conv_check.dual_inf_scale_kappa,
1140 dual_floor_reported: false,
1141 veto_fired: false,
1142 acceptable_veto_fired: false,
1143 masked_acceptable_veto_fired: false,
1144 veto_extra_iters: 0,
1145 rel_infeas_extra_iters: 0,
1146 prev_rel_viol: f64::NAN,
1147 });
1148
1149 let init: Box<dyn crate::init::r#trait::IterateInitializer> = if self.warm_start_init_point
1150 {
1151 Box::new(WarmStartIterateInitializer::with_options(
1152 resolved_warm_options(&self.warm, &self.init),
1153 ))
1154 } else {
1155 let mut d = DefaultIterateInitializer::with_eq_mult_calculator(Box::new(
1156 LeastSquareMults::new(),
1157 ));
1158 d.bound_push = self.init.bound_push;
1159 d.bound_frac = self.init.bound_frac;
1160 d.slack_bound_push = self.init.slack_bound_push;
1161 d.slack_bound_frac = self.init.slack_bound_frac;
1162 d.constr_mult_init_max = self.init.constr_mult_init_max;
1163 d.bound_mult_init_val = self.init.bound_mult_init_val;
1164 d.bound_mult_init_method = self.init.bound_mult_init_method.clone();
1165 d.least_square_init_primal = self.init.least_square_init_primal;
1166 Box::new(d)
1167 };
1168
1169 let eq_mult: Box<dyn crate::eq_mult::r#trait::EqMultCalculator> =
1170 Box::new(LeastSquareMults::new());
1171
1172 let hess: Box<dyn crate::hess::r#trait::HessianUpdater> = match self.hessian_approximation {
1173 HessianApproxChoice::Exact => Box::new(ExactHessianUpdater::new()),
1174 HessianApproxChoice::LimitedMemory => Box::new(LimMemQuasiNewtonUpdater {
1175 update_type: self.limited_memory_update_type,
1176 max_history: self.limited_memory_max_history,
1177 init_val_max: self.limited_memory_init_val_max,
1178 init_val_min: self.limited_memory_init_val_min,
1179 ..LimMemQuasiNewtonUpdater::default()
1180 }),
1181 };
1182
1183 let iter_output: Box<dyn crate::output::r#trait::IterationOutput> = {
1184 use crate::output::orig::{InfPrTag, PrintInfoString};
1185 let mut o = OrigIterationOutput::new();
1186 o.print_frequency_iter = self.output.print_frequency_iter;
1187 o.print_frequency_time = self.output.print_frequency_time;
1188 o.print_info_string = if self.output.print_info_string {
1189 PrintInfoString::Yes
1190 } else {
1191 PrintInfoString::No
1192 };
1193 o.inf_pr_output = if self.output.inf_pr_output_internal {
1194 InfPrTag::Internal
1195 } else {
1196 InfPrTag::Original
1197 };
1198 Box::new(o)
1199 };
1200
1201 AlgorithmBundle {
1202 mu_update,
1203 conv_check,
1204 init,
1205 eq_mult,
1206 hess,
1207 line_search,
1208 iter_output,
1209 search_dir,
1210 }
1211 }
1212}
1213
1214#[cfg(test)]
1215mod tests {
1216 use super::*;
1217
1218 #[test]
1219 fn warm_options_take_the_init_default_not_their_own() {
1220 let mut init = InitOptions::default();
1221 init.bound_mult_init_val = 10.0; // the Mehrotra override value
1222 let mut warm = WarmStartOptions::default();
1223 warm.bound_mult_init_val = 123.0; // stale copy must lose
1224 let resolved = resolved_warm_options(&warm, &init);
1225 assert_eq!(resolved.bound_mult_init_val, 10.0);
1226 // everything else passes through untouched
1227 assert_eq!(resolved.mult_bound_push, warm.mult_bound_push);
1228 assert_eq!(resolved.target_mu, warm.target_mu);
1229 }
1230
1231 #[test]
1232 fn default_builder_assembles() {
1233 let bundle = AlgorithmBuilder::new().build();
1234 // Sanity: the placeholder traits compile and the boxed
1235 // strategies don't panic on construction.
1236 let _ = bundle.line_search.acceptor();
1237 assert!(bundle.search_dir.is_none());
1238 }
1239
1240 #[test]
1241 fn build_with_backend_assembles_search_dir_chain() {
1242 // Drive the builder with the FERAL backend factory; the
1243 // resulting bundle should expose a populated `PdSearchDirCalc`.
1244 let factory: LinearBackendFactory = Box::new(|_| {
1245 Box::new(pounce_feral::FeralSolverInterface::new())
1246 as Box<dyn SparseSymLinearSolverInterface>
1247 });
1248 let bundle = AlgorithmBuilder::new().build_with_backend(factory);
1249 assert!(bundle.search_dir.is_some());
1250 }
1251
1252 #[test]
1253 fn limited_memory_sr1_propagates() {
1254 let b = AlgorithmBuilder {
1255 hessian_approximation: HessianApproxChoice::LimitedMemory,
1256 limited_memory_update_type: UpdateType::Sr1,
1257 ..AlgorithmBuilder::default()
1258 };
1259 let _bundle = b.build();
1260 }
1261
1262 #[test]
1263 fn every_strategy_combination_assembles_without_panic() {
1264 let solvers = [LinearSolverChoice::Ma57, LinearSolverChoice::Feral];
1265 let mu = [MuStrategyChoice::Monotone, MuStrategyChoice::Adaptive];
1266 let hess = [
1267 HessianApproxChoice::Exact,
1268 HessianApproxChoice::LimitedMemory,
1269 ];
1270 let ls = [
1271 LineSearchChoice::Filter,
1272 LineSearchChoice::CgPenalty,
1273 LineSearchChoice::Penalty,
1274 ];
1275 for &linear_solver in &solvers {
1276 for &mu_strategy in &mu {
1277 for &hessian_approximation in &hess {
1278 for &line_search_method in &ls {
1279 let _ = AlgorithmBuilder {
1280 algorithm: AlgorithmChoice::default(),
1281 linear_solver,
1282 linear_system_scaling: LinearSystemScalingChoice::None,
1283 linear_scaling_on_demand: true,
1284 mu_strategy,
1285 mu_oracle: MuOracleKind::QualityFunction,
1286 hessian_approximation,
1287 limited_memory_update_type: UpdateType::Bfgs,
1288 limited_memory_max_history: 6,
1289 limited_memory_init_val_max: 1e8,
1290 limited_memory_init_val_min: 1e-8,
1291 line_search_method,
1292 warm_start_init_point: false,
1293 mehrotra_algorithm: false,
1294 fast_step_computation: false,
1295 kappa_sigma: 1e10,
1296 kappa_d: 1e-5,
1297 tiny_step_tol: 10.0 * Number::EPSILON,
1298 tiny_step_y_tol: 1e-2,
1299 diverging_iterates_tol: 1e20,
1300 dual_diverging_streak: 0,
1301 resto_decline_deferrals: 1,
1302 resto_decline_progress_ratio: 0.5,
1303 kkt_fidelity_tol: 0.0,
1304 conv_check: ConvCheckOptions::default(),
1305 mu: MuOptions::default(),
1306 line_search: LineSearchOptions::default(),
1307 refinement: RefinementOptions::default(),
1308 perturbation: PerturbationOptions::default(),
1309 resto: RestoOptions::default(),
1310 output: OutputOptions::default(),
1311 warm: WarmStartOptions::default(),
1312 sqp: crate::sqp::SqpOptions::default(),
1313 sqp_qp: pounce_qp::QpOptions::default(),
1314 init: InitOptions::default(),
1315 kkt_schur: None,
1316 }
1317 .build();
1318 }
1319 }
1320 }
1321 }
1322 }
1323}