2 * Copyright Nick Thompson, 2019
3 * Use, modification and distribution are subject to the
4 * Boost Software License, Version 1.0. (See accompanying file
5 * LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
8 #include "math_unit_test.hpp"
11 #include <boost/math/interpolators/cardinal_quintic_b_spline.hpp>
12 #ifdef BOOST_HAS_FLOAT128
13 #include <boost/multiprecision/float128.hpp>
14 using boost::multiprecision::float128
;
16 using boost::math::interpolators::cardinal_quintic_b_spline
;
23 Real h
= Real(1)/Real(16);
25 std::vector
<Real
> v(n
, c
);
26 std::pair
<Real
, Real
> left_endpoint_derivatives
{0, 0};
27 std::pair
<Real
, Real
> right_endpoint_derivatives
{0, 0};
28 auto qbs
= cardinal_quintic_b_spline
<Real
>(v
.data(), v
.size(), t0
, h
, left_endpoint_derivatives
, right_endpoint_derivatives
);
33 CHECK_ULP_CLOSE(c
, qbs(t
), 3);
34 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.prime(t
), 400*std::numeric_limits
<Real
>::epsilon());
35 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.double_prime(t
), 60000*std::numeric_limits
<Real
>::epsilon());
41 Real t
= t0
+ i
*h
+ h
/2;
42 CHECK_ULP_CLOSE(c
, qbs(t
), 5);
43 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.prime(t
), 600*std::numeric_limits
<Real
>::epsilon());
44 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.double_prime(t
), 30000*std::numeric_limits
<Real
>::epsilon());
46 CHECK_ULP_CLOSE(c
, qbs(t
), 4);
47 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.prime(t
), 600*std::numeric_limits
<Real
>::epsilon());
48 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.double_prime(t
), 10000*std::numeric_limits
<Real
>::epsilon());
54 void test_constant_estimate_derivatives()
58 Real h
= Real(1)/Real(16);
60 std::vector
<Real
> v(n
, c
);
61 auto qbs
= cardinal_quintic_b_spline
<Real
>(v
.data(), v
.size(), t0
, h
);
66 CHECK_ULP_CLOSE(c
, qbs(t
), 3);
67 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.prime(t
), 1200*std::numeric_limits
<Real
>::epsilon());
68 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.double_prime(t
), 200000*std::numeric_limits
<Real
>::epsilon());
74 Real t
= t0
+ i
*h
+ h
/2;
75 CHECK_ULP_CLOSE(c
, qbs(t
), 8);
76 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.prime(t
), 1200*std::numeric_limits
<Real
>::epsilon());
77 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.double_prime(t
), 80000*std::numeric_limits
<Real
>::epsilon());
79 CHECK_ULP_CLOSE(c
, qbs(t
), 5);
80 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.prime(t
), 1200*std::numeric_limits
<Real
>::epsilon());
81 CHECK_MOLLIFIED_CLOSE(Real(0), qbs
.double_prime(t
), 38000*std::numeric_limits
<Real
>::epsilon());
94 Real h
= Real(1)/Real(16);
96 std::vector
<Real
> y(n
);
97 for (size_t i
= 0; i
< n
; ++i
) {
101 std::pair
<Real
, Real
> left_endpoint_derivatives
{m
, 0};
102 std::pair
<Real
, Real
> right_endpoint_derivatives
{m
, 0};
103 auto qbs
= cardinal_quintic_b_spline
<Real
>(y
.data(), y
.size(), t0
, h
, left_endpoint_derivatives
, right_endpoint_derivatives
);
108 if (!CHECK_ULP_CLOSE(m
*t
+b
, qbs(t
), 3)) {
109 std::cerr
<< " Problem at t = " << t
<< "\n";
111 if(!CHECK_MOLLIFIED_CLOSE(m
, qbs
.prime(t
), 100*abs(m
*t
+b
)*std::numeric_limits
<Real
>::epsilon())) {
112 std::cerr
<< " Problem at t = " << t
<< "\n";
114 if(!CHECK_MOLLIFIED_CLOSE(0, qbs
.double_prime(t
), 10000*abs(m
*t
+b
)*std::numeric_limits
<Real
>::epsilon())) {
115 std::cerr
<< " Problem at t = " << t
<< "\n";
122 Real t
= t0
+ i
*h
+ h
/2;
123 if(!CHECK_ULP_CLOSE(m
*t
+b
, qbs(t
), 4)) {
124 std::cerr
<< " Problem at t = " << t
<< "\n";
126 CHECK_MOLLIFIED_CLOSE(m
, qbs
.prime(t
), 1500*std::numeric_limits
<Real
>::epsilon());
128 if(!CHECK_ULP_CLOSE(m
*t
+b
, qbs(t
), 4)) {
129 std::cerr
<< " Problem at t = " << t
<< "\n";
131 CHECK_MOLLIFIED_CLOSE(m
, qbs
.prime(t
), 3000*std::numeric_limits
<Real
>::epsilon());
137 void test_linear_estimate_derivatives()
143 Real h
= Real(1)/Real(16);
145 std::vector
<Real
> y(n
);
146 for (size_t i
= 0; i
< n
; ++i
) {
151 auto qbs
= cardinal_quintic_b_spline
<Real
>(y
.data(), y
.size(), t0
, h
);
156 if (!CHECK_ULP_CLOSE(m
*t
+b
, qbs(t
), 3)) {
157 std::cerr
<< " Problem at t = " << t
<< "\n";
159 if(!CHECK_MOLLIFIED_CLOSE(m
, qbs
.prime(t
), 100*abs(m
*t
+b
)*std::numeric_limits
<Real
>::epsilon())) {
160 std::cerr
<< " Problem at t = " << t
<< "\n";
162 if(!CHECK_MOLLIFIED_CLOSE(0, qbs
.double_prime(t
), 20000*abs(m
*t
+b
)*std::numeric_limits
<Real
>::epsilon())) {
163 std::cerr
<< " Problem at t = " << t
<< "\n";
170 Real t
= t0
+ i
*h
+ h
/2;
171 if(!CHECK_ULP_CLOSE(m
*t
+b
, qbs(t
), 5)) {
172 std::cerr
<< " Problem at t = " << t
<< "\n";
174 CHECK_MOLLIFIED_CLOSE(m
, qbs
.prime(t
), 1500*std::numeric_limits
<Real
>::epsilon());
176 if(!CHECK_ULP_CLOSE(m
*t
+b
, qbs(t
), 4)) {
177 std::cerr
<< " Problem at t = " << t
<< "\n";
179 CHECK_MOLLIFIED_CLOSE(m
, qbs
.prime(t
), 3000*std::numeric_limits
<Real
>::epsilon());
186 void test_quadratic()
188 Real a
= Real(1)/Real(16);
192 Real h
= Real(1)/Real(16);
194 std::vector
<Real
> y(n
);
195 for (size_t i
= 0; i
< n
; ++i
) {
197 y
[i
] = a
*t
*t
+ b
*t
+ c
;
199 Real t_max
= t0
+ (n
-1)*h
;
200 std::pair
<Real
, Real
> left_endpoint_derivatives
{b
, 2*a
};
201 std::pair
<Real
, Real
> right_endpoint_derivatives
{2*a
*t_max
+ b
, 2*a
};
203 auto qbs
= cardinal_quintic_b_spline
<Real
>(y
, t0
, h
, left_endpoint_derivatives
, right_endpoint_derivatives
);
208 CHECK_ULP_CLOSE(a
*t
*t
+ b
*t
+ c
, qbs(t
), 3);
214 Real t
= t0
+ i
*h
+ h
/2;
215 if(!CHECK_ULP_CLOSE(a
*t
*t
+ b
*t
+ c
, qbs(t
), 5)) {
216 std::cerr
<< " Problem at abscissa t = " << t
<< "\n";
220 if (!CHECK_ULP_CLOSE(a
*t
*t
+ b
*t
+ c
, qbs(t
), 5)) {
221 std::cerr
<< " Problem abscissa t = " << t
<< "\n";
229 void test_quadratic_estimate_derivatives()
231 Real a
= Real(1)/Real(16);
235 Real h
= Real(1)/Real(16);
237 std::vector
<Real
> y(n
);
238 for (size_t i
= 0; i
< n
; ++i
) {
240 y
[i
] = a
*t
*t
+ b
*t
+ c
;
242 auto qbs
= cardinal_quintic_b_spline
<Real
>(y
, t0
, h
);
247 CHECK_ULP_CLOSE(a
*t
*t
+ b
*t
+ c
, qbs(t
), 3);
253 Real t
= t0
+ i
*h
+ h
/2;
254 if(!CHECK_ULP_CLOSE(a
*t
*t
+ b
*t
+ c
, qbs(t
), 10)) {
255 std::cerr
<< " Problem at abscissa t = " << t
<< "\n";
259 if (!CHECK_ULP_CLOSE(a
*t
*t
+ b
*t
+ c
, qbs(t
), 6)) {
260 std::cerr
<< " Problem abscissa t = " << t
<< "\n";
269 test_constant
<double>();
270 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
271 test_constant
<long double>();
274 test_constant_estimate_derivatives
<double>();
275 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
276 test_constant_estimate_derivatives
<long double>();
279 test_linear
<float>();
280 test_linear
<double>();
281 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
282 test_linear
<long double>();
285 test_linear_estimate_derivatives
<double>();
286 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
287 test_linear_estimate_derivatives
<long double>();
290 test_quadratic
<double>();
291 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
292 test_quadratic
<long double>();
295 test_quadratic_estimate_derivatives
<double>();
296 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
297 test_quadratic_estimate_derivatives
<long double>();
300 #ifdef BOOST_HAS_FLOAT128
301 test_constant
<float128
>();
302 test_linear
<float128
>();
303 test_linear_estimate_derivatives
<float128
>();
306 return boost::math::test::report_errors();