scivex_optim/minimize/
mod.rs

1//! Multi-dimensional unconstrained optimization.
2
3mod bfgs;
4mod gradient_descent;
5mod lbfgsb;
6mod nelder_mead;
7
8pub use bfgs::bfgs;
9pub use gradient_descent::gradient_descent;
10pub use lbfgsb::{Bounds, lbfgsb};
11pub use nelder_mead::nelder_mead;
12
13use scivex_core::{Float, Tensor};
14
15/// Result of a multi-dimensional minimization algorithm.
16///
17/// # Examples
18///
19/// ```
20/// # use scivex_core::Tensor;
21/// # use scivex_optim::minimize::{bfgs, MinimizeOptions};
22/// let f = |x: &Tensor<f64>| { let s = x.as_slice(); s[0] * s[0] + s[1] * s[1] };
23/// let g = |x: &Tensor<f64>| {
24///     let s = x.as_slice();
25///     Tensor::from_vec(vec![2.0 * s[0], 2.0 * s[1]], vec![2]).unwrap()
26/// };
27/// let x0 = Tensor::from_vec(vec![3.0_f64, 4.0], vec![2]).unwrap();
28/// let result = bfgs(f, g, &x0, &MinimizeOptions::default()).unwrap();
29/// assert!(result.converged);
30/// ```
31#[cfg_attr(
32    feature = "serde-support",
33    derive(serde::Serialize, serde::Deserialize)
34)]
35#[derive(Debug, Clone)]
36pub struct MinimizeResult<T: Float> {
37    /// The estimated minimizer.
38    pub x: Tensor<T>,
39    /// The function value at the minimizer.
40    pub f_val: T,
41    /// The gradient at the minimizer (if available).
42    pub grad: Option<Tensor<T>>,
43    /// Number of iterations performed.
44    pub iterations: usize,
45    /// Number of function evaluations.
46    pub f_evals: usize,
47    /// Number of gradient evaluations.
48    pub g_evals: usize,
49    /// Whether the algorithm converged.
50    pub converged: bool,
51}
52
53/// Options controlling multi-dimensional minimization.
54///
55/// # Examples
56///
57/// ```
58/// # use scivex_optim::minimize::MinimizeOptions;
59/// let opts = MinimizeOptions::<f64>::default();
60/// assert_eq!(opts.max_iter, 1000);
61/// ```
62#[cfg_attr(
63    feature = "serde-support",
64    derive(serde::Serialize, serde::Deserialize)
65)]
66#[derive(Debug, Clone)]
67pub struct MinimizeOptions<T: Float> {
68    /// Gradient norm tolerance for convergence.
69    pub gtol: T,
70    /// Step-size tolerance.
71    pub xtol: T,
72    /// Function value change tolerance.
73    pub ftol: T,
74    /// Maximum number of iterations.
75    pub max_iter: usize,
76    /// Learning rate (for gradient descent).
77    pub learning_rate: T,
78}
79
80impl<T: Float> Default for MinimizeOptions<T> {
81    fn default() -> Self {
82        Self {
83            gtol: T::from_f64(1e-8),
84            xtol: T::from_f64(1e-12),
85            ftol: T::from_f64(1e-12),
86            max_iter: 1000,
87            learning_rate: T::from_f64(0.01),
88        }
89    }
90}
91
92/// Compute the gradient of `f` at `x` using central finite differences.
93///
94/// `grad_i = (f(x + h*e_i) - f(x - h*e_i)) / (2*h)` where `h` is a
95/// small perturbation and `e_i` is the i-th unit vector.
96///
97/// # Examples
98///
99/// ```
100/// # use scivex_core::Tensor;
101/// # use scivex_optim::minimize::numerical_gradient;
102/// let f = |x: &Tensor<f64>| x.as_slice()[0].powi(2); // f(x) = x²
103/// let x = Tensor::from_vec(vec![3.0_f64], vec![1]).unwrap();
104/// let g = numerical_gradient(&f, &x);
105/// assert!((g.as_slice()[0] - 6.0).abs() < 1e-4); // f'(3) = 6
106/// ```
107pub fn numerical_gradient<T: Float, F: Fn(&Tensor<T>) -> T>(f: &F, x: &Tensor<T>) -> Tensor<T> {
108    let n = x.numel();
109    let h = T::from_f64(1e-8);
110    let two_h = h * T::from_f64(2.0);
111    let mut grad_data = vec![T::zero(); n];
112
113    let x_data = x.as_slice();
114    let mut x_plus = x_data.to_vec();
115    let mut x_minus = x_data.to_vec();
116
117    for i in 0..n {
118        let orig = x_data[i];
119
120        x_plus[i] = orig + h;
121        x_minus[i] = orig - h;
122
123        let t_plus = Tensor::from_vec(x_plus.clone(), vec![n])
124            .expect("perturbed vector length matches original dimension n");
125        let t_minus = Tensor::from_vec(x_minus.clone(), vec![n])
126            .expect("perturbed vector length matches original dimension n");
127
128        grad_data[i] = (f(&t_plus) - f(&t_minus)) / two_h;
129
130        x_plus[i] = orig;
131        x_minus[i] = orig;
132    }
133
134    Tensor::from_vec(grad_data, vec![n]).expect("gradient vector length matches dimension n")
135}
136
137#[cfg(test)]
138mod tests {
139    use super::*;
140
141    #[test]
142    fn test_numerical_gradient_quadratic() {
143        // f(x, y) = x^2 + y^2, grad = (2x, 2y)
144        let f = |x: &Tensor<f64>| {
145            let s = x.as_slice();
146            s[0] * s[0] + s[1] * s[1]
147        };
148
149        let x = Tensor::from_vec(vec![3.0, 4.0], vec![2]).unwrap();
150        let g = numerical_gradient(&f, &x);
151        let gs = g.as_slice();
152
153        assert!((gs[0] - 6.0).abs() < 1e-5, "got {}", gs[0]);
154        assert!((gs[1] - 8.0).abs() < 1e-5, "got {}", gs[1]);
155    }
156}
scivex_optim/minimize/mod.rs

scivex_optim/minimize/
mod.rs
