Skip to main content

pounce_presolve/
coupling.rs

1//! Coupling classification for candidate auxiliary blocks.
2//!
3//! PR 5 of the auxiliary-presolve port (issue #53). After PRs 2–4
4//! find independent square blocks, this module decides whether each
5//! one is **safe** to eliminate.
6//!
7//! A block is unsafe if it interacts with parts of the problem the
8//! IPM still owns:
9//!
10//! - **Inequality coupling** — a block variable appears in any
11//!   inequality row. Eliminating the block would fix that variable's
12//!   value, possibly violating the inequality.
13//! - **Objective coupling** — a block variable contributes to the
14//!   objective gradient. Eliminating fixes the value and removes the
15//!   IPM's freedom to move that variable to improve the objective.
16//!
17//! ripopt anchor: `src/auxiliary_preprocessing.rs:39-59, 1642-1687`.
18
19use std::collections::HashSet;
20
21use pounce_common::types::Number;
22
23use crate::btf::BlockTriangularBlock;
24use crate::incidence::InequalityIncidence;
25
26/// How a candidate block is coupled to the rest of the problem.
27///
28/// Drives the elimination policy: `PureEquality` is always eligible
29/// under `AuxiliaryCouplingPolicy::Safe`, `ObjectiveCoupled` adds in
30/// under `Aggressive`, and the two inequality-coupled variants are
31/// never eliminated in v1 (matches ripopt's conservative default).
32#[derive(Debug, Clone, Copy, PartialEq, Eq)]
33pub enum AuxiliaryCouplingClass {
34    /// Block touches only equality rows and not the objective.
35    PureEquality,
36    /// Block variables appear in the objective gradient.
37    ObjectiveCoupled,
38    /// Block variables appear in at least one inequality row.
39    InequalityCoupled,
40    /// Both of the above.
41    ObjectiveAndInequalityCoupled,
42}
43
44/// Return the set of variable indices where the objective gradient
45/// is non-negligible (`|grad_f[i]| > zero_tol`).
46pub fn objective_gradient_support(grad_f: &[Number], zero_tol: Number) -> HashSet<usize> {
47    grad_f
48        .iter()
49        .enumerate()
50        .filter_map(|(i, &g)| if g.abs() > zero_tol { Some(i) } else { None })
51        .collect()
52}
53
54/// Classify one candidate block's coupling to the rest of the
55/// problem.
56///
57/// # Example
58///
59/// ```
60/// use std::collections::HashSet;
61/// use pounce_presolve::btf::BlockTriangularBlock;
62/// use pounce_presolve::coupling::{classify_block, AuxiliaryCouplingClass};
63/// use pounce_presolve::incidence::{InequalityIncidence, ProbeView};
64///
65/// // Build a problem with one equality row (row 0) and one
66/// // inequality (row 1) touching var 1, and one variable with
67/// // non-zero objective gradient.
68/// let p = ProbeView {
69///     n_vars: 2,
70///     m_rows: 2,
71///     jac_irow: &[0, 0, 1],
72///     jac_jcol: &[0, 1, 1],
73///     jac_values: None,
74///     g_l: &[0.0, 0.0],
75///     g_u: &[0.0, 5.0],
76///     linearity: None,
77///     one_based: false,
78///     eq_tol: 1e-12,
79///     excluded_vars: None,
80///     excluded_rows: None,
81/// };
82/// let ineq = InequalityIncidence::from_probe(&p);
83/// let block = BlockTriangularBlock { eq_rows: vec![0], cols: vec![0, 1] };
84/// let obj: HashSet<usize> = [0usize].into_iter().collect();
85/// let c = classify_block(&block, &ineq, &obj);
86/// assert_eq!(c, AuxiliaryCouplingClass::ObjectiveAndInequalityCoupled);
87/// ```
88pub fn classify_block(
89    block: &BlockTriangularBlock,
90    inequalities: &InequalityIncidence,
91    obj_grad_support: &HashSet<usize>,
92) -> AuxiliaryCouplingClass {
93    let in_obj = block.cols.iter().any(|c| obj_grad_support.contains(c));
94    let in_ineq = block
95        .cols
96        .iter()
97        .any(|&c| inequalities.var_in_inequality(c));
98    match (in_obj, in_ineq) {
99        (false, false) => AuxiliaryCouplingClass::PureEquality,
100        (true, false) => AuxiliaryCouplingClass::ObjectiveCoupled,
101        (false, true) => AuxiliaryCouplingClass::InequalityCoupled,
102        (true, true) => AuxiliaryCouplingClass::ObjectiveAndInequalityCoupled,
103    }
104}
105
106#[cfg(test)]
107mod tests {
108    use super::*;
109    use crate::incidence::ProbeView;
110    use pounce_common::types::{Index, Number};
111
112    fn probe<'a>(
113        n_vars: usize,
114        m_rows: usize,
115        irow: &'a [Index],
116        jcol: &'a [Index],
117        g_l: &'a [Number],
118        g_u: &'a [Number],
119    ) -> ProbeView<'a> {
120        ProbeView {
121            n_vars,
122            m_rows,
123            jac_irow: irow,
124            jac_jcol: jcol,
125            jac_values: None,
126            g_l,
127            g_u,
128            linearity: None,
129            one_based: false,
130            eq_tol: 1e-12,
131            excluded_vars: None,
132            excluded_rows: None,
133        }
134    }
135
136    fn block(cols: Vec<usize>) -> BlockTriangularBlock {
137        BlockTriangularBlock {
138            eq_rows: cols.clone(),
139            cols,
140        }
141    }
142
143    #[test]
144    fn coupling_pure_equality() {
145        // 2 equality rows, no inequalities, zero objective gradient.
146        let p = probe(2, 2, &[0, 1], &[0, 1], &[0.0, 0.0], &[0.0, 0.0]);
147        let ineq = InequalityIncidence::from_probe(&p);
148        let obj = HashSet::new();
149        let b = block(vec![0, 1]);
150        assert_eq!(
151            classify_block(&b, &ineq, &obj),
152            AuxiliaryCouplingClass::PureEquality
153        );
154    }
155
156    #[test]
157    fn coupling_objective_only() {
158        let p = probe(2, 2, &[0, 1], &[0, 1], &[0.0, 0.0], &[0.0, 0.0]);
159        let ineq = InequalityIncidence::from_probe(&p);
160        let obj: HashSet<usize> = [0].into_iter().collect();
161        let b = block(vec![0, 1]);
162        assert_eq!(
163            classify_block(&b, &ineq, &obj),
164            AuxiliaryCouplingClass::ObjectiveCoupled
165        );
166    }
167
168    #[test]
169    fn coupling_inequality_only() {
170        // Row 1 is g ∈ [0, 5] (inequality) and touches col 0.
171        let p = probe(2, 2, &[0, 0, 1], &[0, 1, 0], &[0.0, 0.0], &[0.0, 5.0]);
172        let ineq = InequalityIncidence::from_probe(&p);
173        let obj = HashSet::new();
174        let b = block(vec![0, 1]);
175        assert_eq!(
176            classify_block(&b, &ineq, &obj),
177            AuxiliaryCouplingClass::InequalityCoupled
178        );
179    }
180
181    #[test]
182    fn coupling_both() {
183        let p = probe(2, 2, &[0, 0, 1], &[0, 1, 0], &[0.0, 0.0], &[0.0, 5.0]);
184        let ineq = InequalityIncidence::from_probe(&p);
185        let obj: HashSet<usize> = [1].into_iter().collect();
186        let b = block(vec![0, 1]);
187        assert_eq!(
188            classify_block(&b, &ineq, &obj),
189            AuxiliaryCouplingClass::ObjectiveAndInequalityCoupled
190        );
191    }
192
193    #[test]
194    fn coupling_partial_block_inequality() {
195        // Inequality row touches only var 1; block has both 0 and 1.
196        // Any coupling counts — block is InequalityCoupled.
197        let p = probe(2, 2, &[0, 0, 1], &[0, 1, 1], &[0.0, 0.0], &[0.0, 5.0]);
198        let ineq = InequalityIncidence::from_probe(&p);
199        let obj = HashSet::new();
200        let b = block(vec![0, 1]);
201        assert_eq!(
202            classify_block(&b, &ineq, &obj),
203            AuxiliaryCouplingClass::InequalityCoupled
204        );
205    }
206
207    #[test]
208    fn coupling_zero_grad_tolerance() {
209        // Variable 0 has |grad_f| = 1e-16 < zero_tol. Should not be
210        // flagged as in-objective.
211        let support = objective_gradient_support(&[1e-16, 0.0], 1e-12);
212        assert!(support.is_empty());
213        let p = probe(2, 1, &[0], &[0], &[0.0], &[0.0]);
214        let ineq = InequalityIncidence::from_probe(&p);
215        let b = block(vec![0, 1]);
216        assert_eq!(
217            classify_block(&b, &ineq, &support),
218            AuxiliaryCouplingClass::PureEquality
219        );
220    }
221
222    #[test]
223    fn coupling_empty_block() {
224        let p = probe(2, 1, &[0], &[0], &[0.0], &[0.0]);
225        let ineq = InequalityIncidence::from_probe(&p);
226        let obj: HashSet<usize> = [0, 1].into_iter().collect();
227        let b = BlockTriangularBlock::default();
228        assert_eq!(
229            classify_block(&b, &ineq, &obj),
230            AuxiliaryCouplingClass::PureEquality
231        );
232    }
233
234    #[test]
235    fn coupling_grad_support_helper() {
236        let support = objective_gradient_support(&[2.5, 0.0, -1e-3, 1e-15], 1e-9);
237        let expected: HashSet<usize> = [0, 2].into_iter().collect();
238        assert_eq!(support, expected);
239    }
240}