stats_claw/optimizers/
objectives.rs

1//! Standard test objectives for the optimizer suite.
2//!
3//! [`Quadratic`] is a strictly convex bowl with a known global minimizer, used
4//! to assert convergence and `scipy.optimize` agreement on an easy problem.
5//! [`Rosenbrock`] is the classic non-convex valley with its global minimum at
6//! `(1, 1)`, used to stress the gradient-based, quasi-Newton, and derivative-free
7//! methods.
8//!
9//! # Examples
10//!
11//! ```
12//! use stats_claw::optimizers::objectives::Quadratic;
13//! use stats_claw::optimizers::Objective;
14//!
15//! let q = Quadratic::new(vec![3.0, -2.0]);
16//! // The minimum value is 0, attained at the center.
17//! assert!(q.value(&[3.0, -2.0]).abs() < 1e-15, "f(c) = {}", q.value(&[3.0, -2.0]));
18//! // The gradient is zero at the center.
19//! use stats_claw::optimizers::norm;
20//! assert!(norm(&q.grad(&[3.0, -2.0])) < 1e-15);
21//! ```
22
23use super::Objective;
24
25/// Separable quadratic `f(x) = Σ (xᵢ − cᵢ)²`.
26///
27/// Strictly convex with a single global minimum of `0` at the center `c`. Its
28/// gradient is `2(x − c)` and its Hessian is `2·I`.
29#[derive(Debug, Clone)]
30pub struct Quadratic {
31    /// The minimizer `c`; `f` attains its minimum of `0` here.
32    center: Vec<f64>,
33}
34
35impl Quadratic {
36    /// Builds a quadratic with global minimum `0` at `center`.
37    ///
38    /// # Arguments
39    ///
40    /// * `center` — the minimizer `c`; its length is the problem dimension.
41    ///
42    /// # Returns
43    ///
44    /// A [`Quadratic`] whose [`Objective`] impl minimizes to `center`.
45    #[must_use]
46    pub const fn new(center: Vec<f64>) -> Self {
47        Self { center }
48    }
49}
50
51impl Objective for Quadratic {
52    fn value(&self, x: &[f64]) -> f64 {
53        x.iter()
54            .zip(&self.center)
55            .map(|(a, c)| (a - c) * (a - c))
56            .sum()
57    }
58
59    fn grad(&self, x: &[f64]) -> Vec<f64> {
60        x.iter()
61            .zip(&self.center)
62            .map(|(a, c)| 2.0 * (a - c))
63            .collect()
64    }
65
66    fn hessian(&self, x: &[f64]) -> Vec<Vec<f64>> {
67        let n = x.len();
68        (0..n)
69            .map(|i| (0..n).map(|j| if i == j { 2.0 } else { 0.0 }).collect())
70            .collect()
71    }
72}
73
74/// The two-dimensional Rosenbrock function
75/// `f(x) = (1 − x₀)² + 100·(x₁ − x₀²)²`.
76///
77/// Non-convex, with a narrow curved valley and a global minimum of `0` at
78/// `(1, 1)`. A demanding benchmark for first-order methods.
79#[derive(Debug, Clone, Copy)]
80pub struct Rosenbrock;
81
82impl Objective for Rosenbrock {
83    fn value(&self, x: &[f64]) -> f64 {
84        let x0 = *x.first().unwrap_or(&0.0);
85        let x1 = *x.get(1).unwrap_or(&0.0);
86        let a = 1.0 - x0;
87        let b = x1 - x0 * x0;
88        b.mul_add(100.0 * b, a * a)
89    }
90
91    fn grad(&self, x: &[f64]) -> Vec<f64> {
92        let x0 = *x.first().unwrap_or(&0.0);
93        let x1 = *x.get(1).unwrap_or(&0.0);
94        let d0 = (-400.0 * x0).mul_add(x1 - x0 * x0, -2.0 * (1.0 - x0));
95        let d1 = 200.0 * (x1 - x0 * x0);
96        vec![d0, d1]
97    }
98
99    fn hessian(&self, x: &[f64]) -> Vec<Vec<f64>> {
100        let x0 = *x.first().unwrap_or(&0.0);
101        let x1 = *x.get(1).unwrap_or(&0.0);
102        let h00 = (-400.0_f64).mul_add(x1, (1200.0 * x0).mul_add(x0, 2.0));
103        let h01 = -400.0 * x0;
104        vec![vec![h00, h01], vec![h01, 200.0]]
105    }
106}
107
108#[cfg(test)]
109mod tests {
110    use super::*;
111
112    #[test]
113    fn quadratic_minimum_is_zero_at_center() {
114        let q = Quadratic::new(vec![3.0, -2.0]);
115        assert!(
116            q.value(&[3.0, -2.0]).abs() < 1e-15,
117            "f(c) = {}",
118            q.value(&[3.0, -2.0])
119        );
120        let g = q.grad(&[3.0, -2.0]);
121        assert!(super::super::norm(&g) < 1e-15, "grad at center nonzero");
122    }
123
124    #[test]
125    fn rosenbrock_minimum_is_zero_at_one_one() {
126        let r = Rosenbrock;
127        assert!(
128            r.value(&[1.0, 1.0]).abs() < 1e-15,
129            "f(1,1) = {}",
130            r.value(&[1.0, 1.0])
131        );
132        let g = r.grad(&[1.0, 1.0]);
133        assert!(super::super::norm(&g) < 1e-12, "grad at (1,1) nonzero");
134    }
135
136    #[test]
137    fn rosenbrock_analytic_hessian_matches_finite_difference() {
138        // Override agrees with the trait default (finite-difference) within 1e-4.
139        let r = Rosenbrock;
140        let x = [0.5, 0.3];
141        let analytic = r.hessian(&x);
142        let fd = Objective::hessian(&FdRosen, &x);
143        for (ra, rf) in analytic.iter().zip(&fd) {
144            for (a, f) in ra.iter().zip(rf) {
145                assert!((a - f).abs() < 1e-4, "hessian mismatch: {a} vs {f}");
146            }
147        }
148    }
149
150    /// Rosenbrock without the Hessian override, to exercise the trait default.
151    struct FdRosen;
152    impl Objective for FdRosen {
153        fn value(&self, x: &[f64]) -> f64 {
154            Rosenbrock.value(x)
155        }
156        fn grad(&self, x: &[f64]) -> Vec<f64> {
157            Rosenbrock.grad(x)
158        }
159    }
160}
stats_claw/optimizers/objectives.rs

stats_claw/optimizers/
objectives.rs
