boost/math
2
0
Fork 0
mirror of https://github.com/boostorg/math.git synced 2026-02-10 23:42:23 +00:00
Code Issues Projects Releases Wiki Activity

Fix invalid read in cubic Hermite. [CI SKIP]

Browse Source
This commit is contained in:
Nick
2020-03-22 11:22:36 -04:00
parent bc70b99553
commit 5cbfdb554a
6 changed files with 164 additions and 14 deletions
Download Patch File Download Diff File Expand all files Collapse all files

93
test/cubic_hermite_test.cpp
View File

@@ -12,6 +12,7 @@
#include <array>
#include <vector>
#include <boost/math/interpolators/cubic_hermite.hpp>
#include <boost/math/special_functions/next.hpp>
#include <boost/circular_buffer.hpp>
#ifdef BOOST_HAS_FLOAT128
#include <boost/multiprecision/float128.hpp>
@@ -27,10 +28,11 @@ using boost::math::interpolators::cardinal_cubic_hermite_aos;
template<typename Real>
void test_constant()
{
std::vector<Real> x{0,1,2,3, 9, 22, 81};
Real x0 = 0;
std::vector<Real> x{x0,1,2,3, 9, 22, 81};
std::vector<Real> y(x.size());
for (auto & t : y) {
for (auto & t : y)
{
t = 7;
}
@@ -40,9 +42,19 @@ void test_constant()
auto dydx_copy = dydx;
auto hermite_spline = cubic_hermite(std::move(x_copy), std::move(y_copy), std::move(dydx_copy));
for (Real t = x[0]; t <= x.back(); t += 0.25) {
CHECK_ULP_CLOSE(Real(7), hermite_spline(t), 2);
CHECK_ULP_CLOSE(Real(0), hermite_spline.prime(t), 2);
// Now check the boundaries:
Real tlo = x.front();
Real thi = x.back();
int samples = 5000;
int i = 0;
while (i++ < samples)
{
CHECK_ULP_CLOSE(Real(7), hermite_spline(tlo), 2);
CHECK_ULP_CLOSE(Real(7), hermite_spline(thi), 2);
CHECK_ULP_CLOSE(Real(0), hermite_spline.prime(tlo), 2);
CHECK_ULP_CLOSE(Real(0), hermite_spline.prime(thi), 2);
tlo = boost::math::nextafter(tlo, std::numeric_limits<Real>::max());
thi = boost::math::nextafter(thi, std::numeric_limits<Real>::lowest());
}
boost::circular_buffer<Real> x_buf(x.size());
@@ -226,10 +238,34 @@ void test_cardinal_constant()
auto hermite_spline_aos = cardinal_cubic_hermite_aos(std::move(data), x0, dx);
for (Real t = x0; t <= x0 + 24*dx; t += 0.25) {
CHECK_ULP_CLOSE(Real(7), hermite_spline_aos(t), 2);
CHECK_ULP_CLOSE(Real(0), hermite_spline_aos.prime(t), 2);
if (!CHECK_ULP_CLOSE(Real(7), hermite_spline_aos(t), 2)) {
std::cerr << " Wrong evaluation at t = " << t << "\n";
}
if (!CHECK_ULP_CLOSE(Real(0), hermite_spline_aos.prime(t), 2)) {
std::cerr << " Wrong evaluation at t = " << t << "\n";
}
}
// Now check the boundaries:
Real tlo = x0;
Real thi = x0 + (25-1)*dx;
int samples = 5000;
int i = 0;
while (i++ < samples)
{
CHECK_ULP_CLOSE(Real(7), hermite_spline(tlo), 2);
CHECK_ULP_CLOSE(Real(7), hermite_spline(thi), 2);
CHECK_ULP_CLOSE(Real(7), hermite_spline_aos(tlo), 2);
CHECK_ULP_CLOSE(Real(7), hermite_spline_aos(thi), 2);
CHECK_ULP_CLOSE(Real(0), hermite_spline.prime(tlo), 2);
CHECK_ULP_CLOSE(Real(0), hermite_spline.prime(thi), 2);
CHECK_ULP_CLOSE(Real(0), hermite_spline_aos.prime(tlo), 2);
CHECK_ULP_CLOSE(Real(0), hermite_spline_aos.prime(thi), 2);
tlo = boost::math::nextafter(tlo, std::numeric_limits<Real>::max());
thi = boost::math::nextafter(thi, std::numeric_limits<Real>::lowest());
}
}
@@ -278,6 +314,25 @@ void test_cardinal_linear()
CHECK_ULP_CLOSE(Real(1), hermite_spline_aos.prime(t), 0);
}
Real tlo = x0;
Real thi = x0 + (45-1)*dx;
int samples = 5000;
int i = 0;
while (i++ < samples)
{
CHECK_ULP_CLOSE(Real(tlo), hermite_spline(tlo), 2);
CHECK_ULP_CLOSE(Real(thi), hermite_spline(thi), 2);
CHECK_ULP_CLOSE(Real(1), hermite_spline.prime(tlo), 2);
CHECK_ULP_CLOSE(Real(1), hermite_spline.prime(thi), 2);
CHECK_ULP_CLOSE(Real(tlo), hermite_spline_aos(tlo), 2);
CHECK_ULP_CLOSE(Real(thi), hermite_spline_aos(thi), 2);
CHECK_ULP_CLOSE(Real(1), hermite_spline_aos.prime(tlo), 2);
CHECK_ULP_CLOSE(Real(1), hermite_spline_aos.prime(thi), 2);
tlo = boost::math::nextafter(tlo, std::numeric_limits<Real>::max());
thi = boost::math::nextafter(thi, std::numeric_limits<Real>::lowest());
}
}
@@ -318,6 +373,23 @@ void test_cardinal_quadratic()
CHECK_ULP_CLOSE(t, saos.prime(t), 70);
}
auto [tlo, thi] = s.domain();
int samples = 5000;
int i = 0;
while (i++ < samples)
{
CHECK_ULP_CLOSE(Real(tlo*tlo/2), s(tlo), 3);
CHECK_ULP_CLOSE(Real(thi*thi/2), s(thi), 3);
CHECK_ULP_CLOSE(Real(tlo), s.prime(tlo), 3);
CHECK_ULP_CLOSE(Real(thi), s.prime(thi), 3);
CHECK_ULP_CLOSE(Real(tlo*tlo/2), saos(tlo), 3);
CHECK_ULP_CLOSE(Real(thi*thi/2), saos(thi), 3);
CHECK_ULP_CLOSE(Real(tlo), saos.prime(tlo), 3);
CHECK_ULP_CLOSE(Real(thi), saos.prime(thi), 3);
tlo = boost::math::nextafter(tlo, std::numeric_limits<Real>::max());
thi = boost::math::nextafter(thi, std::numeric_limits<Real>::lowest());
}
}
@@ -359,8 +431,6 @@ int main()
test_cardinal_quadratic<float>();
test_cardinal_interpolation_condition<float>();
test_constant<double>();
test_linear<double>();
test_quadratic<double>();
@@ -370,7 +440,6 @@ int main()
test_cardinal_quadratic<double>();
test_cardinal_interpolation_condition<double>();
test_constant<long double>();
test_linear<long double>();
test_quadratic<long double>();

