spectrafit_models/
polynomial.rs

1use crate::Model;
2
3/// Constant model: `f(x) = c`
4///
5/// Parameters (in order): `[c]`
6pub struct Constant;
7
8impl Model for Constant {
9    fn eval(&self, _x: &[f64], params: &[f64]) -> f64 {
10        params[0]
11    }
12
13    fn jacobian(&self, _x: &[f64], _params: &[f64]) -> Vec<f64> {
14        vec![1.0]
15    }
16
17    #[inline]
18    fn jacobian_into(&self, _x: &[f64], _params: &[f64], out: &mut [f64]) {
19        out[0] = 1.0;
20    }
21
22    fn param_names(&self) -> Vec<std::borrow::Cow<'static, str>> {
23        vec!["c".into()]
24    }
25}
26
27/// Linear model: `f(x) = slope * x₀ + intercept`
28///
29/// Parameters (in order): `[slope, intercept]`
30pub struct Linear;
31
32impl Model for Linear {
33    fn eval(&self, x: &[f64], params: &[f64]) -> f64 {
34        params[0] * x[0] + params[1]
35    }
36
37    fn jacobian(&self, x: &[f64], _params: &[f64]) -> Vec<f64> {
38        vec![x[0], 1.0]
39    }
40
41    #[inline]
42    fn jacobian_into(&self, x: &[f64], _params: &[f64], out: &mut [f64]) {
43        out[0] = x[0];
44        out[1] = 1.0;
45    }
46
47    fn param_names(&self) -> Vec<std::borrow::Cow<'static, str>> {
48        vec!["slope".into(), "intercept".into()]
49    }
50}
51
52/// Quadratic bowl model: `f(x) = amplitude · (x₀ − center)² + offset`
53///
54/// A convex parabola used for the `convex_baseline` benchmark family (clean
55/// quadratic objectives). Summing several of these nodes builds a sum-of-squares
56/// landscape; pairing one with a `Linear` node gives a tilted bowl.
57///
58/// Parameters (in order): `[amplitude, center, offset]`
59pub struct Quadratic;
60
61impl Model for Quadratic {
62    fn eval(&self, x: &[f64], params: &[f64]) -> f64 {
63        let d = x[0] - params[1];
64        params[0] * d * d + params[2]
65    }
66
67    fn jacobian(&self, x: &[f64], params: &[f64]) -> Vec<f64> {
68        let d = x[0] - params[1];
69        // ∂/∂amplitude = d²; ∂/∂center = −2·A·d; ∂/∂offset = 1
70        vec![d * d, -2.0 * params[0] * d, 1.0]
71    }
72
73    #[inline]
74    fn jacobian_into(&self, x: &[f64], params: &[f64], out: &mut [f64]) {
75        let d = x[0] - params[1];
76        out[0] = d * d;
77        out[1] = -2.0 * params[0] * d;
78        out[2] = 1.0;
79    }
80
81    fn param_names(&self) -> Vec<std::borrow::Cow<'static, str>> {
82        vec!["amplitude".into(), "center".into(), "offset".into()]
83    }
84}
85
86#[cfg(test)]
87mod tests {
88    use super::*;
89    use approx::assert_relative_eq;
90
91    #[test]
92    fn constant_eval() {
93        let c = Constant;
94        assert_relative_eq!(c.eval(&[99.0], &[5.0]), 5.0, epsilon = 1e-12);
95    }
96
97    #[test]
98    fn constant_jacobian() {
99        let c = Constant;
100        let j = c.jacobian(&[1.0], &[5.0]);
101        assert_eq!(j.len(), 1);
102        assert_relative_eq!(j[0], 1.0, epsilon = 1e-12);
103    }
104
105    #[test]
106    fn linear_eval() {
107        let l = Linear;
108        // slope=3, intercept=1, x=2 → 3*2+1 = 7
109        assert_relative_eq!(l.eval(&[2.0], &[3.0, 1.0]), 7.0, epsilon = 1e-12);
110    }
111
112    #[test]
113    fn linear_jacobian_shape() {
114        let l = Linear;
115        let j = l.jacobian(&[2.0], &[3.0, 1.0]);
116        assert_eq!(j.len(), 2);
117    }
118
119    #[test]
120    fn linear_jacobian_values() {
121        let l = Linear;
122        // jac[0] = x₀ = 2.0; jac[1] = 1.0
123        let j = l.jacobian(&[2.0], &[3.0, 1.0]);
124        assert_relative_eq!(j[0], 2.0, epsilon = 1e-12);
125        assert_relative_eq!(j[1], 1.0, epsilon = 1e-12);
126    }
127
128    #[test]
129    fn linear_eval_at_zero() {
130        let l = Linear;
131        // x=0 → intercept only
132        assert_relative_eq!(l.eval(&[0.0], &[5.0, 3.0]), 3.0, epsilon = 1e-12);
133    }
134
135    #[test]
136    fn quadratic_eval() {
137        let q = Quadratic;
138        // A=2, c=0.5, b=0.3, x=2 → 2*(1.5)^2 + 0.3 = 4.8
139        assert_relative_eq!(q.eval(&[2.0], &[2.0, 0.5, 0.3]), 4.8, epsilon = 1e-12);
140        // at the vertex x=c the bowl equals the offset
141        assert_relative_eq!(q.eval(&[0.5], &[2.0, 0.5, 0.3]), 0.3, epsilon = 1e-12);
142    }
143
144    #[test]
145    fn quadratic_jacobian_matches_numerical() {
146        // Already a central difference; widened from one regime/x-point to
147        // three regimes and several x-points relative to the vertex, and
148        // tightened epsilon from 1e-6 to 1e-7.
149        let q = Quadratic;
150        let param_sets = [
151            [1.8, 0.4, -0.2], // nominal
152            [1e-3, 0.0, 0.0], // small amplitude, vertex at origin
153            [5.0, -2.0, 3.0], // offset vertex, large amplitude
154        ];
155        for params in param_sets {
156            let center = params[1];
157            for &mult in &[-3.0_f64, -1.0, -0.1, 0.0, 0.1, 1.0, 3.0] {
158                let x = [center + mult];
159                let analytic = q.jacobian(&x, &params);
160                for i in 0..params.len() {
161                    let h = 1e-6 * params[i].abs().max(1.0);
162                    let mut pp = params;
163                    let mut pm = params;
164                    pp[i] += h;
165                    pm[i] -= h;
166                    let num = (q.eval(&x, &pp) - q.eval(&x, &pm)) / (2.0 * h);
167                    assert_relative_eq!(analytic[i], num, epsilon = 1e-7, max_relative = 1e-6);
168                }
169            }
170        }
171    }
172}
spectrafit_models/polynomial.rs

spectrafit_models/
polynomial.rs
