argmin_testfunctions/
threehumpcamel.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//! # Three-hump camel test function
9//!
10//! Defined as
11//!
12//! $$
13//! f(x_1,\\,x_2) = 2x_1^2 - 1.05x_1^4 + \frac{x_1^6}{6} + x_1x_2 + x_2^2
14//! $$
15//!
16//! where $x_i \in [-5, 5]$.
17//!
18//! The global minimum is at $f(x_1,\\,x_2) = f(0,\\,0) = 0$.
19
20use num::{Float, FromPrimitive};
21
22/// Three-hump camel test function.
23///
24/// Defined as
25///
26/// $$
27/// f(x_1,\\,x_2) = 2x_1^2 - 1.05x_1^4 + \frac{x_1^6}{6} + x_1x_2 + x_2^2
28/// $$
29///
30/// where $x_i \in [-5, 5]$.
31///
32/// The global minimum is at $f(x_1,\\,x_2) = f(0,\\,0) = 0$.
33pub fn threehumpcamel<T>(param: &[T; 2]) -> T
34where
35    T: Float + FromPrimitive,
36{
37    let [x1, x2] = *param;
38
39    T::from_f64(2.0).unwrap() * x1.powi(2) - T::from_f64(1.05).unwrap() * x1.powi(4)
40        + x1.powi(6) / T::from_f64(6.0).unwrap()
41        + x1 * x2
42        + x2.powi(2)
43}
44
45/// Derivative of Three-hump camel test function.
46pub fn threehumpcamel_derivative<T>(param: &[T; 2]) -> [T; 2]
47where
48    T: Float + FromPrimitive,
49{
50    let [x1, x2] = *param;
51
52    let n2 = T::from_f64(2.0).unwrap();
53    let n4 = T::from_f64(4.0).unwrap();
54    let n4_2 = T::from_f64(4.2).unwrap();
55
56    [x1.powi(5) - n4_2 * x1.powi(3) + n4 * x1 + x2, n2 * x2 + x1]
57}
58
59/// Hessian of Three-hump camel test function.
60pub fn threehumpcamel_hessian<T>(param: &[T; 2]) -> [[T; 2]; 2]
61where
62    T: Float + FromPrimitive,
63{
64    let [x1, _] = *param;
65
66    let n1 = T::from_f64(1.0).unwrap();
67    let n2 = T::from_f64(2.0).unwrap();
68    let n4 = T::from_f64(4.0).unwrap();
69    let n5 = T::from_f64(5.0).unwrap();
70    let n12_6 = T::from_f64(12.6).unwrap();
71
72    let a = n5 * x1.powi(4) - n12_6 * x1.powi(2) + n4;
73    let b = n2;
74    let offdiag = n1;
75
76    [[a, offdiag], [offdiag, b]]
77}
78
79#[cfg(test)]
80mod tests {
81    use super::*;
82    use approx::assert_relative_eq;
83    use finitediff::FiniteDiff;
84    use proptest::prelude::*;
85    use std::{f32, f64};
86
87    #[test]
88    fn test_threehumpcamel_optimum() {
89        assert_relative_eq!(
90            threehumpcamel(&[0.0_f32, 0.0_f32]),
91            0.0,
92            epsilon = f32::EPSILON
93        );
94        assert_relative_eq!(
95            threehumpcamel(&[0.0_f64, 0.0_f64]),
96            0.0,
97            epsilon = f64::EPSILON
98        );
99
100        let deriv = threehumpcamel_derivative(&[0.0, 0.0]);
101        for i in 0..2 {
102            assert_relative_eq!(deriv[i], 0.0, epsilon = f64::EPSILON);
103        }
104    }
105
106    proptest! {
107        #[test]
108        fn test_threehumpcamel_derivative(a in -5.0..5.0, b in -5.0..5.0) {
109            let param = [a, b];
110            let derivative = threehumpcamel_derivative(&param);
111            let derivative_fd = Vec::from(param).central_diff(&|x| threehumpcamel(&[x[0], x[1]]));
112            for i in 0..derivative.len() {
113                assert_relative_eq!(
114                    derivative[i],
115                    derivative_fd[i],
116                    epsilon = 1e-5,
117                    max_relative = 1e-2
118                );
119            }
120        }
121    }
122
123    proptest! {
124        #[test]
125        fn test_threehumpcamel_hessian_finitediff(a in -5.0..5.0, b in -5.0..5.0) {
126            let param = [a, b];
127            let hessian = threehumpcamel_hessian(&param);
128            let hessian_fd =
129                Vec::from(param).central_hessian(&|x| threehumpcamel_derivative(&[x[0], x[1]]).to_vec());
130            let n = hessian.len();
131            // println!("1: {hessian:?} at {a}/{b}");
132            // println!("2: {hessian_fd:?} at {a}/{b}");
133            for i in 0..n {
134                assert_eq!(hessian[i].len(), n);
135                for j in 0..n {
136                    if hessian_fd[i][j].is_finite() {
137                        assert_relative_eq!(
138                            hessian[i][j],
139                            hessian_fd[i][j],
140                            epsilon = 1e-5,
141                            max_relative = 1e-2
142                        );
143                    }
144                }
145            }
146        }
147    }
148}