8
test/daubechies_scaling_test.cpp
View File

@@ -427,6 +427,14 @@ void test_quadratures()
{
CHECK_ULP_CLOSE(Real(0), phi(xlo), 0);
CHECK_ULP_CLOSE(Real(0), phi(xhi), 0);
if constexpr (p > 2) {
assert(abs(phi.prime(xlo)) <= 5);
assert(abs(phi.prime(xhi)) <= 5);
if constexpr (p > 5) {
assert(abs(phi.double_prime(xlo)) <= 5);
assert(abs(phi.double_prime(xhi)) <= 5);
}
}
xlo = std::nextafter(xlo, std::numeric_limits<Real>::lowest());
xhi = std::nextafter(xhi, std::numeric_limits<Real>::max());
}

21
test/daubechies_wavelet_test.cpp
View File

@@ -11,6 +11,7 @@
#include <iomanip>
#include <iostream>
#include <random>
#include <cmath>
#include <boost/core/demangle.hpp>
#include <boost/hana/for_each.hpp>
#include <boost/hana/ext/std/integer_sequence.hpp>
@@ -81,6 +82,15 @@ void test_quadratures()
{
CHECK_ULP_CLOSE(Real(0), psi(xlo), 0);
CHECK_ULP_CLOSE(Real(0), psi(xhi), 0);
if constexpr (p > 2)
{
CHECK_ULP_CLOSE(Real(0), psi.prime(xlo), 0);
CHECK_ULP_CLOSE(Real(0), psi.prime(xhi), 0);
if constexpr (p >= 6) {
CHECK_ULP_CLOSE(Real(0), psi.double_prime(xlo), 0);
CHECK_ULP_CLOSE(Real(0), psi.double_prime(xhi), 0);
}
}
xlo = std::nextafter(xlo, std::numeric_limits<Real>::lowest());
xhi = std::nextafter(xhi, std::numeric_limits<Real>::max());
}
@@ -88,8 +98,19 @@ void test_quadratures()
xlo = a;
xhi = b;
for (int i = 0; i < samples; ++i) {
std::cout << std::setprecision(std::numeric_limits<Real>::max_digits10);
assert(abs(psi(xlo)) <= 5);
assert(abs(psi(xhi)) <= 5);
if constexpr (p > 2)
{
assert(abs(psi.prime(xlo)) <= 5);
assert(abs(psi.prime(xhi)) <= 5);
if constexpr (p >= 6)
{
assert(abs(psi.double_prime(xlo)) <= 5);
assert(abs(psi.double_prime(xhi)) <= 5);
}
}
xlo = std::nextafter(xlo, std::numeric_limits<Real>::max());
xhi = std::nextafter(xhi, std::numeric_limits<Real>::lowest());
}