Skip to main content

pounce_algorithm/mu/oracle/
probing.rs

1//! Probing (Mehrotra) mu oracle — port of
2//! `IpProbingMuOracle.{hpp,cpp}`. Phase 10.
3//!
4//! Given the affine-step complementarity `mu_aff` and current `mu`,
5//! the probing rule chooses
6//!
7//! ```text
8//!   sigma = min((mu_aff / mu)^3, sigma_max)
9//!   mu_new = sigma * mu
10//! ```
11//!
12//! per `IpProbingMuOracle.cpp:117-121` and the standard Mehrotra
13//! predictor-corrector formula.
14
15use crate::ipopt_cq::IpoptCqHandle;
16use crate::ipopt_data::IpoptDataHandle;
17use crate::ipopt_nlp::IpoptNlp;
18use crate::iterates_vector::IteratesVector;
19use crate::kkt::pd_search_dir_calc::PdSearchDirCalc;
20use crate::mu::oracle::r#trait::MuOracle;
21use pounce_common::types::Number;
22use std::cell::RefCell;
23use std::rc::Rc;
24
25pub struct ProbingMuOracle {
26    pub sigma_max: Number,
27    pub mu_min: Number,
28    pub mu_max: Number,
29    pub mu_curr: Number,
30    pub mu_aff: Number,
31}
32
33impl Default for ProbingMuOracle {
34    fn default() -> Self {
35        Self {
36            // Default `sigma_max = 100` per `RegisterOptions`.
37            sigma_max: 100.0,
38            mu_min: 1e-11,
39            mu_max: 1e5,
40            mu_curr: 1.0,
41            mu_aff: 1.0,
42        }
43    }
44}
45
46impl ProbingMuOracle {
47    pub fn new() -> Self {
48        Self::default()
49    }
50
51    /// Pure-arithmetic probing formula from `IpProbingMuOracle.cpp:117`.
52    pub fn probing_mu(mu_curr: Number, mu_aff: Number, sigma_max: Number) -> Number {
53        let sigma = (mu_aff / mu_curr).powi(3).min(sigma_max);
54        sigma * mu_curr
55    }
56
57    /// Drive the affine (predictor) step through `pd_search_dir`,
58    /// snapshot `μ_aff` from the projected complementarity, then return
59    /// the probed `μ_new ∈ [mu_min, mu_max]`. Returns `None` if the
60    /// affine linear solve fails (caller falls back to LOQO).
61    ///
62    /// `tau` is the fraction-to-the-boundary parameter (typically
63    /// `1.0` for the affine probe per upstream's
64    /// `IpProbingMuOracle.cpp:107`).
65    pub fn calculate_mu_with_affine_step(
66        &mut self,
67        data: &IpoptDataHandle,
68        cq: &IpoptCqHandle,
69        nlp: &Rc<RefCell<dyn IpoptNlp>>,
70        pd_search_dir: &mut PdSearchDirCalc,
71        tau: Number,
72    ) -> Option<Number> {
73        if !pd_search_dir.compute_affine_step(data, cq, nlp) {
74            return None;
75        }
76        let delta_aff: IteratesVector = data.borrow().delta_aff.clone()?;
77        let cq_ref = cq.borrow();
78        let mu_curr = data.borrow().curr_mu;
79        let alpha_pri = cq_ref.aff_step_alpha_primal_max(&delta_aff, tau);
80        let alpha_du = cq_ref.aff_step_alpha_dual_max(&delta_aff, tau);
81        let mu_aff = cq_ref.aff_step_compl_avrg(&delta_aff, alpha_pri, alpha_du);
82        self.mu_curr = mu_curr;
83        self.mu_aff = mu_aff;
84        let raw = Self::probing_mu(mu_curr, mu_aff, self.sigma_max);
85        Some(raw.clamp(self.mu_min, self.mu_max))
86    }
87}
88
89impl MuOracle for ProbingMuOracle {
90    fn calculate_mu(&mut self) -> Option<Number> {
91        let raw = Self::probing_mu(self.mu_curr, self.mu_aff, self.sigma_max);
92        Some(raw.clamp(self.mu_min, self.mu_max))
93    }
94}
95
96#[cfg(test)]
97mod tests {
98    use super::*;
99
100    #[test]
101    fn probing_aff_equals_curr_keeps_mu() {
102        // ratio=1 → sigma=1 → mu_new=mu_curr.
103        assert_eq!(ProbingMuOracle::probing_mu(1.0, 1.0, 100.0), 1.0);
104    }
105
106    #[test]
107    fn probing_aff_half_curr() {
108        // ratio=0.5 → sigma=0.125 → mu_new=0.125 * 1.0.
109        let m = ProbingMuOracle::probing_mu(1.0, 0.5, 100.0);
110        assert!((m - 0.125).abs() < 1e-15);
111    }
112
113    #[test]
114    fn probing_caps_at_sigma_max() {
115        // ratio=10 → sigma=1000, capped at 100. mu_new=100.
116        let m = ProbingMuOracle::probing_mu(1.0, 10.0, 100.0);
117        assert!((m - 100.0).abs() < 1e-13);
118    }
119
120    #[test]
121    fn calculate_mu_via_trait_clamped() {
122        let mut o = ProbingMuOracle {
123            sigma_max: 100.0,
124            mu_min: 0.5,
125            mu_max: 10.0,
126            mu_curr: 1.0,
127            mu_aff: 0.001, // sigma ≈ 1e-9, mu_new ≈ 1e-9, clamps to 0.5.
128        };
129        assert_eq!(o.calculate_mu(), Some(0.5));
130    }
131}