argmin_testfunctions/
easom.rs

1// Copyright 2018-2024 argmin developers
2//
3// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
4// http://apache.org/licenses/LICENSE-2.0> or the MIT license <LICENSE-MIT or
5// http://opensource.org/licenses/MIT>, at your option. This file may not be
6// copied, modified, or distributed except according to those terms.
7
8//! # Easom test function
9//!
10//! Defined as
11//!
12//! $$
13//! f(x_1, x_2) = - \cos(x_1)\cos(x_2)\exp\left(-(x_1 - \pi)^2 - (x_2 - \pi)^2\right)
14//! $$
15//!
16//! where $x_i \in [-100,\\,100]$.
17//!
18//! The global minimum is at $f(x_1, x_2) = f(\pi, \pi) = -1$.
19
20use num::{Float, FromPrimitive};
21use std::f64::consts::PI;
22
23/// Easom test function
24///
25/// Defined as
26///
27/// $$
28/// f(x_1, x_2) = - \cos(x_1)\cos(x_2)\exp\left(-(x_1 - \pi)^2 - (x_2 - \pi)^2\right)
29/// $$
30///
31/// where $x_i \in [-100,\\,100]$.
32///
33/// The global minimum is at $f(x_1, x_2) = f(\pi, \pi) = -1$.
34pub fn easom<T>(param: &[T; 2]) -> T
35where
36    T: Float + FromPrimitive,
37{
38    let [x1, x2] = *param;
39    let pi = T::from_f64(PI).unwrap();
40    -x1.cos() * x2.cos() * (-(x1 - pi).powi(2) - (x2 - pi).powi(2)).exp()
41}
42
43/// Derivative of Easom test function
44pub fn easom_derivative<T>(param: &[T; 2]) -> [T; 2]
45where
46    T: Float + FromPrimitive,
47{
48    let [x1, x2] = *param;
49
50    let pi = T::from_f64(PI).unwrap();
51    let n2 = T::from_f64(2.0).unwrap();
52
53    let factor = (-(x1 - pi).powi(2) - (x2 - pi).powi(2)).exp();
54
55    [
56        factor * x2.cos() * (x1.sin() + n2 * x1 * x1.cos() - n2 * pi * x1.cos()),
57        factor * x1.cos() * (x2.sin() + n2 * x2 * x2.cos() - n2 * pi * x2.cos()),
58    ]
59}
60
61/// Hessian of Easom test function
62pub fn easom_hessian<T>(param: &[T; 2]) -> [[T; 2]; 2]
63where
64    T: Float + FromPrimitive,
65{
66    let [x1, x2] = *param;
67
68    let pi = T::from_f64(PI).unwrap();
69    let n2 = T::from_f64(2.0).unwrap();
70    let n4 = T::from_f64(4.0).unwrap();
71    let n3 = T::from_f64(3.0).unwrap();
72    let n8 = T::from_f64(8.0).unwrap();
73
74    let x1cos = x1.cos();
75    let x1sin = x1.sin();
76    let x2cos = x2.cos();
77    let x2sin = x2.sin();
78    let factor = (-(x1 - pi).powi(2) - (x2 - pi).powi(2)).exp();
79    let offdiag = factor * (x1sin + n2 * (x1 - pi) * x1cos) * (n2 * (pi - x2) * x2cos - x2sin);
80
81    [
82        [
83            factor
84                * x2cos
85                * (n4 * (pi - x1) * x1sin
86                    + (-n4 * x1.powi(2) + n8 * pi * x1 - n4 * pi.powi(2) + n3) * x1cos),
87            offdiag,
88        ],
89        [
90            offdiag,
91            factor
92                * x1cos
93                * (n4 * (pi - x2) * x2sin
94                    + (-n4 * x2.powi(2) + n8 * pi * x2 - n4 * pi.powi(2) + n3) * x2cos),
95        ],
96    ]
97}
98
99#[cfg(test)]
100mod tests {
101    use super::*;
102    use approx::assert_relative_eq;
103    use finitediff::FiniteDiff;
104    use proptest::prelude::*;
105    use std::{f32, f32::consts::PI as PI32, f64, f64::consts::PI as PI64};
106
107    #[test]
108    fn test_easom_optimum() {
109        assert_relative_eq!(easom(&[PI32, PI32]), -1.0_f32, epsilon = f32::EPSILON);
110        assert_relative_eq!(easom(&[PI64, PI64]), -1.0_f64, epsilon = f64::EPSILON);
111
112        let deriv = easom_derivative(&[PI64, PI64]);
113        for i in 0..2 {
114            assert_relative_eq!(deriv[i], 0.0_f64, epsilon = f64::EPSILON);
115        }
116    }
117
118    proptest! {
119        #[test]
120        fn test_easom_derivative_finitediff(a in -100.0..100.0, b in -100.0..100.0) {
121            let param = [a, b];
122            let derivative = easom_derivative(&param);
123            let derivative_fd = Vec::from(param).central_diff(&|x| easom(&[x[0], x[1]]));
124            for i in 0..derivative.len() {
125                assert_relative_eq!(
126                    derivative[i],
127                    derivative_fd[i],
128                    epsilon = 1e-5,
129                    max_relative = 1e-2
130                );
131            }
132        }
133    }
134
135    proptest! {
136        #[test]
137        fn test_easom_hessian_finitediff(a in -100.0..100.0, b in -100.0..100.0) {
138            let param = [a, b];
139            let hessian = easom_hessian(&param);
140            let hessian_fd =
141                Vec::from(param).forward_hessian(&|x| easom_derivative(&[x[0], x[1]]).to_vec());
142            let n = hessian.len();
143            for i in 0..n {
144                assert_eq!(hessian[i].len(), n);
145                for j in 0..n {
146                    assert_relative_eq!(
147                        hessian[i][j],
148                        hessian_fd[i][j],
149                        epsilon = 1e-5,
150                        max_relative = 1e-2
151                    );
152                }
153            }
154        }
155    }
156}