Skip to main content

pounce_cli/
builtin.rs

1//! Built-in TNLP test problems for the CLI. Each problem is a
2//! self-contained `impl TNLP` so the CLI can run end-to-end without
3//! parsing an `.nl` file.
4//!
5//! Currently shipped:
6//!
7//! * `quadratic` — `min (x[0]-3)^2 + (x[1]-4)^2`, unconstrained,
8//!   exact Hessian = `2I`. Optimum at `(3, 4)`, `f* = 0`.
9//! * `rosenbrock` — `min 100*(x[1]-x[0]^2)^2 + (1-x[0])^2`,
10//!   unconstrained, exact Hessian. Optimum at `(1, 1)`, `f* = 0`.
11
12use pounce_common::types::{Index, Number};
13use pounce_nlp::tnlp::{
14    BoundsInfo, IndexStyle, IpoptCq, IpoptData, NlpInfo, Solution, SparsityRequest, StartingPoint,
15    TNLP,
16};
17use std::cell::RefCell;
18use std::rc::Rc;
19
20pub fn list() -> Vec<&'static str> {
21    vec![
22        "quadratic",
23        "rosenbrock",
24        "bounded-quadratic",
25        "eq-quadratic",
26        "circle",
27        "infeasible-eq",
28    ]
29}
30
31pub fn lookup(name: &str) -> Option<Rc<RefCell<dyn TNLP>>> {
32    match name {
33        "quadratic" => Some(Rc::new(RefCell::new(Quadratic::default()))),
34        "rosenbrock" => Some(Rc::new(RefCell::new(Rosenbrock::default()))),
35        "bounded-quadratic" => Some(Rc::new(RefCell::new(BoundedQuadratic::default()))),
36        "eq-quadratic" => Some(Rc::new(RefCell::new(EqQuadratic::default()))),
37        "circle" => Some(Rc::new(RefCell::new(Circle::default()))),
38        "infeasible-eq" => Some(Rc::new(RefCell::new(InfeasibleEq::default()))),
39        _ => None,
40    }
41}
42
43// --------------------------------------------------------------------
44// Quadratic: min (x0 - 3)^2 + (x1 - 4)^2
45// --------------------------------------------------------------------
46
47#[derive(Debug, Default)]
48pub struct Quadratic {
49    pub final_x: Option<[Number; 2]>,
50    pub final_obj: Number,
51}
52
53impl TNLP for Quadratic {
54    fn get_nlp_info(&mut self) -> Option<NlpInfo> {
55        Some(NlpInfo {
56            n: 2,
57            m: 0,
58            nnz_jac_g: 0,
59            nnz_h_lag: 2, // diagonal Hessian, lower triangle
60            index_style: IndexStyle::C,
61        })
62    }
63
64    fn get_bounds_info(&mut self, b: BoundsInfo<'_>) -> bool {
65        b.x_l.iter_mut().for_each(|v| *v = -2e19);
66        b.x_u.iter_mut().for_each(|v| *v = 2e19);
67        true
68    }
69
70    fn get_starting_point(&mut self, sp: StartingPoint<'_>) -> bool {
71        sp.x.copy_from_slice(&[0.0, 0.0]);
72        true
73    }
74
75    fn eval_f(&mut self, x: &[Number], _new_x: bool) -> Option<Number> {
76        Some((x[0] - 3.0).powi(2) + (x[1] - 4.0).powi(2))
77    }
78
79    fn eval_grad_f(&mut self, x: &[Number], _new_x: bool, grad: &mut [Number]) -> bool {
80        grad[0] = 2.0 * (x[0] - 3.0);
81        grad[1] = 2.0 * (x[1] - 4.0);
82        true
83    }
84
85    fn eval_g(&mut self, _x: &[Number], _new_x: bool, _g: &mut [Number]) -> bool {
86        true
87    }
88
89    fn eval_jac_g(
90        &mut self,
91        _x: Option<&[Number]>,
92        _new_x: bool,
93        _mode: SparsityRequest<'_>,
94    ) -> bool {
95        true
96    }
97
98    fn eval_h(
99        &mut self,
100        _x: Option<&[Number]>,
101        _new_x: bool,
102        obj_factor: Number,
103        _lambda: Option<&[Number]>,
104        _new_lambda: bool,
105        mode: SparsityRequest<'_>,
106    ) -> bool {
107        match mode {
108            SparsityRequest::Structure { irow, jcol } => {
109                irow.copy_from_slice(&[0, 1]);
110                jcol.copy_from_slice(&[0, 1]);
111            }
112            SparsityRequest::Values { values } => {
113                values[0] = 2.0 * obj_factor;
114                values[1] = 2.0 * obj_factor;
115            }
116        }
117        true
118    }
119
120    fn finalize_solution(&mut self, sol: Solution<'_>, _d: &IpoptData, _q: &IpoptCq) {
121        self.final_x = Some([sol.x[0], sol.x[1]]);
122        self.final_obj = sol.obj_value;
123    }
124}
125
126// --------------------------------------------------------------------
127// Rosenbrock: min 100 (x1 - x0^2)^2 + (1 - x0)^2
128// --------------------------------------------------------------------
129
130#[derive(Debug, Default)]
131pub struct Rosenbrock {
132    pub final_x: Option<[Number; 2]>,
133    pub final_obj: Number,
134}
135
136impl TNLP for Rosenbrock {
137    fn get_nlp_info(&mut self) -> Option<NlpInfo> {
138        Some(NlpInfo {
139            n: 2,
140            m: 0,
141            nnz_jac_g: 0,
142            nnz_h_lag: 3, // dense 2x2 lower triangle: (0,0), (1,0), (1,1)
143            index_style: IndexStyle::C,
144        })
145    }
146
147    fn get_bounds_info(&mut self, b: BoundsInfo<'_>) -> bool {
148        b.x_l.iter_mut().for_each(|v| *v = -2e19);
149        b.x_u.iter_mut().for_each(|v| *v = 2e19);
150        true
151    }
152
153    fn get_starting_point(&mut self, sp: StartingPoint<'_>) -> bool {
154        sp.x.copy_from_slice(&[-1.2, 1.0]);
155        true
156    }
157
158    fn eval_f(&mut self, x: &[Number], _new_x: bool) -> Option<Number> {
159        let a = x[1] - x[0] * x[0];
160        let b = 1.0 - x[0];
161        Some(100.0 * a * a + b * b)
162    }
163
164    fn eval_grad_f(&mut self, x: &[Number], _new_x: bool, grad: &mut [Number]) -> bool {
165        // d/dx0 = -400 x0 (x1 - x0^2) - 2 (1 - x0)
166        // d/dx1 =  200 (x1 - x0^2)
167        grad[0] = -400.0 * x[0] * (x[1] - x[0] * x[0]) - 2.0 * (1.0 - x[0]);
168        grad[1] = 200.0 * (x[1] - x[0] * x[0]);
169        true
170    }
171
172    fn eval_g(&mut self, _x: &[Number], _new_x: bool, _g: &mut [Number]) -> bool {
173        true
174    }
175
176    fn eval_jac_g(
177        &mut self,
178        _x: Option<&[Number]>,
179        _new_x: bool,
180        _mode: SparsityRequest<'_>,
181    ) -> bool {
182        true
183    }
184
185    fn eval_h(
186        &mut self,
187        x: Option<&[Number]>,
188        _new_x: bool,
189        obj_factor: Number,
190        _lambda: Option<&[Number]>,
191        _new_lambda: bool,
192        mode: SparsityRequest<'_>,
193    ) -> bool {
194        match mode {
195            SparsityRequest::Structure { irow, jcol } => {
196                // Lower triangle (row >= col): (0,0), (1,0), (1,1).
197                irow.copy_from_slice(&[0, 1, 1]);
198                jcol.copy_from_slice(&[0, 0, 1]);
199            }
200            SparsityRequest::Values { values } => {
201                let x = x.unwrap_or(&[0.0, 0.0]);
202                let h00 = -400.0 * (x[1] - 3.0 * x[0] * x[0]) + 2.0;
203                let h10 = -400.0 * x[0];
204                let h11 = 200.0;
205                values[0] = obj_factor * h00;
206                values[1] = obj_factor * h10;
207                values[2] = obj_factor * h11;
208            }
209        }
210        true
211    }
212
213    fn finalize_solution(&mut self, sol: Solution<'_>, _d: &IpoptData, _q: &IpoptCq) {
214        self.final_x = Some([sol.x[0], sol.x[1]]);
215        self.final_obj = sol.obj_value;
216    }
217}
218
219// Small helper for doc references
220#[allow(dead_code)]
221fn _ix(_: Index) {}
222
223// --------------------------------------------------------------------
224// BoundedQuadratic: min (x0-3)^2 + (x1-4)^2 s.t. 0 <= x0 <= 2, 0 <= x1 <= 2
225// Optimum is at the corner (2, 2) with f* = 1 + 4 = 5.
226// --------------------------------------------------------------------
227
228#[derive(Debug, Default)]
229pub struct BoundedQuadratic {
230    pub final_x: Option<[Number; 2]>,
231    pub final_obj: Number,
232}
233
234impl TNLP for BoundedQuadratic {
235    fn get_nlp_info(&mut self) -> Option<NlpInfo> {
236        Some(NlpInfo {
237            n: 2,
238            m: 0,
239            nnz_jac_g: 0,
240            nnz_h_lag: 2,
241            index_style: IndexStyle::C,
242        })
243    }
244
245    fn get_bounds_info(&mut self, b: BoundsInfo<'_>) -> bool {
246        b.x_l.copy_from_slice(&[0.0, 0.0]);
247        b.x_u.copy_from_slice(&[2.0, 2.0]);
248        true
249    }
250
251    fn get_starting_point(&mut self, sp: StartingPoint<'_>) -> bool {
252        sp.x.copy_from_slice(&[1.0, 1.0]);
253        true
254    }
255
256    fn eval_f(&mut self, x: &[Number], _new_x: bool) -> Option<Number> {
257        Some((x[0] - 3.0).powi(2) + (x[1] - 4.0).powi(2))
258    }
259
260    fn eval_grad_f(&mut self, x: &[Number], _new_x: bool, grad: &mut [Number]) -> bool {
261        grad[0] = 2.0 * (x[0] - 3.0);
262        grad[1] = 2.0 * (x[1] - 4.0);
263        true
264    }
265
266    fn eval_g(&mut self, _x: &[Number], _new_x: bool, _g: &mut [Number]) -> bool {
267        true
268    }
269
270    fn eval_jac_g(
271        &mut self,
272        _x: Option<&[Number]>,
273        _new_x: bool,
274        _mode: SparsityRequest<'_>,
275    ) -> bool {
276        true
277    }
278
279    fn eval_h(
280        &mut self,
281        _x: Option<&[Number]>,
282        _new_x: bool,
283        obj_factor: Number,
284        _lambda: Option<&[Number]>,
285        _new_lambda: bool,
286        mode: SparsityRequest<'_>,
287    ) -> bool {
288        match mode {
289            SparsityRequest::Structure { irow, jcol } => {
290                irow.copy_from_slice(&[0, 1]);
291                jcol.copy_from_slice(&[0, 1]);
292            }
293            SparsityRequest::Values { values } => {
294                values[0] = 2.0 * obj_factor;
295                values[1] = 2.0 * obj_factor;
296            }
297        }
298        true
299    }
300
301    fn finalize_solution(&mut self, sol: Solution<'_>, _d: &IpoptData, _q: &IpoptCq) {
302        self.final_x = Some([sol.x[0], sol.x[1]]);
303        self.final_obj = sol.obj_value;
304    }
305}
306
307// --------------------------------------------------------------------
308// EqQuadratic: min x0^2 + x1^2  s.t.  x0 + x1 = 1
309// Optimum at (1/2, 1/2), f* = 1/2, multiplier y = -1.
310// --------------------------------------------------------------------
311
312#[derive(Debug, Default)]
313pub struct EqQuadratic {
314    pub final_x: Option<[Number; 2]>,
315    pub final_obj: Number,
316}
317
318impl TNLP for EqQuadratic {
319    fn get_nlp_info(&mut self) -> Option<NlpInfo> {
320        Some(NlpInfo {
321            n: 2,
322            m: 1,
323            nnz_jac_g: 2,
324            nnz_h_lag: 2,
325            index_style: IndexStyle::C,
326        })
327    }
328
329    fn get_bounds_info(&mut self, b: BoundsInfo<'_>) -> bool {
330        b.x_l.iter_mut().for_each(|v| *v = -2e19);
331        b.x_u.iter_mut().for_each(|v| *v = 2e19);
332        b.g_l[0] = 1.0;
333        b.g_u[0] = 1.0;
334        true
335    }
336
337    fn get_starting_point(&mut self, sp: StartingPoint<'_>) -> bool {
338        sp.x.copy_from_slice(&[0.0, 0.0]);
339        true
340    }
341
342    fn eval_f(&mut self, x: &[Number], _new_x: bool) -> Option<Number> {
343        Some(x[0] * x[0] + x[1] * x[1])
344    }
345
346    fn eval_grad_f(&mut self, x: &[Number], _new_x: bool, grad: &mut [Number]) -> bool {
347        grad[0] = 2.0 * x[0];
348        grad[1] = 2.0 * x[1];
349        true
350    }
351
352    fn eval_g(&mut self, x: &[Number], _new_x: bool, g: &mut [Number]) -> bool {
353        g[0] = x[0] + x[1];
354        true
355    }
356
357    fn eval_jac_g(
358        &mut self,
359        _x: Option<&[Number]>,
360        _new_x: bool,
361        mode: SparsityRequest<'_>,
362    ) -> bool {
363        match mode {
364            SparsityRequest::Structure { irow, jcol } => {
365                irow.copy_from_slice(&[0, 0]);
366                jcol.copy_from_slice(&[0, 1]);
367            }
368            SparsityRequest::Values { values } => {
369                values[0] = 1.0;
370                values[1] = 1.0;
371            }
372        }
373        true
374    }
375
376    fn eval_h(
377        &mut self,
378        _x: Option<&[Number]>,
379        _new_x: bool,
380        obj_factor: Number,
381        _lambda: Option<&[Number]>,
382        _new_lambda: bool,
383        mode: SparsityRequest<'_>,
384    ) -> bool {
385        match mode {
386            SparsityRequest::Structure { irow, jcol } => {
387                irow.copy_from_slice(&[0, 1]);
388                jcol.copy_from_slice(&[0, 1]);
389            }
390            SparsityRequest::Values { values } => {
391                values[0] = 2.0 * obj_factor;
392                values[1] = 2.0 * obj_factor;
393            }
394        }
395        true
396    }
397
398    fn finalize_solution(&mut self, sol: Solution<'_>, _d: &IpoptData, _q: &IpoptCq) {
399        self.final_x = Some([sol.x[0], sol.x[1]]);
400        self.final_obj = sol.obj_value;
401    }
402}
403
404// --------------------------------------------------------------------
405// Circle: min  x0  s.t.  x0^2 + x1^2 = 1
406// Optimum at (-1, 0), f* = -1, multiplier y = 1/2.
407// Tests nonlinear equality constraint with non-trivial Hessian
408// contribution from the constraint (∇²g_0 = 2I) into the Lagrangian.
409// --------------------------------------------------------------------
410
411#[derive(Debug, Default)]
412pub struct Circle {
413    pub final_x: Option<[Number; 2]>,
414    pub final_obj: Number,
415}
416
417impl TNLP for Circle {
418    fn get_nlp_info(&mut self) -> Option<NlpInfo> {
419        Some(NlpInfo {
420            n: 2,
421            m: 1,
422            nnz_jac_g: 2,
423            nnz_h_lag: 2,
424            index_style: IndexStyle::C,
425        })
426    }
427
428    fn get_bounds_info(&mut self, b: BoundsInfo<'_>) -> bool {
429        b.x_l.iter_mut().for_each(|v| *v = -2e19);
430        b.x_u.iter_mut().for_each(|v| *v = 2e19);
431        b.g_l[0] = 1.0;
432        b.g_u[0] = 1.0;
433        true
434    }
435
436    fn get_starting_point(&mut self, sp: StartingPoint<'_>) -> bool {
437        sp.x.copy_from_slice(&[-0.5, 0.5]);
438        true
439    }
440
441    fn eval_f(&mut self, x: &[Number], _new_x: bool) -> Option<Number> {
442        Some(x[0])
443    }
444
445    fn eval_grad_f(&mut self, _x: &[Number], _new_x: bool, grad: &mut [Number]) -> bool {
446        grad[0] = 1.0;
447        grad[1] = 0.0;
448        true
449    }
450
451    fn eval_g(&mut self, x: &[Number], _new_x: bool, g: &mut [Number]) -> bool {
452        g[0] = x[0] * x[0] + x[1] * x[1];
453        true
454    }
455
456    fn eval_jac_g(
457        &mut self,
458        x: Option<&[Number]>,
459        _new_x: bool,
460        mode: SparsityRequest<'_>,
461    ) -> bool {
462        match mode {
463            SparsityRequest::Structure { irow, jcol } => {
464                irow.copy_from_slice(&[0, 0]);
465                jcol.copy_from_slice(&[0, 1]);
466            }
467            SparsityRequest::Values { values } => {
468                let x = x.unwrap_or(&[0.0, 0.0]);
469                values[0] = 2.0 * x[0];
470                values[1] = 2.0 * x[1];
471            }
472        }
473        true
474    }
475
476    fn eval_h(
477        &mut self,
478        _x: Option<&[Number]>,
479        _new_x: bool,
480        _obj_factor: Number,
481        lambda: Option<&[Number]>,
482        _new_lambda: bool,
483        mode: SparsityRequest<'_>,
484    ) -> bool {
485        match mode {
486            SparsityRequest::Structure { irow, jcol } => {
487                irow.copy_from_slice(&[0, 1]);
488                jcol.copy_from_slice(&[0, 1]);
489            }
490            SparsityRequest::Values { values } => {
491                // ∇²L = obj_factor * ∇²f + λ * ∇²g_0 = 0 + λ * 2I
492                let lam = lambda.map(|l| l[0]).unwrap_or(0.0);
493                values[0] = 2.0 * lam;
494                values[1] = 2.0 * lam;
495            }
496        }
497        true
498    }
499
500    fn finalize_solution(&mut self, sol: Solution<'_>, _d: &IpoptData, _q: &IpoptCq) {
501        self.final_x = Some([sol.x[0], sol.x[1]]);
502        self.final_obj = sol.obj_value;
503    }
504}
505
506// --------------------------------------------------------------------
507// InfeasibleEq: min x0^2 + x1^2
508//   s.t.  x0 + x1 = 1   (g_0)
509//         x0 + x1 = 2   (g_1)
510// The two equalities are mutually contradictory, so no feasible point
511// exists. The standard solve drives the restoration phase, which also
512// cannot achieve feasibility, returning Restoration_Failed. With
513// `l1_fallback_on_restoration_failure=yes` (or
514// `l1_exact_penalty_barrier=yes`), the CLI then performs a second
515// inner solve via the ℓ₁-exact penalty-barrier wrapper. That second
516// pass is what exercises the multi-pass restoration factory provider
517// path — the very path that previously panicked with
518// "restoration factory invoked more than once".
519// --------------------------------------------------------------------
520
521#[derive(Debug, Default)]
522pub struct InfeasibleEq {
523    pub final_x: Option<[Number; 2]>,
524    pub final_obj: Number,
525}
526
527impl TNLP for InfeasibleEq {
528    fn get_nlp_info(&mut self) -> Option<NlpInfo> {
529        Some(NlpInfo {
530            n: 2,
531            m: 2,
532            nnz_jac_g: 4,
533            nnz_h_lag: 2,
534            index_style: IndexStyle::C,
535        })
536    }
537
538    fn get_bounds_info(&mut self, b: BoundsInfo<'_>) -> bool {
539        b.x_l.iter_mut().for_each(|v| *v = -2e19);
540        b.x_u.iter_mut().for_each(|v| *v = 2e19);
541        b.g_l[0] = 1.0;
542        b.g_u[0] = 1.0;
543        b.g_l[1] = 2.0;
544        b.g_u[1] = 2.0;
545        true
546    }
547
548    fn get_starting_point(&mut self, sp: StartingPoint<'_>) -> bool {
549        sp.x.copy_from_slice(&[0.0, 0.0]);
550        true
551    }
552
553    fn eval_f(&mut self, x: &[Number], _new_x: bool) -> Option<Number> {
554        Some(x[0] * x[0] + x[1] * x[1])
555    }
556
557    fn eval_grad_f(&mut self, x: &[Number], _new_x: bool, grad: &mut [Number]) -> bool {
558        grad[0] = 2.0 * x[0];
559        grad[1] = 2.0 * x[1];
560        true
561    }
562
563    fn eval_g(&mut self, x: &[Number], _new_x: bool, g: &mut [Number]) -> bool {
564        g[0] = x[0] + x[1];
565        g[1] = x[0] + x[1];
566        true
567    }
568
569    fn eval_jac_g(
570        &mut self,
571        _x: Option<&[Number]>,
572        _new_x: bool,
573        mode: SparsityRequest<'_>,
574    ) -> bool {
575        match mode {
576            SparsityRequest::Structure { irow, jcol } => {
577                irow.copy_from_slice(&[0, 0, 1, 1]);
578                jcol.copy_from_slice(&[0, 1, 0, 1]);
579            }
580            SparsityRequest::Values { values } => {
581                values.copy_from_slice(&[1.0, 1.0, 1.0, 1.0]);
582            }
583        }
584        true
585    }
586
587    fn eval_h(
588        &mut self,
589        _x: Option<&[Number]>,
590        _new_x: bool,
591        obj_factor: Number,
592        _lambda: Option<&[Number]>,
593        _new_lambda: bool,
594        mode: SparsityRequest<'_>,
595    ) -> bool {
596        match mode {
597            SparsityRequest::Structure { irow, jcol } => {
598                irow.copy_from_slice(&[0, 1]);
599                jcol.copy_from_slice(&[0, 1]);
600            }
601            SparsityRequest::Values { values } => {
602                values[0] = 2.0 * obj_factor;
603                values[1] = 2.0 * obj_factor;
604            }
605        }
606        true
607    }
608
609    fn finalize_solution(&mut self, sol: Solution<'_>, _d: &IpoptData, _q: &IpoptCq) {
610        self.final_x = Some([sol.x[0], sol.x[1]]);
611        self.final_obj = sol.obj_value;
612    }
613}
614
615#[cfg(test)]
616mod tests {
617    use super::*;
618
619    #[test]
620    fn list_contains_known_problems() {
621        let l = list();
622        assert!(l.contains(&"quadratic"));
623        assert!(l.contains(&"rosenbrock"));
624    }
625
626    #[test]
627    fn quadratic_evaluates_correctly() {
628        let mut q = Quadratic::default();
629        let f = q.eval_f(&[3.0, 4.0], true).unwrap();
630        assert_eq!(f, 0.0);
631        let mut g = [0.0; 2];
632        q.eval_grad_f(&[0.0, 0.0], true, &mut g);
633        assert_eq!(g, [-6.0, -8.0]);
634    }
635
636    #[test]
637    fn rosenbrock_grad_zero_at_optimum() {
638        let mut r = Rosenbrock::default();
639        let f = r.eval_f(&[1.0, 1.0], true).unwrap();
640        assert!(f.abs() < 1e-15);
641        let mut g = [0.0; 2];
642        r.eval_grad_f(&[1.0, 1.0], true, &mut g);
643        assert!(g[0].abs() < 1e-12);
644        assert!(g[1].abs() < 1e-12);
645    }
646
647    #[test]
648    fn lookup_returns_known_and_rejects_unknown() {
649        assert!(lookup("quadratic").is_some());
650        assert!(lookup("rosenbrock").is_some());
651        assert!(lookup("nonsense").is_none());
652    }
653}