pounce_qp/problem.rs
1//! QP problem definition and the solution / warm-start types it
2//! pairs with. Storage borrows from `pounce-linalg` types directly;
3//! the QP solver never copies the Hessian or Jacobian.
4
5use crate::error::{QpError, QpStatus};
6use crate::working_set::WorkingSet;
7use pounce_common::Number;
8use pounce_linalg::triplet::{GenTMatrix, SymTMatrix};
9use std::time::Duration;
10
11/// Caller-supplied hint about the inertia of `H`. Lets a strictly-
12/// convex problem skip the inertia-correction probe of §4.5.
13/// `Unknown` is always safe (the solver detects indefiniteness from
14/// the LDLᵀ factor of the KKT block) and is the default.
15#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
16pub enum HessianInertia {
17 /// `H` is symmetric positive semi-definite.
18 Psd,
19 /// `H` has (potentially) negative eigenvalues; the inertia-
20 /// control path is required.
21 Indefinite,
22 /// Caller offers no claim; solver probes via factor inertia.
23 #[default]
24 Unknown,
25}
26
27/// A convex-or-nonconvex sparse QP:
28/// ```text
29/// min ½ xᵀ H x + gᵀ x
30/// s.t. bl ≤ A x ≤ bu
31/// xl ≤ x ≤ xu
32/// ```
33/// Two-sided general bounds (equality is `bl = bu`); two-sided
34/// variable bounds (fixed is `xl = xu`; free is `±NLP_*_BOUND_INF`).
35/// `H` is symmetric — only the upper triangle is stored.
36///
37/// The lifetime parameter ties the borrowed problem data to the
38/// caller's storage; the solver never copies any of these fields.
39pub struct QpProblem<'a> {
40 pub n: usize,
41 pub m: usize,
42 pub h: &'a SymTMatrix,
43 pub g: &'a [Number],
44 pub a: &'a GenTMatrix,
45 pub bl: &'a [Number],
46 pub bu: &'a [Number],
47 pub xl: &'a [Number],
48 pub xu: &'a [Number],
49 pub hessian_inertia: HessianInertia,
50}
51
52impl<'a> QpProblem<'a> {
53 /// Validate every dimension and bound-ordering invariant the
54 /// solver relies on. Called once at the top of `solve` before
55 /// any work happens.
56 pub fn validate(&self) -> Result<(), QpError> {
57 if self.h.space().dim() as usize != self.n {
58 return Err(QpError::DimensionMismatch(format!(
59 "H is {}×{} but n = {}",
60 self.h.space().dim(),
61 self.h.space().dim(),
62 self.n
63 )));
64 }
65 if self.g.len() != self.n {
66 return Err(QpError::DimensionMismatch(format!(
67 "g.len() = {} but n = {}",
68 self.g.len(),
69 self.n
70 )));
71 }
72 if self.a.space().n_rows() as usize != self.m || self.a.space().n_cols() as usize != self.n
73 {
74 return Err(QpError::DimensionMismatch(format!(
75 "A is {}×{} but expected {}×{}",
76 self.a.space().n_rows(),
77 self.a.space().n_cols(),
78 self.m,
79 self.n
80 )));
81 }
82 if self.bl.len() != self.m || self.bu.len() != self.m {
83 return Err(QpError::DimensionMismatch(format!(
84 "bl.len() = {}, bu.len() = {}, but m = {}",
85 self.bl.len(),
86 self.bu.len(),
87 self.m
88 )));
89 }
90 if self.xl.len() != self.n || self.xu.len() != self.n {
91 return Err(QpError::DimensionMismatch(format!(
92 "xl.len() = {}, xu.len() = {}, but n = {}",
93 self.xl.len(),
94 self.xu.len(),
95 self.n
96 )));
97 }
98 for (i, (&l, &u)) in self.bl.iter().zip(self.bu.iter()).enumerate() {
99 if l > u {
100 return Err(QpError::InvertedBounds(format!(
101 "constraint row {i}: bl = {l} > bu = {u}"
102 )));
103 }
104 }
105 for (i, (&l, &u)) in self.xl.iter().zip(self.xu.iter()).enumerate() {
106 if l > u {
107 return Err(QpError::InvertedBounds(format!(
108 "variable {i}: xl = {l} > xu = {u}"
109 )));
110 }
111 }
112 Ok(())
113 }
114}
115
116/// Warm-start seed: previous primal-dual iterate plus working set.
117/// Passed to [`crate::QpSolver::solve`] as `Some(ws)`; `None` is the
118/// cold-start path through phase-1 elastic mode (§4.3).
119#[derive(Debug, Clone)]
120pub struct QpWarmStart {
121 pub x: Vec<Number>,
122 /// Lagrange multipliers for the general constraints, length `m`.
123 pub lambda_g: Vec<Number>,
124 /// Bound multipliers, length `n`, packed signed
125 /// (`z_l − z_u`). Positive ⇒ lower-bound active, negative ⇒
126 /// upper-bound active.
127 pub lambda_x: Vec<Number>,
128 pub working: WorkingSet,
129}
130
131/// Solver output. `working` is the new working set, suitable for
132/// passing as the next solve's warm start.
133#[derive(Debug, Clone)]
134pub struct QpSolution {
135 pub x: Vec<Number>,
136 pub lambda_g: Vec<Number>,
137 pub lambda_x: Vec<Number>,
138 pub working: WorkingSet,
139 pub obj: Number,
140 pub status: QpStatus,
141 pub stats: QpStats,
142}
143
144/// Per-solve counters and timers reported alongside the solution.
145/// Phase 5a uses these for the §8.2 scaling-sweep plots and the
146/// §8.5 warm-start sweep.
147#[derive(Debug, Clone, Default)]
148pub struct QpStats {
149 /// Total active-set changes (adds + drops) across the solve.
150 pub n_working_set_changes: u32,
151 /// Number of times the cached factorization was discarded and
152 /// the base KKT was refactored (the §4.2 reset cycle).
153 pub n_refactor: u32,
154 /// Number of Schur-complement rank-1 updates applied (the
155 /// bounded-cost path between refactors).
156 pub n_schur_updates: u32,
157 /// Whether the solve passed through phase-1 elastic mode.
158 pub used_phase1: bool,
159 /// Wall-clock time spent inside `solve`.
160 pub time: Duration,
161}