otspot_core/problem/
mod.rs

1//! LP問題定義モジュール
2//!
3//! 線形計画問題（LP）の構造定義・制約種別・ソルバー結果の表現を提供する。
4//! 問題は標準形 `min c^T x  s.t.  Ax {<=,>=,=} b,  x in [lb, ub]` で定義される。
5
6pub mod certificate;
7use certificate::{BoundGapCertificate, OptimalCertificate};
8
9use crate::error::SolverError;
10use crate::options::WarmStartBasis;
11use crate::sparse::CscMatrix;
12use std::fmt;
13
14/// Returns `true` if `(lb, ub)` is a valid bound pair.
15///
16/// Valid means: neither endpoint is NaN, `lb <= ub`, `lb < +∞`, and `ub > -∞`.
17/// The last two conditions reject empty intervals that `lb > ub` misses when both
18/// endpoints are the same infinity — e.g., `(+∞, +∞)` or `(-∞, -∞)`.
19pub(crate) fn is_valid_bound_pair(lb: f64, ub: f64) -> bool {
20    !lb.is_nan() && !ub.is_nan() && lb <= ub && lb < f64::INFINITY && ub > f64::NEG_INFINITY
21}
22
23/// LP問題における制約条件の種別
24#[non_exhaustive]
25#[derive(Debug, Clone, Copy, PartialEq)]
26pub enum ConstraintType {
27    /// 以下（<=）
28    Le,
29    /// 以上（>=）
30    Ge,
31    /// 等式（==）
32    Eq,
33}
34
35/// Route taken by a solve call (populated per-result, race-free).
36#[non_exhaustive]
37#[derive(Debug, Clone, Copy, PartialEq, Default)]
38pub enum SolveRoute {
39    /// Route not yet set (default for uninitialized results).
40    #[default]
41    Unknown,
42    /// Called directly via `crate::lp::solve_lp_with`.
43    LpDirect,
44    /// LP forwarded from `solve_qp_with(Q=0)`.
45    LpForwardedFromQp,
46    /// QP solved via IPM (Q≠0).
47    QpIpm,
48}
49
50/// Per-solve routing and warm-start statistics (race-free, per-result).
51///
52/// Replaces process-global `AtomicU64` counters so parallel tests observe
53/// independent stats without reset/race issues.
54#[derive(Debug, Clone, Default)]
55pub struct SolveStats {
56    /// Route taken for this solve.
57    pub route: SolveRoute,
58    /// Whether the solver stopped because the deadline (timeout_secs / deadline) was reached.
59    ///
60    /// `true` iff `result.status == SolveStatus::Timeout`. Deterministic sentinel for
61    /// deadline-enforcement tests: assert this field instead of measuring wall time.
62    pub deadline_triggered: bool,
63    /// Whether the postsolve saddle-point Krylov IR was skipped because the solution
64    /// already met the user tolerance (`kkt_already_pass`). Deterministic sentinel for
65    /// the gate: removing the gate (always refine) flips this to `false`.
66    pub postsolve_krylov_ir_skipped: bool,
67    /// For LP solves via the QP dispatch (`LpForwardedFromQp`): `true` if the
68    /// returned result came from the IPM, `false` if from simplex. Lets callers
69    /// (e.g. benchmarks) label the actual route instead of a static size guess.
70    pub lp_ipm_path: bool,
71    /// LP simplex used the bounded-variable Eq+finite-UB Phase I path.
72    ///
73    /// This is a load-bearing route sentinel for large Netlib LPs whose upper
74    /// bounds must stay as variable bounds instead of expanded UB rows.
75    pub bounded_eq_ub_path: bool,
76}
77
78/// ソルバーの求解結果ステータス
79#[non_exhaustive]
80#[derive(Debug, Clone, PartialEq)]
81pub enum SolveStatus {
82    /// 最適解が求まった
83    Optimal,
84    /// 局所的最適解（非凸QP: 慣性修正付きIPMが収束したKKT点）
85    ///
86    /// Q行列が不定（indefinite）の場合、大域的最適性は保証されないが、
87    /// KKT条件を満たす局所最適解またはサドル点が返される。
88    /// 慣性修正（Gershgorin 評価から導出した δI 加算）により IPM を収束させた。
89    LocallyOptimal,
90    /// 問題が実行不可能（infeasible）
91    Infeasible,
92    /// 問題が非有界（unbounded）
93    Unbounded,
94    /// 反復回数上限に到達した（最適性未確認）
95    MaxIterations,
96    /// 解は見つかったが精度基準未達（偽Optimal検出: スケール解除後の残差超過）
97    SuboptimalSolution,
98    /// タイムアウト（timeout_secs を超過した）
99    Timeout,
100    /// 数値エラー（LDL分解失敗等、問題が数値的に解けない）
101    NumericalError,
102    /// Q行列が不定（非凸QP）。IPMはQ正半定値を前提とする。
103    NonConvex(String),
104    /// 非凸 QP の局所最適解 (= `solve_qp_global` 経由で incumbent あり、ε-global 証明なし)。
105    ///
106    /// BB driver が deadline / max_nodes / max_depth で打ち切られ、incumbent ある状態。
107    /// `LocallyOptimal` (= IPM inertia 補正後の単発解) と区別して、caller が「探索打切」
108    /// vs「単発 KKT 収束」を識別できる。`Optimal` には**含めない** (= global proof なし)。
109    NonconvexLocal,
110    /// 非凸 QP の大域 ε-最適解 (= `solve_qp_global` で gap_tol まで証明済み + Q が indefinite)。
111    ///
112    /// `Optimal` は「Q が PSD で IPM/BB が global 達成」専用に維持し、indefinite Q の場合は
113    /// 本 variant で明示分離する (caller が「global 証明済」かを fact で判別)。
114    NonconvexGlobal,
115    /// Problem type not supported by this solver.
116    ///
117    /// Returned when the caller passes a problem that cannot be handled, e.g.
118    /// a QCQP (quadratic constraints present) submitted to the QP/LP entry.
119    NotSupported(String),
120}
121
122impl fmt::Display for SolveStatus {
123    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
124        match self {
125            SolveStatus::Optimal => write!(f, "Optimal"),
126            SolveStatus::LocallyOptimal => write!(f, "LocallyOptimal"),
127            SolveStatus::Infeasible => write!(f, "Infeasible"),
128            SolveStatus::Unbounded => write!(f, "Unbounded"),
129            SolveStatus::MaxIterations => write!(f, "MaxIterations"),
130            SolveStatus::SuboptimalSolution => write!(f, "SuboptimalSolution"),
131            SolveStatus::Timeout => write!(f, "Timeout"),
132            SolveStatus::NumericalError => write!(f, "NumericalError"),
133            SolveStatus::NonConvex(msg) => write!(f, "NonConvex({})", msg),
134            SolveStatus::NonconvexLocal => write!(f, "NonconvexLocal"),
135            SolveStatus::NonconvexGlobal => write!(f, "NonconvexGlobal"),
136            SolveStatus::NotSupported(msg) => write!(f, "NotSupported({})", msg),
137        }
138    }
139}
140
141/// LP/QP共通求解結果型
142///
143/// LP求解（Simplex等）と QP求解（AS/IPM/Concurrent）の両方で使用できる統一結果型。
144/// LP固有フィールド（`reduced_costs`, `slack`, `warm_start_basis`）は QP求解時は空/None。
145/// QP固有フィールド（`bound_duals`, `iterations`）は LP求解時は空/0。
146#[derive(Debug, Clone)]
147pub struct SolverResult {
148    /// 求解ステータス
149    pub status: SolveStatus,
150    /// 最適目的関数値（最適解が存在する場合）
151    pub objective: f64,
152    /// 解ベクトル（最適解が存在する場合）
153    pub solution: Vec<f64>,
154    /// 双対変数ベクトル（各制約の影価格、最適解が存在する場合）
155    pub dual_solution: Vec<f64>,
156    // --- LP固有フィールド ---
157    /// 被縮小費用ベクトル（各決定変数に対して、最適解が存在する場合）
158    pub reduced_costs: Vec<f64>,
159    /// スラック変数ベクトル（各制約のスラック b_i - a_i^T x、最適解が存在する場合）
160    pub slack: Vec<f64>,
161    /// warm-start用の基底情報（Optimal時のみ Some）
162    pub warm_start_basis: Option<WarmStartBasis>,
163    // --- QP固有フィールド ---
164    /// Bound dual values (shadow prices for variable bounds).
165    ///
166    /// Maps to original variable indices via col_map.
167    /// Empty if no bound constraints are active.
168    ///
169    /// - 除去変数 (presolveで固定された変数) の bound_dual = 0.0 (近似)
170    /// - presolve tightening で追加された境界の dual は報告しない（元問題基準）
171    /// - 配列順: `[lb_dual(j0), ..., lb_dual(j_{n_lb-1}), ub_dual(j0), ..., ub_dual(j_{n_ub-1})]`
172    pub bound_duals: Vec<f64>,
173    /// 反復回数（WSR実績回数）
174    pub iterations: usize,
175    /// 最終反復の残差実値 (pfeas, dfeas, duality_gap)。Optimal/MaxIterations時のみ Some。
176    pub final_residuals: Option<(f64, f64, f64)>,
177    /// 相対双対ギャップ (|p_obj - d_obj| / max(|p|,|d|,1))。
178    /// IPPMM 内部の best-so-far に紐づく値。unscale_ipm_result の Suboptimal→Optimal 昇格ゲート用。
179    /// None = 未計測（LP simplex 等 gap を持たない経路）。
180    pub duality_gap_rel: Option<f64>,
181    /// 各 phase の所要時間 (LP simplex 経路のみ、None なら未計測)。
182    /// 「どこに時間が掛かっているか」事実観測用 (CLAUDE.md「順調に収束に向けて探索」)。
183    pub timing_breakdown: Option<TimingBreakdown>,
184    /// Postsolve が最終的に採用した y_orig の dfeas violation (bound-aware sup ノルム).
185    /// LP simplex + presolve 経路のみ Some。caller (solve_with) が値が `PIVOT_TOL` を
186    /// 超えるとき presolve=off で再解する fallback gate に使う (greenbea-class 問題対策)。
187    pub postsolve_dfeas: Option<f64>,
188    /// Per-solve routing and warm-start statistics (race-free).
189    pub stats: SolveStats,
190    /// Branch-and-bound gap certificate.
191    ///
192    /// Present iff the solver completed a B&B search with a fully authenticated gap
193    /// (no `proof_uncertain` region, `within_gap` satisfied). `None` for direct LP/QP solves.
194    pub bound_gap_cert: Option<BoundGapCertificate>,
195    /// KKT optimality certificate — minted by `prove_optimal` on the B&B incumbent.
196    ///
197    /// Set when `finalize_proven` verifies all KKT conditions on the returned point.
198    /// `None` for demoted (LocallyOptimal/NonconvexLocal) or non-B&B results.
199    pub opt_cert: Option<OptimalCertificate>,
200}
201
202/// 各 phase 所要時間 (μs精度)。LP simplex と QP IPM の両経路で共用。
203#[derive(Debug, Clone, Copy, Default, PartialEq)]
204pub struct TimingBreakdown {
205    // ── LP simplex 経路 ───────────────────────────────────────────────────────
206    /// Presolve 全体 (run_presolve)
207    pub presolve_us: u64,
208    /// 縮約後 simplex 本体 (solve_without_presolve)
209    pub solve_us: u64,
210    /// Postsolve 総計 (run_postsolve / QP 経路では下記 4 field の合計)
211    pub postsolve_us: u64,
212
213    // ── QP IPM 経路: IPM 反復内部 ────────────────────────────────────────────
214    /// KKT 行列 LDL 数値因子化の累計 (全 iteration 合計)
215    pub ipm_factorize_us: u64,
216    /// KKT solve (predictor/corrector/Gondzio) の累計
217    pub ipm_solve_us: u64,
218    /// LDL regularization retry の累計回数 (健全性プローブ失敗含む)
219    pub ipm_reg_retries: u32,
220    /// MINRES (iterative) backend が 1 回以上使われたか
221    pub ipm_used_iterative: bool,
222
223    // ── QP IPM 経路: postsolve 段階別 ────────────────────────────────────────
224    /// postsolve_qp_with_dual_recovery (reduced → orig 空間写像)
225    pub postsolve_map_us: u64,
226    /// refine_postsolve_dual_lsq (元空間 y LSQ refine)
227    pub postsolve_lsq_us: u64,
228    /// refine_postsolve_recovery (Stage 0: SingletonRow 後退代入)
229    pub postsolve_recovery_us: u64,
230    /// refine_post_processing (Stage 1+2: primal projection / y-z refit)
231    pub postsolve_refine_us: u64,
232    /// refine_krylov_and_projection (saddle-point Krylov IR)
233    pub postsolve_krylov_ir_us: u64,
234}
235
236impl Default for SolverResult {
237    fn default() -> Self {
238        SolverResult {
239            status: SolveStatus::NumericalError,
240            objective: 0.0,
241            solution: vec![],
242            dual_solution: vec![],
243            reduced_costs: vec![],
244            slack: vec![],
245            warm_start_basis: None,
246            bound_duals: vec![],
247            iterations: 0,
248            final_residuals: None,
249            duality_gap_rel: None,
250            timing_breakdown: None,
251            postsolve_dfeas: None,
252            stats: SolveStats::default(),
253            bound_gap_cert: None,
254            opt_cert: None,
255        }
256    }
257}
258
259impl fmt::Display for SolverResult {
260    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
261        write!(f, "Status: {}, Objective: {}", self.status, self.objective)
262    }
263}
264
265/// 線形計画問題: min c^T x + obj_offset  s.t.  Ax {op} b,  x in [lb, ub]
266///
267/// 目的関数・制約行列・右辺ベクトル・変数上下限をまとめて保持する。
268/// 制約種別（`<=`, `>=`, `=`）と変数ごとの上下限を個別に指定できる。
269///
270/// `obj_offset` は MPS/QPS の N 行 RHS (objective constant) に対応する。
271/// `LpProblem::new` / `new_general` では 0.0 に初期化される。
272#[derive(Debug, Clone)]
273pub struct LpProblem {
274    /// 目的関数係数ベクトル（長さ: `num_vars`）
275    pub c: Vec<f64>,
276    /// 制約行列（CSC形式、サイズ: `num_constraints` x `num_vars`）
277    pub a: CscMatrix,
278    /// 制約右辺ベクトル（長さ: `num_constraints`）
279    pub b: Vec<f64>,
280    /// 決定変数の数
281    pub num_vars: usize,
282    /// 制約式の数
283    pub num_constraints: usize,
284    /// 各制約の種別（長さ: `num_constraints`）
285    pub constraint_types: Vec<ConstraintType>,
286    /// 各変数の上下限 `(lower, upper)`（長さ: `num_vars`）
287    pub bounds: Vec<(f64, f64)>,
288    /// 問題名（オプション）
289    pub name: Option<String>,
290    /// 目的定数項 (MPS N-row RHS); 最適値 = c^T x* + obj_offset
291    pub obj_offset: f64,
292}
293
294impl LpProblem {
295    /// 新しいLP問題を検証付きで生成する（後方互換版）
296    ///
297    /// 標準形 `min c^T x  s.t.  Ax <= b,  x >= 0` を作成する。
298    /// 全制約を `<=`、全変数の下限を 0・上限を `+∞` とする。
299    ///
300    /// # 引数
301    /// * `c` - 目的関数係数ベクトル
302    /// * `a` - 制約行列（CSC形式）
303    /// * `b` - 制約右辺ベクトル
304    ///
305    /// # 戻り値
306    /// * `Ok(LpProblem)` - 次元が有効な場合
307    /// * `Err(String)` - 次元不一致などの検証エラー時
308    pub fn new(c: Vec<f64>, a: CscMatrix, b: Vec<f64>) -> Result<Self, SolverError> {
309        let num_vars = c.len();
310        let num_constraints = b.len();
311
312        // Set defaults for backward compatibility
313        let constraint_types = vec![ConstraintType::Le; num_constraints];
314        let bounds = vec![(0.0, f64::INFINITY); num_vars];
315        let name = None;
316
317        Self::new_general(c, a, b, constraint_types, bounds, name)
318    }
319
320    /// 制約種別と変数上下限を完全指定して新しいLP問題を生成する
321    ///
322    /// # 引数
323    /// * `c` - 目的関数係数ベクトル
324    /// * `a` - 制約行列（CSC形式）
325    /// * `b` - 制約右辺ベクトル
326    /// * `constraint_types` - 各制約の種別（`Le` / `Ge` / `Eq`）
327    /// * `bounds` - 各変数の上下限 `(lower, upper)`
328    /// * `name` - 問題名（オプション）
329    ///
330    /// # 戻り値
331    /// * `Ok(LpProblem)` - 次元が有効な場合
332    /// * `Err(String)` - 次元不一致などの検証エラー時
333    pub fn new_general(
334        c: Vec<f64>,
335        a: CscMatrix,
336        b: Vec<f64>,
337        constraint_types: Vec<ConstraintType>,
338        bounds: Vec<(f64, f64)>,
339        name: Option<String>,
340    ) -> Result<Self, SolverError> {
341        // Validate dimensions
342        if c.len() != a.ncols {
343            return Err(SolverError::DimensionMismatch {
344                field: "c",
345                expected: a.ncols,
346                got: c.len(),
347            });
348        }
349        if b.len() != a.nrows {
350            return Err(SolverError::DimensionMismatch {
351                field: "b",
352                expected: a.nrows,
353                got: b.len(),
354            });
355        }
356        if constraint_types.len() != b.len() {
357            return Err(SolverError::DimensionMismatch {
358                field: "constraint_types",
359                expected: b.len(),
360                got: constraint_types.len(),
361            });
362        }
363        if bounds.len() != c.len() {
364            return Err(SolverError::DimensionMismatch {
365                field: "bounds",
366                expected: c.len(),
367                got: bounds.len(),
368            });
369        }
370        for (i, &v) in c.iter().enumerate() {
371            if !v.is_finite() {
372                return Err(SolverError::NonFiniteCoefficient {
373                    field: "c",
374                    index: i,
375                });
376            }
377        }
378        for (i, &v) in b.iter().enumerate() {
379            if !v.is_finite() {
380                return Err(SolverError::NonFiniteCoefficient {
381                    field: "b",
382                    index: i,
383                });
384            }
385        }
386        for (i, &v) in a.values.iter().enumerate() {
387            if !v.is_finite() {
388                return Err(SolverError::NonFiniteCoefficient {
389                    field: "A",
390                    index: i,
391                });
392            }
393        }
394        for (i, &(lb, ub)) in bounds.iter().enumerate() {
395            if !is_valid_bound_pair(lb, ub) {
396                return Err(SolverError::InvalidBounds { index: i, lb, ub });
397            }
398        }
399
400        Ok(LpProblem {
401            num_vars: c.len(),
402            num_constraints: b.len(),
403            c,
404            a,
405            b,
406            constraint_types,
407            bounds,
408            name,
409            obj_offset: 0.0,
410        })
411    }
412}
413
414impl fmt::Display for LpProblem {
415    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
416        write!(
417            f,
418            "LP: min c^T x, {} vars, {} constraints",
419            self.num_vars, self.num_constraints
420        )
421    }
422}
423
424#[cfg(test)]
425mod tests {
426    use super::*;
427    use crate::error::SolverError;
428
429    #[test]
430    fn test_lp_problem_new_valid() {
431        // 2 variables, 2 constraints
432        let c = vec![1.0, 2.0];
433        let a = CscMatrix::new(2, 2);
434        let b = vec![5.0, 6.0];
435
436        let lp = LpProblem::new(c, a, b).unwrap();
437        assert_eq!(lp.num_vars, 2);
438        assert_eq!(lp.num_constraints, 2);
439    }
440
441    #[test]
442    fn test_lp_problem_new_invalid_c_dimension() {
443        // c.len() = 3, but a.ncols = 2
444        let c = vec![1.0, 2.0, 3.0];
445        let a = CscMatrix::new(2, 2);
446        let b = vec![5.0, 6.0];
447
448        let result = LpProblem::new(c, a, b);
449        assert!(result.is_err());
450        assert!(matches!(
451            result.unwrap_err(),
452            SolverError::DimensionMismatch { field: "c", .. }
453        ));
454    }
455
456    #[test]
457    fn test_lp_problem_new_invalid_b_dimension() {
458        // b.len() = 3, but a.nrows = 2
459        let c = vec![1.0, 2.0];
460        let a = CscMatrix::new(2, 2);
461        let b = vec![5.0, 6.0, 7.0];
462
463        let result = LpProblem::new(c, a, b);
464        assert!(result.is_err());
465        assert!(matches!(
466            result.unwrap_err(),
467            SolverError::DimensionMismatch { field: "b", .. }
468        ));
469    }
470
471    #[test]
472    fn test_lp_problem_display() {
473        let c = vec![1.0, 2.0];
474        let a = CscMatrix::new(2, 2);
475        let b = vec![5.0, 6.0];
476        let lp = LpProblem::new(c, a, b).unwrap();
477
478        let display = format!("{}", lp);
479        assert_eq!(display, "LP: min c^T x, 2 vars, 2 constraints");
480    }
481
482    #[test]
483    fn test_solve_status_display() {
484        assert_eq!(format!("{}", SolveStatus::Optimal), "Optimal");
485        assert_eq!(format!("{}", SolveStatus::Infeasible), "Infeasible");
486        assert_eq!(format!("{}", SolveStatus::Unbounded), "Unbounded");
487    }
488
489    #[test]
490    fn test_solver_result_display() {
491        let result = SolverResult {
492            status: SolveStatus::Optimal,
493            objective: 42.5,
494            solution: vec![1.0, 2.0],
495            dual_solution: vec![],
496            reduced_costs: vec![],
497            slack: vec![],
498            warm_start_basis: None,
499            ..Default::default()
500        };
501        let display = format!("{}", result);
502        assert_eq!(display, "Status: Optimal, Objective: 42.5");
503    }
504
505    #[test]
506    fn solver_result_default_is_not_success() {
507        let result = SolverResult::default();
508        assert_eq!(result.status, SolveStatus::NumericalError);
509        assert!(result.solution.is_empty());
510    }
511
512    fn make_lp(
513        c: Vec<f64>,
514        b: Vec<f64>,
515        a_vals: Vec<f64>,
516        bounds: Vec<(f64, f64)>,
517    ) -> Result<LpProblem, SolverError> {
518        let n = c.len();
519        let m = b.len();
520        let a = if a_vals.is_empty() {
521            CscMatrix::new(m, n)
522        } else {
523            let rows = vec![0usize; n];
524            let cols: Vec<usize> = (0..n).collect();
525            CscMatrix::from_triplets(&rows, &cols, &a_vals, m, n).unwrap()
526        };
527        let ct = vec![ConstraintType::Le; m];
528        LpProblem::new_general(c, a, b, ct, bounds, None)
529    }
530
531    #[test]
532    fn lp_valid_accepted() {
533        let res = make_lp(
534            vec![1.0, 2.0],
535            vec![5.0],
536            vec![1.0, 1.0],
537            vec![(0.0, f64::INFINITY), (0.0, 10.0)],
538        );
539        assert!(res.is_ok());
540    }
541
542    #[test]
543    fn lp_nan_in_c_rejected() {
544        let bad_vals = [f64::NAN, f64::INFINITY, f64::NEG_INFINITY];
545        for bad in bad_vals {
546            let res = make_lp(
547                vec![bad, 1.0],
548                vec![5.0],
549                vec![1.0, 1.0],
550                vec![(0.0, f64::INFINITY); 2],
551            );
552            assert!(
553                matches!(
554                    res,
555                    Err(SolverError::NonFiniteCoefficient { field: "c", .. })
556                ),
557                "expected NonFiniteCoefficient for c={bad}"
558            );
559        }
560    }
561
562    #[test]
563    fn lp_nan_in_b_rejected() {
564        let bad_vals = [f64::NAN, f64::INFINITY, f64::NEG_INFINITY];
565        for bad in bad_vals {
566            let res = make_lp(
567                vec![1.0, 2.0],
568                vec![bad],
569                vec![1.0, 1.0],
570                vec![(0.0, f64::INFINITY); 2],
571            );
572            assert!(
573                matches!(
574                    res,
575                    Err(SolverError::NonFiniteCoefficient { field: "b", .. })
576                ),
577                "expected NonFiniteCoefficient for b={bad}"
578            );
579        }
580    }
581
582    #[test]
583    fn lp_nan_in_a_rejected() {
584        let n = 2;
585        let bad_vals = [f64::NAN, f64::INFINITY, f64::NEG_INFINITY];
586        for bad in bad_vals {
587            // from_triplets drops NaN via DROP_TOL; inject bad value directly.
588            let mut a = CscMatrix::from_triplets(&[0], &[0], &[1.0], 1, n).unwrap();
589            a.values[0] = bad;
590            let res = LpProblem::new_general(
591                vec![1.0, 2.0],
592                a,
593                vec![5.0],
594                vec![ConstraintType::Le],
595                vec![(0.0, f64::INFINITY); n],
596                None,
597            );
598            assert!(
599                matches!(
600                    res,
601                    Err(SolverError::NonFiniteCoefficient { field: "A", .. })
602                ),
603                "expected NonFiniteCoefficient for A val={bad}"
604            );
605        }
606    }
607
608    #[test]
609    fn lp_nan_in_bounds_rejected() {
610        let cases: Vec<(f64, f64)> = vec![(f64::NAN, 1.0), (0.0, f64::NAN), (f64::NAN, f64::NAN)];
611        for (lb, ub) in cases {
612            let res = make_lp(
613                vec![1.0, 2.0],
614                vec![5.0],
615                vec![1.0, 1.0],
616                vec![(lb, ub), (0.0, f64::INFINITY)],
617            );
618            assert!(
619                matches!(res, Err(SolverError::InvalidBounds { index: 0, .. })),
620                "expected InvalidBounds for ({lb},{ub})"
621            );
622        }
623    }
624
625    #[test]
626    fn lp_lb_gt_ub_rejected() {
627        let cases: Vec<(f64, f64)> = vec![
628            (5.0, 1.0),
629            (1.0, 0.0),
630            (f64::INFINITY, f64::NEG_INFINITY),
631            (0.1, 0.0),
632        ];
633        for (lb, ub) in cases {
634            let res = make_lp(
635                vec![1.0, 2.0],
636                vec![5.0],
637                vec![1.0, 1.0],
638                vec![(lb, ub), (0.0, f64::INFINITY)],
639            );
640            assert!(
641                matches!(res, Err(SolverError::InvalidBounds { .. })),
642                "expected InvalidBounds for lb={lb} ub={ub}"
643            );
644        }
645    }
646
647    #[test]
648    fn lp_inf_bounds_accepted() {
649        let res = make_lp(
650            vec![1.0, 2.0],
651            vec![5.0],
652            vec![1.0, 1.0],
653            vec![(f64::NEG_INFINITY, f64::INFINITY), (0.0, f64::INFINITY)],
654        );
655        assert!(res.is_ok(), "±inf bounds should be valid");
656    }
657
658    // Bug C/D: (inf,inf) and (-inf,-inf) form empty intervals that lb>ub misses
659    // because inf==inf and -inf==-inf. These must be rejected.
660    #[test]
661    fn lp_empty_interval_same_inf_rejected() {
662        let cases: Vec<(f64, f64)> = vec![
663            (f64::INFINITY, f64::INFINITY),
664            (f64::NEG_INFINITY, f64::NEG_INFINITY),
665        ];
666        for (lb, ub) in cases {
667            let res = make_lp(
668                vec![1.0, 2.0],
669                vec![5.0],
670                vec![1.0, 1.0],
671                vec![(lb, ub), (0.0, f64::INFINITY)],
672            );
673            assert!(
674                matches!(res, Err(SolverError::InvalidBounds { index: 0, .. })),
675                "expected InvalidBounds for ({lb},{ub})"
676            );
677        }
678    }
679
680    // Sentinel: valid half-infinite and fixed-var bounds must be accepted.
681    // If bound validation is replaced by a no-op, fixed-var or half-inf cases
682    // would still be Ok (no rejection) — these are *positive* sentinels that
683    // confirm the validator doesn't over-reject.  The reject sentinels above
684    // ensure no-op → FAIL.
685    #[test]
686    fn lp_valid_bound_variants_accepted() {
687        let valid_cases: Vec<Vec<(f64, f64)>> = vec![
688            vec![(f64::NEG_INFINITY, 5.0), (0.0, f64::INFINITY)],
689            vec![(0.0, f64::INFINITY), (f64::NEG_INFINITY, f64::INFINITY)],
690            vec![(3.0, 3.0), (0.0, 10.0)],
691            vec![(0.0, 0.0), (0.0, 0.0)],
692        ];
693        for bounds in valid_cases {
694            let res = make_lp(vec![1.0, 2.0], vec![5.0], vec![1.0, 1.0], bounds.clone());
695            assert!(res.is_ok(), "expected Ok for bounds={bounds:?}");
696        }
697    }
698
699    // Sentinel: lb=+inf (alone) and ub=-inf (alone) must each be rejected.
700    #[test]
701    fn lp_lone_inf_bound_rejected() {
702        let cases: Vec<(f64, f64)> = vec![
703            (f64::INFINITY, 5.0),
704            (f64::INFINITY, f64::INFINITY),
705            (0.0, f64::NEG_INFINITY),
706            (f64::NEG_INFINITY, f64::NEG_INFINITY),
707        ];
708        for (lb, ub) in cases {
709            let res = make_lp(
710                vec![1.0, 2.0],
711                vec![5.0],
712                vec![1.0, 1.0],
713                vec![(lb, ub), (0.0, f64::INFINITY)],
714            );
715            assert!(
716                matches!(res, Err(SolverError::InvalidBounds { index: 0, .. })),
717                "expected InvalidBounds for ({lb},{ub})"
718            );
719        }
720    }
721}
otspot_core/problem/mod.rs

otspot_core/problem/
mod.rs
