Skip to main content

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}
149
150impl Default for ConvCheckOptions {
151    fn default() -> Self {
152        Self {
153            tol: 1e-8,
154            dual_inf_tol: 1.0,
155            constr_viol_tol: 1e-4,
156            compl_inf_tol: 1e-4,
157            acceptable_tol: 1e-6,
158            acceptable_dual_inf_tol: 1e10,
159            acceptable_constr_viol_tol: 1e-2,
160            acceptable_compl_inf_tol: 1e-2,
161            acceptable_obj_change_tol: 1e20,
162            acceptable_iter: 15,
163            max_iter: 3000,
164            max_cpu_time: 1e6,
165            max_wall_time: 1e6,
166            infeas_stationarity_tol: 1e-8,
167            infeas_viol_kappa: 1e2,
168            infeas_max_streak: 5,
169        }
170    }
171}
172
173#[derive(Debug, Clone)]
174pub struct AlgorithmBuilder {
175    /// Top-level algorithm dispatch. Default `InteriorPoint` ⇒
176    /// `build_with_backend` returns the existing `AlgorithmBundle`
177    /// (consumed by `IpoptAlgorithm`). `ActiveSetSqp` ⇒ caller
178    /// must use `build_sqp_with_backend` to assemble the Phase 5b
179    /// `SqpAlgorithm`. The two builder methods sit side by side
180    /// because the assembled algorithm shape differs (IPM bundle
181    /// vs SQP struct).
182    pub algorithm: AlgorithmChoice,
183    pub linear_solver: LinearSolverChoice,
184    /// Symmetric scaling method for the augmented KKT system. Wired
185    /// into [`TSymLinearSolver`] by [`Self::build_with_backend`].
186    /// Mirrors upstream `linear_system_scaling` (`IpAlgBuilder.cpp:538-560`).
187    pub linear_system_scaling: LinearSystemScalingChoice,
188    /// Lazy-vs-eager scaling toggle (`linear_scaling_on_demand`,
189    /// `IpTSymLinearSolver.cpp:50-58`). Only consulted when
190    /// `linear_system_scaling != None`. Upstream default is `true`
191    /// (compute scaling only on the first solve that fails / shows
192    /// poor conditioning); pounce mirrors that. Set to `false` to
193    /// scale every factorization.
194    pub linear_scaling_on_demand: bool,
195    pub mu_strategy: MuStrategyChoice,
196    /// Selector forwarded to [`AdaptiveMuUpdate`] when
197    /// `mu_strategy = Adaptive`. Ignored for `Monotone`. Defaults to
198    /// `QualityFunction` per upstream's `RegisterOptions` default.
199    pub mu_oracle: MuOracleKind,
200    pub hessian_approximation: HessianApproxChoice,
201    pub limited_memory_update_type: UpdateType,
202    pub line_search_method: LineSearchChoice,
203    pub warm_start_init_point: bool,
204    /// `mehrotra_algorithm` — when true, [`PdSearchDirCalc`] folds
205    /// the Mehrotra second-order complementarity term into the
206    /// search-direction RHS. Mirrors upstream's
207    /// `IpAlgBuilder.cpp:Mehrotra` flag. Requires `mu_strategy =
208    /// Adaptive` so that an affine step is computed each iteration;
209    /// [`Self::build_with_backend`] does not enforce this — the
210    /// option-parser in `application.rs` is responsible for the
211    /// cascading defaults (`mu_oracle = probing` etc.).
212    pub mehrotra_algorithm: bool,
213    pub conv_check: ConvCheckOptions,
214    pub mu: MuOptions,
215    pub line_search: LineSearchOptions,
216    pub output: OutputOptions,
217    pub warm: WarmStartOptions,
218    /// SQP-specific options (consulted only when
219    /// `algorithm = ActiveSetSqp`).
220    pub sqp: crate::sqp::SqpOptions,
221    pub init: InitOptions,
222}
223
224/// Knobs read off `OptionsList` and baked into
225/// [`DefaultIterateInitializer`]. Defaults mirror
226/// `IpDefaultIterateInitializer.cpp:RegisterOptions`. The Mehrotra
227/// cascade in `application.rs` overrides `bound_push`, `bound_frac`,
228/// and `bound_mult_init_val` to upstream's more-aggressive values
229/// (`10`, `0.2`, `1.0`).
230#[derive(Debug, Clone)]
231pub struct InitOptions {
232    pub bound_push: Number,
233    pub bound_frac: Number,
234    pub slack_bound_push: Number,
235    pub slack_bound_frac: Number,
236    pub constr_mult_init_max: Number,
237    pub bound_mult_init_val: Number,
238    /// `bound_mult_init_method`: `"constant"` (default) or `"mu-based"`
239    /// (matches upstream's `IpDefaultIterateInitializer.cpp`).
240    pub bound_mult_init_method: String,
241    /// `least_square_init_primal` — replace the user's starting `x`
242    /// with the min-norm primal that satisfies the linearized
243    /// constraints. Used by the Mehrotra cascade in `application.rs`
244    /// to drop iter-0 primal infeasibility on LP-shaped problems.
245    /// Mirrors upstream `IpDefaultIterateInitializer.cpp:200-222`.
246    pub least_square_init_primal: bool,
247}
248
249impl Default for InitOptions {
250    fn default() -> Self {
251        Self {
252            bound_push: 1e-2,
253            bound_frac: 1e-2,
254            slack_bound_push: 1e-2,
255            slack_bound_frac: 1e-2,
256            constr_mult_init_max: 1e3,
257            bound_mult_init_val: 1.0,
258            bound_mult_init_method: "constant".into(),
259            least_square_init_primal: false,
260        }
261    }
262}
263
264/// Knobs read off `OptionsList` and baked into
265/// [`WarmStartIterateInitializer`]. Defaults mirror
266/// `IpWarmStartIterateInitializer.cpp:RegisterOptions`.
267///
268/// Wired today: `mult_init_max` (clamps |y_c|, |y_d| and caps z/v
269/// blocks) and `target_mu` (overrides `data.curr_mu` at iter 0).
270/// The remaining knobs (`bound_push`, `bound_frac`, `slack_bound_push`,
271/// `slack_bound_frac`, `mult_bound_push`, `entire_iterate`,
272/// `same_structure`) are stored on the initializer but not yet
273/// consumed — `WarmStartIterateInitializer::set_initial_iterates`
274/// currently trusts the caller-populated `data.curr` rather than
275/// re-running the upstream `push_variables` machinery.
276#[derive(Debug, Clone)]
277pub struct WarmStartOptions {
278    pub bound_push: Number,
279    pub bound_frac: Number,
280    pub slack_bound_push: Number,
281    pub slack_bound_frac: Number,
282    pub mult_bound_push: Number,
283    pub mult_init_max: Number,
284    pub target_mu: Number,
285    pub entire_iterate: bool,
286    pub same_structure: bool,
287}
288
289impl Default for WarmStartOptions {
290    fn default() -> Self {
291        Self {
292            bound_push: 1e-3,
293            bound_frac: 1e-3,
294            slack_bound_push: 1e-3,
295            slack_bound_frac: 1e-3,
296            mult_bound_push: 1e-3,
297            mult_init_max: 1e6,
298            target_mu: 0.0,
299            entire_iterate: false,
300            same_structure: false,
301        }
302    }
303}
304
305/// Knobs read off `OptionsList` and baked into the assembled
306/// `MonotoneMuUpdate` or `AdaptiveMuUpdate`. Defaults mirror
307/// `IpMonotoneMuUpdate.cpp` / `IpAdaptiveMuUpdate.cpp:RegisterOptions`.
308/// `mu_max` defaults to the sentinel `-1`; positive values are baked
309/// into both updaters at build time (adaptive interprets `-1` as
310/// "lazy-init from `mu_max_fact * avrg_compl`").
311#[derive(Debug, Clone)]
312pub struct MuOptions {
313    pub mu_init: Number,
314    pub mu_max: Number,
315    pub mu_max_fact: Number,
316    pub mu_min: Number,
317    pub mu_target: Number,
318    pub mu_linear_decrease_factor: Number,
319    pub mu_superlinear_decrease_power: Number,
320    pub mu_allow_fast_monotone_decrease: bool,
321    pub barrier_tol_factor: Number,
322    /// `sigma_max` / `sigma_min` — clamp on the centering parameter σ
323    /// chosen by `QualityFunctionMuOracle`. Only consumed when
324    /// `mu_strategy=adaptive` and `mu_oracle=quality-function`.
325    /// Defaults from `IpQualityFunctionMuOracle.cpp:RegisterOptions`.
326    pub sigma_max: Number,
327    pub sigma_min: Number,
328    /// `adaptive_mu_globalization` — globalization strategy for the
329    /// adaptive μ-selection mode. Mirrors
330    /// `IpAdaptiveMuUpdate.cpp:RegisterOptions`. Default is
331    /// `ObjConstrFilter`; the Mehrotra cascade switches to
332    /// `NeverMonotoneMode` to disable globalization entirely.
333    pub adaptive_mu_globalization: crate::mu::adaptive::AdaptiveMuGlobalization,
334    /// `quality_function_norm_type` — norm used inside the quality
335    /// function to aggregate the three KKT components. Forwarded to
336    /// `QualityFunctionMuOracle` when `mu_oracle=quality-function`.
337    pub quality_function_norm_type: crate::mu::oracle::quality_function::NormType,
338    /// `quality_function_centrality` — centrality penalty term added
339    /// to the quality function.
340    pub quality_function_centrality: crate::mu::oracle::quality_function::CentralityType,
341    /// `quality_function_balancing_term` — balancing penalty term in
342    /// the quality function (kicks in when complementarity is far
343    /// below infeasibilities).
344    pub quality_function_balancing_term: crate::mu::oracle::quality_function::BalancingTermType,
345    /// `quality_function_max_section_steps` — cap on golden-section
346    /// iterations when picking σ. Default 8.
347    pub quality_function_max_section_steps: i32,
348    /// `quality_function_section_sigma_tol` — width tolerance in
349    /// σ-space for golden section. Default 1e-2.
350    pub quality_function_section_sigma_tol: Number,
351    /// `quality_function_section_qf_tol` — relative flatness
352    /// tolerance for golden section. Default 0.0.
353    pub quality_function_section_qf_tol: Number,
354    /// `adaptive_mu_safeguard_factor` — guard for the LOQO fallback
355    /// in adaptive mode. Default 0.0.
356    pub adaptive_mu_safeguard_factor: Number,
357    /// `adaptive_mu_monotone_init_factor` — multiplier on the
358    /// average complementarity when seeding monotone mode after a
359    /// free-mode bailout. Default 0.8.
360    pub adaptive_mu_monotone_init_factor: Number,
361    /// `adaptive_mu_restore_previous_iterate` — restore the most
362    /// recent free-mode iterate when switching to fixed mode.
363    /// Default `false`.
364    pub adaptive_mu_restore_previous_iterate: bool,
365    /// `adaptive_mu_kkterror_red_iters` — window length for the
366    /// `KKT_ERROR` globalization history. Default 4.
367    pub adaptive_mu_kkterror_red_iters: usize,
368    /// `adaptive_mu_kkterror_red_fact` — required relative reduction
369    /// of the KKT error over the window. Default 0.9999.
370    pub adaptive_mu_kkterror_red_fact: Number,
371    /// `adaptive_mu_kkt_norm_type` — norm used to score the iterate
372    /// in adaptive globalization decisions.
373    pub adaptive_mu_kkt_norm_type: crate::mu::adaptive::AdaptiveMuKktNorm,
374    /// `probing_iterate_quality_factor` (default 1e4, pounce-specific
375    /// — see pounce#58). When the probing (Mehrotra) μ-oracle is
376    /// about to read `curr_avrg_compl()` for its `mu_curr` input, a
377    /// single imbalanced `(s_i, z_i)` pair can inflate the average
378    /// 5+ orders above the stored `data.curr_mu`. The oracle then
379    /// returns `σ · mu_curr` ≫ previous μ, throwing the iterate out
380    /// of the convergence neighborhood. This guard short-circuits
381    /// that case by signalling restoration when the ratio
382    /// `curr_avrg_compl / curr_mu` exceeds the factor. Set to 0 or
383    /// any non-positive value to disable.
384    pub probing_iterate_quality_factor: Number,
385}
386
387impl Default for MuOptions {
388    fn default() -> Self {
389        Self {
390            mu_init: 0.1,
391            mu_max: -1.0,
392            mu_max_fact: 1e3,
393            mu_min: 1e-11,
394            mu_target: 0.0,
395            mu_linear_decrease_factor: 0.2,
396            mu_superlinear_decrease_power: 1.5,
397            mu_allow_fast_monotone_decrease: true,
398            barrier_tol_factor: 10.0,
399            sigma_max: 1e2,
400            sigma_min: 1e-6,
401            adaptive_mu_globalization:
402                crate::mu::adaptive::AdaptiveMuGlobalization::ObjConstrFilter,
403            quality_function_norm_type:
404                crate::mu::oracle::quality_function::NormType::TwoNormSquared,
405            quality_function_centrality: crate::mu::oracle::quality_function::CentralityType::None,
406            quality_function_balancing_term:
407                crate::mu::oracle::quality_function::BalancingTermType::None,
408            quality_function_max_section_steps: 8,
409            quality_function_section_sigma_tol: 1e-2,
410            quality_function_section_qf_tol: 0.0,
411            adaptive_mu_safeguard_factor: 0.0,
412            adaptive_mu_monotone_init_factor: 0.8,
413            adaptive_mu_restore_previous_iterate: false,
414            adaptive_mu_kkterror_red_iters: 4,
415            adaptive_mu_kkterror_red_fact: 0.9999,
416            adaptive_mu_kkt_norm_type: crate::mu::adaptive::AdaptiveMuKktNorm::TwoNormSquared,
417            probing_iterate_quality_factor: 1e4,
418        }
419    }
420}
421
422/// Knobs baked into the assembled [`BacktrackingLineSearch`]. Defaults
423/// mirror `IpBacktrackingLineSearch.cpp:RegisterOptions`.
424#[derive(Debug, Clone)]
425pub struct LineSearchOptions {
426    pub watchdog_shortened_iter_trigger: Index,
427    pub watchdog_trial_iter_max: Index,
428    /// `soft_resto_pderror_reduction_factor` — required relative
429    /// reduction in the primal-dual error for a soft-resto step.
430    /// `0` disables the soft restoration phase.
431    pub soft_resto_pderror_reduction_factor: Number,
432    /// `max_soft_resto_iters` — cap on consecutive soft-resto
433    /// iterations before full restoration is forced.
434    pub max_soft_resto_iters: Index,
435    /// `accept_every_trial_step` — short-circuits the filter / alpha
436    /// loop and accepts the full fraction-to-the-boundary step every
437    /// outer iteration. Mirrors upstream's
438    /// `IpBacktrackingLineSearch::accept_every_trial_step_`. Drops
439    /// global convergence guarantees; only safe for problems where the
440    /// Newton step is already a descent step (LPs, convex QPs). The
441    /// Mehrotra cascade in `application.rs` flips this on.
442    pub accept_every_trial_step: bool,
443    /// `alpha_for_y` — policy for the equality-multiplier (y_c / y_d)
444    /// step length. Upstream default is `Primal`; the Mehrotra cascade
445    /// switches to `BoundMult`.
446    pub alpha_for_y: crate::line_search::backtracking::AlphaForY,
447}
448
449impl Default for LineSearchOptions {
450    fn default() -> Self {
451        Self {
452            watchdog_shortened_iter_trigger: 10,
453            watchdog_trial_iter_max: 3,
454            soft_resto_pderror_reduction_factor: 1.0 - 1e-4,
455            max_soft_resto_iters: 10,
456            accept_every_trial_step: false,
457            alpha_for_y: crate::line_search::backtracking::AlphaForY::Primal,
458        }
459    }
460}
461
462/// Knobs baked into the assembled [`OrigIterationOutput`]. Defaults
463/// mirror `IpOrigIterationOutput.cpp:RegisterOptions` /
464/// `IpAlgorithmRegOp.cpp`.
465#[derive(Debug, Clone)]
466pub struct OutputOptions {
467    pub print_frequency_iter: Index,
468    pub print_frequency_time: Number,
469    /// `print_info_string` (default `false`). When on, the iter row
470    /// ends with the contents of `IpoptData::info_string` so users
471    /// can read the per-iteration diagnostic tags.
472    pub print_info_string: bool,
473    /// `inf_pr_output` — `"original"` (default) prints the unscaled
474    /// NLP primal infeasibility; `"internal"` prints the internal
475    /// reformulated violation. Only meaningful once NLP-side scaling
476    /// is in play; until then both modes produce the same number.
477    pub inf_pr_output_internal: bool,
478}
479
480impl Default for OutputOptions {
481    fn default() -> Self {
482        Self {
483            print_frequency_iter: 1,
484            print_frequency_time: 0.0,
485            print_info_string: false,
486            inf_pr_output_internal: false,
487        }
488    }
489}
490
491impl Default for AlgorithmBuilder {
492    fn default() -> Self {
493        Self {
494            algorithm: AlgorithmChoice::default(),
495            linear_solver: LinearSolverChoice::Feral,
496            linear_system_scaling: LinearSystemScalingChoice::None,
497            linear_scaling_on_demand: true,
498            mu_strategy: MuStrategyChoice::Monotone,
499            mu_oracle: MuOracleKind::QualityFunction,
500            hessian_approximation: HessianApproxChoice::Exact,
501            limited_memory_update_type: UpdateType::Bfgs,
502            line_search_method: LineSearchChoice::Filter,
503            warm_start_init_point: false,
504            mehrotra_algorithm: false,
505            conv_check: ConvCheckOptions::default(),
506            mu: MuOptions::default(),
507            line_search: LineSearchOptions::default(),
508            output: OutputOptions::default(),
509            warm: WarmStartOptions::default(),
510            sqp: crate::sqp::SqpOptions::default(),
511            init: InitOptions::default(),
512        }
513    }
514}
515
516impl AlgorithmBuilder {
517    pub fn new() -> Self {
518        Self::default()
519    }
520
521    /// Assemble the strategy bundle without a search-direction
522    /// calculator. Used by structural unit tests that don't want to
523    /// pull in a linear-solver backend.
524    pub fn build(&self) -> AlgorithmBundle {
525        self.build_inner(None)
526    }
527
528    /// Same as [`Self::build`] but also constructs the
529    /// `SymLinearSolver → AugSystemSolver → PdFullSpaceSolver →
530    /// PdSearchDirCalc` chain via the supplied `factory`.
531    pub fn build_with_backend(&self, mut factory: LinearBackendFactory) -> AlgorithmBundle {
532        let backend = factory(self.linear_solver);
533        let scaling: Option<Box<dyn pounce_linsol::TSymScalingMethod>> =
534            match self.linear_system_scaling {
535                LinearSystemScalingChoice::None => None,
536                LinearSystemScalingChoice::Ruiz => {
537                    Some(Box::new(pounce_linsol::RuizTSymScalingMethod::new()))
538                }
539                LinearSystemScalingChoice::Mc19 => {
540                    tracing::warn!(target: "pounce::algorithm",
541                        "pounce: linear_system_scaling=mc19 not yet implemented; using no scaling"
542                    );
543                    None
544                }
545            };
546        let linsol = TSymLinearSolver::new(backend, scaling, self.linear_scaling_on_demand);
547        let inner_aug = StdAugSystemSolver::new(linsol);
548        // Limited-memory mode publishes the Hessian as a
549        // `LowRankUpdateSymMatrix`; wrap the standard solver in the
550        // Sherman-Morrison-Woodbury low-rank solver so the augmented
551        // system factorizes only the diagonal `B0` and the quasi-Newton
552        // update is applied as a rank-`m` correction (`O(n·m)` memory).
553        let aug_solver: Box<dyn AugSystemSolver> = if matches!(
554            self.hessian_approximation,
555            HessianApproxChoice::LimitedMemory
556        ) {
557            Box::new(LowRankAugSystemSolver::new(Box::new(inner_aug)))
558        } else {
559            Box::new(inner_aug)
560        };
561        let perturb = Rc::new(RefCell::new(PdPerturbationHandler::new()));
562        let pd_solver = PdFullSpaceSolver::new(aug_solver, perturb);
563        let mut search_dir = PdSearchDirCalc::new(pd_solver);
564        search_dir.mehrotra_algorithm = self.mehrotra_algorithm;
565        self.build_inner(Some(search_dir))
566    }
567
568    /// Phase 5b assembly path for the SQP algorithm. Consults
569    /// `self.algorithm`: when `ActiveSetSqp`, constructs an
570    /// `SqpAlgorithm` using the supplied backend factory for the
571    /// QP subproblem solver; otherwise returns `None` so the
572    /// caller can fall back to the IPM `build_with_backend`.
573    ///
574    /// Sister to `build_with_backend`: the SQP algorithm doesn't
575    /// share `AlgorithmBundle`'s shape (no mu_update / no IPM
576    /// line search), so the two paths return different types.
577    pub fn build_sqp_with_backend(
578        &self,
579        mut factory: LinearBackendFactory,
580    ) -> Option<crate::sqp::SqpAlgorithm> {
581        if !matches!(self.algorithm, AlgorithmChoice::ActiveSetSqp) {
582            return None;
583        }
584        let backend = factory(self.linear_solver);
585        let qp_solver = pounce_qp::ParametricActiveSetSolver::new(backend);
586        Some(crate::sqp::SqpAlgorithm::new(qp_solver, self.sqp.clone()))
587    }
588
589    fn build_inner(&self, search_dir: Option<PdSearchDirCalc>) -> AlgorithmBundle {
590        let mu_update: Box<dyn crate::mu::r#trait::MuUpdate> = match self.mu_strategy {
591            MuStrategyChoice::Monotone => {
592                let mut m = MonotoneMuUpdate::new();
593                m.mu_init = self.mu.mu_init;
594                // `mu_max` sentinel `-1` keeps the monotone default
595                // (1e5); only override on a user-supplied positive.
596                if self.mu.mu_max > 0.0 {
597                    m.mu_max = self.mu.mu_max;
598                }
599                m.mu_min = self.mu.mu_min;
600                m.mu_target = self.mu.mu_target;
601                m.mu_linear_decrease_factor = self.mu.mu_linear_decrease_factor;
602                m.mu_superlinear_decrease_power = self.mu.mu_superlinear_decrease_power;
603                m.mu_allow_fast_monotone_decrease = self.mu.mu_allow_fast_monotone_decrease;
604                m.barrier_tol_factor = self.mu.barrier_tol_factor;
605                m.compl_inf_tol = self.conv_check.compl_inf_tol;
606                Box::new(m)
607            }
608            MuStrategyChoice::Adaptive => {
609                let mut adaptive = AdaptiveMuUpdate::new();
610                adaptive.mu_oracle = self.mu_oracle;
611                adaptive.mu_init = self.mu.mu_init;
612                // Adaptive treats `mu_max == -1` as "lazy init from
613                // `mu_max_fact * curr_avrg_compl`" — forward the
614                // sentinel as-is.
615                adaptive.mu_max = self.mu.mu_max;
616                adaptive.mu_max_fact = self.mu.mu_max_fact;
617                adaptive.mu_min = self.mu.mu_min;
618                adaptive.mu_linear_decrease_factor = self.mu.mu_linear_decrease_factor;
619                adaptive.mu_superlinear_decrease_power = self.mu.mu_superlinear_decrease_power;
620                adaptive.barrier_tol_factor = self.mu.barrier_tol_factor;
621                adaptive.sigma_min = self.mu.sigma_min;
622                adaptive.sigma_max = self.mu.sigma_max;
623                adaptive.adaptive_mu_globalization = self.mu.adaptive_mu_globalization;
624                adaptive.qf_norm_type = self.mu.quality_function_norm_type;
625                adaptive.qf_centrality_type = self.mu.quality_function_centrality;
626                adaptive.qf_balancing_term = self.mu.quality_function_balancing_term;
627                adaptive.qf_max_section_steps = self.mu.quality_function_max_section_steps;
628                adaptive.qf_section_sigma_tol = self.mu.quality_function_section_sigma_tol;
629                adaptive.qf_section_qf_tol = self.mu.quality_function_section_qf_tol;
630                adaptive.probing_iterate_quality_factor = self.mu.probing_iterate_quality_factor;
631                adaptive.adaptive_mu_safeguard_factor = self.mu.adaptive_mu_safeguard_factor;
632                adaptive.adaptive_mu_monotone_init_factor =
633                    self.mu.adaptive_mu_monotone_init_factor;
634                adaptive.restore_accepted_iterate = self.mu.adaptive_mu_restore_previous_iterate;
635                adaptive.adaptive_mu_kkterror_red_iters = self.mu.adaptive_mu_kkterror_red_iters;
636                adaptive.adaptive_mu_kkterror_red_fact = self.mu.adaptive_mu_kkterror_red_fact;
637                adaptive.adaptive_mu_kkt_norm = self.mu.adaptive_mu_kkt_norm_type;
638                Box::new(adaptive)
639            }
640        };
641
642        let acceptor: Box<dyn BacktrackingLsAcceptor> = match self.line_search_method {
643            LineSearchChoice::Filter => Box::new(FilterLsAcceptor::default()),
644            LineSearchChoice::Penalty => Box::new(PenaltyLsAcceptor::default()),
645            // CG-penalty acceptor lands with the rest of the
646            // CG-penalty path; fall back to the penalty acceptor's
647            // surface for now.
648            LineSearchChoice::CgPenalty => Box::new(PenaltyLsAcceptor::default()),
649        };
650        let mut line_search = BacktrackingLineSearch::new(acceptor);
651        line_search.watchdog_shortened_iter_trigger =
652            self.line_search.watchdog_shortened_iter_trigger;
653        line_search.watchdog_trial_iter_max = self.line_search.watchdog_trial_iter_max;
654        line_search.soft_resto_pderror_reduction_factor =
655            self.line_search.soft_resto_pderror_reduction_factor;
656        line_search.max_soft_resto_iters = self.line_search.max_soft_resto_iters;
657        line_search.accept_every_trial_step = self.line_search.accept_every_trial_step;
658        line_search.alpha_for_y = self.line_search.alpha_for_y;
659
660        let conv_check: Box<dyn crate::conv_check::r#trait::ConvCheck> =
661            Box::new(OptErrorConvCheck {
662                tol: self.conv_check.tol,
663                dual_inf_tol: self.conv_check.dual_inf_tol,
664                constr_viol_tol: self.conv_check.constr_viol_tol,
665                compl_inf_tol: self.conv_check.compl_inf_tol,
666                acceptable_tol: self.conv_check.acceptable_tol,
667                acceptable_dual_inf_tol: self.conv_check.acceptable_dual_inf_tol,
668                acceptable_constr_viol_tol: self.conv_check.acceptable_constr_viol_tol,
669                acceptable_compl_inf_tol: self.conv_check.acceptable_compl_inf_tol,
670                acceptable_obj_change_tol: self.conv_check.acceptable_obj_change_tol,
671                acceptable_iter: self.conv_check.acceptable_iter,
672                max_iter: self.conv_check.max_iter,
673                max_cpu_time: self.conv_check.max_cpu_time,
674                max_wall_time: self.conv_check.max_wall_time,
675                acceptable_count: 0,
676                last_acceptable_obj: None,
677                infeas_stationarity_tol: self.conv_check.infeas_stationarity_tol,
678                infeas_viol_kappa: self.conv_check.infeas_viol_kappa,
679                infeas_max_streak: self.conv_check.infeas_max_streak,
680                infeas_streak: 0,
681            });
682
683        let init: Box<dyn crate::init::r#trait::IterateInitializer> = if self.warm_start_init_point
684        {
685            Box::new(WarmStartIterateInitializer::with_options(self.warm.clone()))
686        } else {
687            let mut d = DefaultIterateInitializer::with_eq_mult_calculator(Box::new(
688                LeastSquareMults::new(),
689            ));
690            d.bound_push = self.init.bound_push;
691            d.bound_frac = self.init.bound_frac;
692            d.slack_bound_push = self.init.slack_bound_push;
693            d.slack_bound_frac = self.init.slack_bound_frac;
694            d.constr_mult_init_max = self.init.constr_mult_init_max;
695            d.bound_mult_init_val = self.init.bound_mult_init_val;
696            d.bound_mult_init_method = self.init.bound_mult_init_method.clone();
697            d.least_square_init_primal = self.init.least_square_init_primal;
698            Box::new(d)
699        };
700
701        let eq_mult: Box<dyn crate::eq_mult::r#trait::EqMultCalculator> =
702            Box::new(LeastSquareMults::new());
703
704        let hess: Box<dyn crate::hess::r#trait::HessianUpdater> = match self.hessian_approximation {
705            HessianApproxChoice::Exact => Box::new(ExactHessianUpdater::new()),
706            HessianApproxChoice::LimitedMemory => Box::new(LimMemQuasiNewtonUpdater {
707                update_type: self.limited_memory_update_type,
708                ..LimMemQuasiNewtonUpdater::default()
709            }),
710        };
711
712        let iter_output: Box<dyn crate::output::r#trait::IterationOutput> = {
713            use crate::output::orig::{InfPrTag, PrintInfoString};
714            let mut o = OrigIterationOutput::new();
715            o.print_frequency_iter = self.output.print_frequency_iter;
716            o.print_frequency_time = self.output.print_frequency_time;
717            o.print_info_string = if self.output.print_info_string {
718                PrintInfoString::Yes
719            } else {
720                PrintInfoString::No
721            };
722            o.inf_pr_output = if self.output.inf_pr_output_internal {
723                InfPrTag::Internal
724            } else {
725                InfPrTag::Original
726            };
727            Box::new(o)
728        };
729
730        AlgorithmBundle {
731            mu_update,
732            conv_check,
733            init,
734            eq_mult,
735            hess,
736            line_search,
737            iter_output,
738            search_dir,
739        }
740    }
741}
742
743#[cfg(test)]
744mod tests {
745    use super::*;
746
747    #[test]
748    fn default_builder_assembles() {
749        let bundle = AlgorithmBuilder::new().build();
750        // Sanity: the placeholder traits compile and the boxed
751        // strategies don't panic on construction.
752        let _ = bundle.line_search.acceptor();
753        assert!(bundle.search_dir.is_none());
754    }
755
756    #[test]
757    fn build_with_backend_assembles_search_dir_chain() {
758        // Drive the builder with the FERAL backend factory; the
759        // resulting bundle should expose a populated `PdSearchDirCalc`.
760        let factory: LinearBackendFactory = Box::new(|_| {
761            Box::new(pounce_feral::FeralSolverInterface::new())
762                as Box<dyn SparseSymLinearSolverInterface>
763        });
764        let bundle = AlgorithmBuilder::new().build_with_backend(factory);
765        assert!(bundle.search_dir.is_some());
766    }
767
768    #[test]
769    fn limited_memory_sr1_propagates() {
770        let b = AlgorithmBuilder {
771            hessian_approximation: HessianApproxChoice::LimitedMemory,
772            limited_memory_update_type: UpdateType::Sr1,
773            ..AlgorithmBuilder::default()
774        };
775        let _bundle = b.build();
776    }
777
778    #[test]
779    fn every_strategy_combination_assembles_without_panic() {
780        let solvers = [LinearSolverChoice::Ma57, LinearSolverChoice::Feral];
781        let mu = [MuStrategyChoice::Monotone, MuStrategyChoice::Adaptive];
782        let hess = [
783            HessianApproxChoice::Exact,
784            HessianApproxChoice::LimitedMemory,
785        ];
786        let ls = [
787            LineSearchChoice::Filter,
788            LineSearchChoice::CgPenalty,
789            LineSearchChoice::Penalty,
790        ];
791        for &linear_solver in &solvers {
792            for &mu_strategy in &mu {
793                for &hessian_approximation in &hess {
794                    for &line_search_method in &ls {
795                        let _ = AlgorithmBuilder {
796                            algorithm: AlgorithmChoice::default(),
797                            linear_solver,
798                            linear_system_scaling: LinearSystemScalingChoice::None,
799                            linear_scaling_on_demand: true,
800                            mu_strategy,
801                            mu_oracle: MuOracleKind::QualityFunction,
802                            hessian_approximation,
803                            limited_memory_update_type: UpdateType::Bfgs,
804                            line_search_method,
805                            warm_start_init_point: false,
806                            mehrotra_algorithm: false,
807                            conv_check: ConvCheckOptions::default(),
808                            mu: MuOptions::default(),
809                            line_search: LineSearchOptions::default(),
810                            output: OutputOptions::default(),
811                            warm: WarmStartOptions::default(),
812                            sqp: crate::sqp::SqpOptions::default(),
813                            init: InitOptions::default(),
814                        }
815                        .build();
816                    }
817                }
818            }
819        }
820    }
821}