summaryrefslogtreecommitdiffstats
path: root/src/boost/libs/math/test/test_autodiff_1.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'src/boost/libs/math/test/test_autodiff_1.cpp')
-rw-r--r--src/boost/libs/math/test/test_autodiff_1.cpp695
1 files changed, 695 insertions, 0 deletions
diff --git a/src/boost/libs/math/test/test_autodiff_1.cpp b/src/boost/libs/math/test/test_autodiff_1.cpp
new file mode 100644
index 00000000..9417d78b
--- /dev/null
+++ b/src/boost/libs/math/test/test_autodiff_1.cpp
@@ -0,0 +1,695 @@
+// Copyright Matthew Pulver 2018 - 2019.
+// Distributed under the Boost Software License, Version 1.0.
+// (See accompanying file LICENSE_1_0.txt or copy at
+// https://www.boost.org/LICENSE_1_0.txt)
+
+#include "test_autodiff.hpp"
+
+BOOST_AUTO_TEST_SUITE(test_autodiff_1)
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(constructors, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ // Verify value-initialized instance has all 0 entries.
+ const autodiff_fvar<T, m> empty1 = autodiff_fvar<T, m>();
+ for (auto i : boost::irange(m + 1)) {
+ BOOST_CHECK_EQUAL(empty1.derivative(i), 0);
+ }
+ const auto empty2 = autodiff_fvar<T, m, n>();
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ BOOST_CHECK_EQUAL(empty2.derivative(i, j), 0);
+ }
+ }
+ // Single variable
+ const T cx = 10.0;
+ const auto x = make_fvar<T, m>(cx);
+ for (auto i : boost::irange(m + 1)) {
+ if (i == 0u) {
+ BOOST_CHECK_EQUAL(x.derivative(i), cx);
+ } else if (i == 1) {
+ BOOST_CHECK_EQUAL(x.derivative(i), 1);
+ } else {
+ BOOST_CHECK_EQUAL(x.derivative(i), 0);
+ }
+ }
+ const autodiff_fvar<T, n> xn = x;
+ for (auto i : boost::irange(n + 1)) {
+ if (i == 0) {
+ BOOST_CHECK_EQUAL(xn.derivative(i), cx);
+ } else if (i == 1) {
+ BOOST_CHECK_EQUAL(xn.derivative(i), 1);
+ } else {
+ BOOST_CHECK_EQUAL(xn.derivative(i), 0);
+ }
+ }
+ // Second independent variable
+ const T cy = 100.0;
+ const auto y = make_fvar<T, m, n>(cy);
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(y.derivative(i, j), cy);
+ } else if (i == 0 && j == 1) {
+ BOOST_CHECK_EQUAL(y.derivative(i, j), 1.0);
+ } else {
+ BOOST_CHECK_EQUAL(y.derivative(i, j), 0.0);
+ }
+ }
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(implicit_constructors, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ const autodiff_fvar<T, m> x = 3;
+ const autodiff_fvar<T, m> one = uncast_return(x);
+ const autodiff_fvar<T, m> two_and_a_half = 2.5;
+ BOOST_CHECK_EQUAL(static_cast<T>(x), 3.0);
+ BOOST_CHECK_EQUAL(static_cast<T>(one), 1.0);
+ BOOST_CHECK_EQUAL(static_cast<T>(two_and_a_half), 2.5);
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(assignment, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 10.0;
+ const T cy = 10.0;
+ autodiff_fvar<T, m, n>
+ empty; // Uninitialized variable<> may have non-zero values.
+ // Single variable
+ auto x = make_fvar<T, m>(cx);
+ empty = static_cast<decltype(empty)>(
+ x); // Test static_cast of single-variable to double-variable type.
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(empty.derivative(i, j), cx);
+ } else if (i == 1 && j == 0) {
+ BOOST_CHECK_EQUAL(empty.derivative(i, j), 1.0);
+ } else {
+ BOOST_CHECK_EQUAL(empty.derivative(i, j), 0.0);
+ }
+ }
+ }
+ auto y = make_fvar<T, m, n>(cy);
+ empty = y; // default assignment operator
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(empty.derivative(i, j), cy);
+ } else if (i == 0 && j == 1) {
+ BOOST_CHECK_EQUAL(empty.derivative(i, j), 1.0);
+ } else {
+ BOOST_CHECK_EQUAL(empty.derivative(i, j), 0.0);
+ }
+ }
+ }
+ empty = cx; // set a constant
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(empty.derivative(i, j), cx);
+ } else {
+ BOOST_CHECK_EQUAL(empty.derivative(i, j), 0.0);
+ }
+ }
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(ostream, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ const T cx = 10;
+ const auto x = make_fvar<T, m>(cx);
+ std::ostringstream ss;
+ ss << "x = " << x;
+ BOOST_CHECK_EQUAL(ss.str(), "x = depth(1)(10,1,0,0)");
+ ss.str(std::string());
+ const auto scalar = make_fvar<T,0>(cx);
+ ss << "scalar = " << scalar;
+ BOOST_CHECK_EQUAL(ss.str(), "scalar = depth(1)(10)");
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(addition_assignment, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 10.0;
+ auto sum = autodiff_fvar<T, m, n>(); // zero-initialized
+ // Single variable
+ const auto x = make_fvar<T, m>(cx);
+ sum += x;
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(sum.derivative(i, j), cx);
+ } else if (i == 1 && j == 0) {
+ BOOST_CHECK_EQUAL(sum.derivative(i, j), 1.0);
+ } else {
+ BOOST_CHECK_EQUAL(sum.derivative(i, j), 0.0);
+ }
+ }
+ }
+ // Arithmetic constant
+ const T cy = 11.0;
+ sum = 0;
+ sum += cy;
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(sum.derivative(i, j), cy);
+ } else {
+ BOOST_CHECK_EQUAL(sum.derivative(i, j), 0.0);
+ }
+ }
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(subtraction_assignment, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 10.0;
+ auto sum = autodiff_fvar<T, m, n>(); // zero-initialized
+ // Single variable
+ const auto x = make_fvar<T, m>(cx);
+ sum -= x;
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(sum.derivative(i, j), -cx);
+ } else if (i == 1 && j == 0) {
+ BOOST_CHECK_EQUAL(sum.derivative(i, j), -1.0);
+ } else {
+ BOOST_CHECK_EQUAL(sum.derivative(i, j), 0.0);
+ }
+ }
+ }
+ // Arithmetic constant
+ const T cy = 11.0;
+ sum = 0;
+ sum -= cy;
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(sum.derivative(i, j), -cy);
+ } else {
+ BOOST_CHECK_EQUAL(sum.derivative(i, j), 0.0);
+ }
+ }
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(multiplication_assignment, T, all_float_types) {
+ // Try explicit bracing based on feedback. Doesn't add very much except 26
+ // extra lines.
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 10.0;
+ auto product = autodiff_fvar<T, m, n>(1); // unit constant
+ // Single variable
+ auto x = make_fvar<T, m>(cx);
+ product *= x;
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(product.derivative(i, j), cx);
+ } else if (i == 1 && j == 0) {
+ BOOST_CHECK_EQUAL(product.derivative(i, j), 1.0);
+ } else {
+ BOOST_CHECK_EQUAL(product.derivative(i, j), 0.0);
+ }
+ }
+ }
+ // Arithmetic constant
+ const T cy = 11.0;
+ product = 1;
+ product *= cy;
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(product.derivative(i, j), cy);
+ } else {
+ BOOST_CHECK_EQUAL(product.derivative(i, j), 0.0);
+ }
+ }
+ }
+ // 0 * inf = nan
+ x = make_fvar<T, m>(0.0);
+ x *= std::numeric_limits<T>::infinity();
+ // std::cout << "x = " << x << std::endl;
+ for (auto i : boost::irange(m + 1)) {
+ if (i == 0) {
+ BOOST_CHECK(boost::math::isnan(static_cast<T>(x))); // Correct
+ // BOOST_CHECK_EQUAL(x.derivative(i) == 0.0); // Wrong. See
+ // multiply_assign_by_root_type().
+ } else if (i == 1) {
+ BOOST_CHECK(boost::math::isinf(x.derivative(i)));
+ } else {
+ BOOST_CHECK_EQUAL(x.derivative(i), 0.0);
+ }
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(division_assignment, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 16.0;
+ auto quotient = autodiff_fvar<T, m, n>(1); // unit constant
+ // Single variable
+ const auto x = make_fvar<T, m>(cx);
+ quotient /= x;
+ BOOST_CHECK_EQUAL(quotient.derivative(0, 0), 1 / cx);
+ BOOST_CHECK_EQUAL(quotient.derivative(1, 0), -1 / pow(cx, 2));
+ BOOST_CHECK_EQUAL(quotient.derivative(2, 0), 2 / pow(cx, 3));
+ BOOST_CHECK_EQUAL(quotient.derivative(3, 0), -6 / pow(cx, 4));
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(std::size_t(1), n + 1)) {
+ BOOST_CHECK_EQUAL(quotient.derivative(i, j), 0.0);
+ }
+ }
+ // Arithmetic constant
+ const T cy = 32.0;
+ quotient = 1;
+ quotient /= cy;
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(quotient.derivative(i, j), 1 / cy);
+ } else {
+ BOOST_CHECK_EQUAL(quotient.derivative(i, j), 0.0);
+ }
+ }
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(unary_signs, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 16.0;
+ autodiff_fvar<T, m, n> lhs;
+ // Single variable
+ const auto x = make_fvar<T, m>(cx);
+ lhs = static_cast<decltype(lhs)>(-x);
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(lhs.derivative(i, j), -cx);
+ } else if (i == 1 && j == 0) {
+ BOOST_CHECK_EQUAL(lhs.derivative(i, j), -1.0);
+ } else {
+ BOOST_CHECK_EQUAL(lhs.derivative(i, j), 0.0);
+ }
+ }
+ }
+ lhs = static_cast<decltype(lhs)>(+x);
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ if (i == 0 && j == 0) {
+ BOOST_CHECK_EQUAL(lhs.derivative(i, j), cx);
+ } else if (i == 1 && j == 0) {
+ BOOST_CHECK_EQUAL(lhs.derivative(i, j), 1.0);
+ } else {
+ BOOST_CHECK_EQUAL(lhs.derivative(i, j), 0.0);
+ }
+ }
+ }
+}
+
+// TODO 3 tests for 3 operator+() definitions.
+BOOST_AUTO_TEST_CASE_TEMPLATE(cast_double, T, all_float_types) {
+ const T ca(13);
+ const T i(12);
+ constexpr std::size_t m = 3;
+ const auto x = make_fvar<T, m>(ca);
+ BOOST_CHECK_LT(i, x);
+ BOOST_CHECK_EQUAL(i * x, i * ca);
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(int_double_casting, T, all_float_types) {
+ const T ca = 3.0;
+ const auto x0 = make_fvar<T, 0>(ca);
+ BOOST_CHECK_EQUAL(static_cast<T>(x0), ca);
+ const auto x1 = make_fvar<T, 1>(ca);
+ BOOST_CHECK_EQUAL(static_cast<T>(x1), ca);
+ const auto x2 = make_fvar<T, 2>(ca);
+ BOOST_CHECK_EQUAL(static_cast<T>(x2), ca);
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(scalar_addition, T, all_float_types) {
+ const T ca = 3.0;
+ const T cb = 4.0;
+ const auto sum0 = autodiff_fvar<T, 0>(ca) + autodiff_fvar<T, 0>(cb);
+ BOOST_CHECK_EQUAL(ca + cb, static_cast<T>(sum0));
+ const auto sum1 = autodiff_fvar<T, 0>(ca) + cb;
+ BOOST_CHECK_EQUAL(ca + cb, static_cast<T>(sum1));
+ const auto sum2 = ca + autodiff_fvar<T, 0>(cb);
+ BOOST_CHECK_EQUAL(ca + cb, static_cast<T>(sum2));
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(power8, T, all_float_types) {
+ constexpr std::size_t n = 8u;
+ const T ca = 3.0;
+ auto x = make_fvar<T, n>(ca);
+ // Test operator*=()
+ x *= x;
+ x *= x;
+ x *= x;
+ const T power_factorial = boost::math::factorial<T>(n);
+ for (auto i : boost::irange(n + 1)) {
+ BOOST_CHECK_CLOSE(
+ static_cast<T>(x.derivative(i)),
+ static_cast<T>(power_factorial /
+ boost::math::factorial<T>(static_cast<unsigned>(n - i)) *
+ pow(ca, n - i)),
+ std::numeric_limits<T>::epsilon());
+ }
+ x = make_fvar<T, n>(ca);
+ // Test operator*()
+ x = x * x * x * x * x * x * x * x;
+ for (auto i : boost::irange(n + 1)) {
+ BOOST_CHECK_CLOSE(
+ x.derivative(i),
+ power_factorial /
+ boost::math::factorial<T>(static_cast<unsigned>(n - i)) *
+ pow(ca, n - i),
+ std::numeric_limits<T>::epsilon());
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(dim1_multiplication, T, all_float_types) {
+ constexpr std::size_t m = 2;
+ constexpr std::size_t n = 3;
+ const T cy = 4.0;
+ auto y0 = make_fvar<T, m>(cy);
+ auto y = make_fvar<T, n>(cy);
+ y *= y0;
+ BOOST_CHECK_EQUAL(y.derivative(0), cy * cy);
+ BOOST_CHECK_EQUAL(y.derivative(1), 2 * cy);
+ BOOST_CHECK_EQUAL(y.derivative(2), 2.0);
+ BOOST_CHECK_EQUAL(y.derivative(3), 0.0);
+ y = y * cy;
+ BOOST_CHECK_EQUAL(y.derivative(0), cy * cy * cy);
+ BOOST_CHECK_EQUAL(y.derivative(1), 2 * cy * cy);
+ BOOST_CHECK_EQUAL(y.derivative(2), 2.0 * cy);
+ BOOST_CHECK_EQUAL(y.derivative(3), 0.0);
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(dim1and2_multiplication, T, all_float_types) {
+ constexpr std::size_t m = 2;
+ constexpr std::size_t n = 3;
+ const T cx = 3.0;
+ const T cy = 4.0;
+ auto x = make_fvar<T, m>(cx);
+ auto y = make_fvar<T, m, n>(cy);
+ y *= x;
+ BOOST_CHECK_EQUAL(y.derivative(0, 0), cx * cy);
+ BOOST_CHECK_EQUAL(y.derivative(0, 1), cx);
+ BOOST_CHECK_EQUAL(y.derivative(1, 0), cy);
+ BOOST_CHECK_EQUAL(y.derivative(1, 1), 1.0);
+ for (auto i : boost::irange(std::size_t(1), m)) {
+ for (auto j : boost::irange(std::size_t(1), n)) {
+ if (i == 1 && j == 1) {
+ BOOST_CHECK_EQUAL(y.derivative(i, j), 1.0);
+ } else {
+ BOOST_CHECK_EQUAL(y.derivative(i, j), 0.0);
+ }
+ }
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(dim2_addition, T, all_float_types) {
+ constexpr std::size_t m = 2;
+ constexpr std::size_t n = 3;
+ const T cx = 3.0;
+ const auto x = make_fvar<T, m>(cx);
+ BOOST_CHECK_EQUAL(x.derivative(0), cx);
+ BOOST_CHECK_EQUAL(x.derivative(1), 1.0);
+ BOOST_CHECK_EQUAL(x.derivative(2), 0.0);
+ const T cy = 4.0;
+ const auto y = make_fvar<T, m, n>(cy);
+ BOOST_CHECK_EQUAL(static_cast<T>(y.derivative(0)), cy);
+ BOOST_CHECK_EQUAL(static_cast<T>(y.derivative(1)),
+ 0.0); // partial of y w.r.t. x.
+
+ BOOST_CHECK_EQUAL(y.derivative(0, 0), cy);
+ BOOST_CHECK_EQUAL(y.derivative(0, 1), 1.0);
+ BOOST_CHECK_EQUAL(y.derivative(1, 0), 0.0);
+ BOOST_CHECK_EQUAL(y.derivative(1, 1), 0.0);
+ const auto z = x + y;
+ BOOST_CHECK_EQUAL(z.derivative(0, 0), cx + cy);
+ BOOST_CHECK_EQUAL(z.derivative(0, 1), 1.0);
+ BOOST_CHECK_EQUAL(z.derivative(1, 0), 1.0);
+ BOOST_CHECK_EQUAL(z.derivative(1, 1), 0.0);
+ // The following 4 are unnecessarily more expensive than the previous 4.
+ BOOST_CHECK_EQUAL(z.derivative(0).derivative(0), cx + cy);
+ BOOST_CHECK_EQUAL(z.derivative(0).derivative(1), 1.0);
+ BOOST_CHECK_EQUAL(z.derivative(1).derivative(0), 1.0);
+ BOOST_CHECK_EQUAL(z.derivative(1).derivative(1), 0.0);
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(dim2_multiplication, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 6.0;
+ const auto x = make_fvar<T, m>(cx);
+ const T cy = 5.0;
+ const auto y = make_fvar<T, 0, n>(cy);
+ const auto z = x * x * y * y * y;
+ BOOST_CHECK_EQUAL(z.derivative(0, 0), cx * cx * cy * cy * cy); // x^2 * y^3
+ BOOST_CHECK_EQUAL(z.derivative(0, 1), cx * cx * 3 * cy * cy); // x^2 * 3y^2
+ BOOST_CHECK_EQUAL(z.derivative(0, 2), cx * cx * 6 * cy); // x^2 * 6y
+ BOOST_CHECK_EQUAL(z.derivative(0, 3), cx * cx * 6); // x^2 * 6
+ BOOST_CHECK_EQUAL(z.derivative(0, 4), 0.0); // x^2 * 0
+ BOOST_CHECK_EQUAL(z.derivative(1, 0), 2 * cx * cy * cy * cy); // 2x * y^3
+ BOOST_CHECK_EQUAL(z.derivative(1, 1), 2 * cx * 3 * cy * cy); // 2x * 3y^2
+ BOOST_CHECK_EQUAL(z.derivative(1, 2), 2 * cx * 6 * cy); // 2x * 6y
+ BOOST_CHECK_EQUAL(z.derivative(1, 3), 2 * cx * 6); // 2x * 6
+ BOOST_CHECK_EQUAL(z.derivative(1, 4), 0.0); // 2x * 0
+ BOOST_CHECK_EQUAL(z.derivative(2, 0), 2 * cy * cy * cy); // 2 * y^3
+ BOOST_CHECK_EQUAL(z.derivative(2, 1), 2 * 3 * cy * cy); // 2 * 3y^2
+ BOOST_CHECK_EQUAL(z.derivative(2, 2), 2 * 6 * cy); // 2 * 6y
+ BOOST_CHECK_EQUAL(z.derivative(2, 3), 2 * 6); // 2 * 6
+ BOOST_CHECK_EQUAL(z.derivative(2, 4), 0.0); // 2 * 0
+ BOOST_CHECK_EQUAL(z.derivative(3, 0), 0.0); // 0 * y^3
+ BOOST_CHECK_EQUAL(z.derivative(3, 1), 0.0); // 0 * 3y^2
+ BOOST_CHECK_EQUAL(z.derivative(3, 2), 0.0); // 0 * 6y
+ BOOST_CHECK_EQUAL(z.derivative(3, 3), 0.0); // 0 * 6
+ BOOST_CHECK_EQUAL(z.derivative(3, 4), 0.0); // 0 * 0
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(dim2_multiplication_and_subtraction, T,
+ all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 6.0;
+ const auto x = make_fvar<T, m>(cx);
+ const T cy = 5.0;
+ const auto y = make_fvar<T, 0, n>(cy);
+ const auto z = x * x - y * y;
+ BOOST_CHECK_EQUAL(z.derivative(0, 0), cx * cx - cy * cy);
+ BOOST_CHECK_EQUAL(z.derivative(0, 1), -2 * cy);
+ BOOST_CHECK_EQUAL(z.derivative(0, 2), -2.0);
+ BOOST_CHECK_EQUAL(z.derivative(0, 3), 0.0);
+ BOOST_CHECK_EQUAL(z.derivative(0, 4), 0.0);
+ BOOST_CHECK_EQUAL(z.derivative(1, 0), 2 * cx);
+ BOOST_CHECK_EQUAL(z.derivative(2, 0), 2.0);
+ for (auto i : boost::irange(std::size_t(1), m + 1)) {
+ for (auto j : boost::irange(std::size_t(1), n + 1)) {
+ BOOST_CHECK_EQUAL(z.derivative(i, j), 0.0);
+ }
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(inverse, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ const T cx = 4.0;
+ const auto x = make_fvar<T, m>(cx);
+ const auto xinv = x.inverse();
+ BOOST_CHECK_EQUAL(xinv.derivative(0), 1 / cx);
+ BOOST_CHECK_EQUAL(xinv.derivative(1), -1 / pow(cx, 2));
+ BOOST_CHECK_EQUAL(xinv.derivative(2), 2 / pow(cx, 3));
+ BOOST_CHECK_EQUAL(xinv.derivative(3), -6 / pow(cx, 4));
+ const auto zero = make_fvar<T, m>(0);
+ const auto inf = zero.inverse();
+ for (auto i : boost::irange(m + 1)) {
+ BOOST_CHECK_EQUAL(inf.derivative(i),
+ (i % 2 == 1 ? -1 : 1) *
+ std::numeric_limits<T>::infinity());
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(division, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 16.0;
+ auto x = make_fvar<T, m>(cx);
+ const T cy = 4.0;
+ auto y = make_fvar<T, 1, n>(cy);
+ auto z = x * x / (y * y);
+ BOOST_CHECK_EQUAL(z.derivative(0, 0), cx * cx / (cy * cy)); // x^2 * y^-2
+ BOOST_CHECK_EQUAL(z.derivative(0, 1), cx * cx * (-2) * pow(cy, -3));
+ BOOST_CHECK_EQUAL(z.derivative(0, 2), cx * cx * (6) * pow(cy, -4));
+ BOOST_CHECK_EQUAL(z.derivative(0, 3), cx * cx * (-24) * pow(cy, -5));
+ BOOST_CHECK_EQUAL(z.derivative(0, 4), cx * cx * (120) * pow(cy, -6));
+ BOOST_CHECK_EQUAL(z.derivative(1, 0), 2 * cx / (cy * cy));
+ BOOST_CHECK_EQUAL(z.derivative(1, 1), 2 * cx * (-2) * pow(cy, -3));
+ BOOST_CHECK_EQUAL(z.derivative(1, 2), 2 * cx * (6) * pow(cy, -4));
+ BOOST_CHECK_EQUAL(z.derivative(1, 3), 2 * cx * (-24) * pow(cy, -5));
+ BOOST_CHECK_EQUAL(z.derivative(1, 4), 2 * cx * (120) * pow(cy, -6));
+ BOOST_CHECK_EQUAL(z.derivative(2, 0), 2 / (cy * cy));
+ BOOST_CHECK_EQUAL(z.derivative(2, 1), 2 * (-2) * pow(cy, -3));
+ BOOST_CHECK_EQUAL(z.derivative(2, 2), 2 * (6) * pow(cy, -4));
+ BOOST_CHECK_EQUAL(z.derivative(2, 3), 2 * (-24) * pow(cy, -5));
+ BOOST_CHECK_EQUAL(z.derivative(2, 4), 2 * (120) * pow(cy, -6));
+ for (auto j : boost::irange(n + 1)) {
+ BOOST_CHECK_EQUAL(z.derivative(3, j), 0.0);
+ }
+
+ auto x1 = make_fvar<T, m>(cx);
+ auto z1 = x1 / cy;
+ BOOST_CHECK_EQUAL(z1.derivative(0), cx / cy);
+ BOOST_CHECK_EQUAL(z1.derivative(1), 1 / cy);
+ BOOST_CHECK_EQUAL(z1.derivative(2), 0.0);
+ BOOST_CHECK_EQUAL(z1.derivative(3), 0.0);
+ auto y2 = make_fvar<T, m, n>(cy);
+ auto z2 = cx / y2;
+ BOOST_CHECK_EQUAL(z2.derivative(0, 0), cx / cy);
+ BOOST_CHECK_EQUAL(z2.derivative(0, 1), -cx / pow(cy, 2));
+ BOOST_CHECK_EQUAL(z2.derivative(0, 2), 2 * cx / pow(cy, 3));
+ BOOST_CHECK_EQUAL(z2.derivative(0, 3), -6 * cx / pow(cy, 4));
+ BOOST_CHECK_EQUAL(z2.derivative(0, 4), 24 * cx / pow(cy, 5));
+ for (auto i : boost::irange(std::size_t(1), m + 1)) {
+ for (auto j : boost::irange(n + 1)) {
+ BOOST_CHECK_EQUAL(z2.derivative(i, j), 0.0);
+ }
+ }
+
+ const auto z3 = y / x;
+ BOOST_CHECK_EQUAL(z3.derivative(0, 0), cy / cx);
+ BOOST_CHECK_EQUAL(z3.derivative(0, 1), 1 / cx);
+ BOOST_CHECK_EQUAL(z3.derivative(1, 0), -cy / pow(cx, 2));
+ BOOST_CHECK_EQUAL(z3.derivative(1, 1), -1 / pow(cx, 2));
+ BOOST_CHECK_EQUAL(z3.derivative(2, 0), 2 * cy / pow(cx, 3));
+ BOOST_CHECK_EQUAL(z3.derivative(2, 1), 2 / pow(cx, 3));
+ BOOST_CHECK_EQUAL(z3.derivative(3, 0), -6 * cy / pow(cx, 4));
+ BOOST_CHECK_EQUAL(z3.derivative(3, 1), -6 / pow(cx, 4));
+ for (auto i : boost::irange(m + 1)) {
+ for (auto j : boost::irange(std::size_t(2), n + 1)) {
+ BOOST_CHECK_EQUAL(z3.derivative(i, j), 0.0);
+ }
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(equality, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 10.0;
+ const T cy = 10.0;
+ const auto x = make_fvar<T, m>(cx);
+ const auto y = make_fvar<T, 0, n>(cy);
+ BOOST_CHECK_EQUAL(x, y);
+ BOOST_CHECK_EQUAL(x, cy);
+ BOOST_CHECK_EQUAL(cx, y);
+ BOOST_CHECK_EQUAL(cy, x);
+ BOOST_CHECK_EQUAL(y, cx);
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(inequality, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 10.0;
+ const T cy = 11.0;
+ const auto x = make_fvar<T, m>(cx);
+ const auto y = make_fvar<T, 0, n>(cy);
+ BOOST_CHECK_NE(x, y);
+ BOOST_CHECK_NE(x, cy);
+ BOOST_CHECK_NE(cx, y);
+ BOOST_CHECK_NE(cy, x);
+ BOOST_CHECK_NE(y, cx);
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(less_than_or_equal_to, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 10.0;
+ const T cy = 11.0;
+ const auto x = make_fvar<T, m>(cx);
+ const auto y = make_fvar<T, 0, n>(cy);
+ BOOST_CHECK_LE(x, y);
+ BOOST_CHECK_LE(x, y - 1);
+ BOOST_CHECK_LT(x, y);
+ BOOST_CHECK_LE(x, cy);
+ BOOST_CHECK_LE(x, cy - 1);
+ BOOST_CHECK_LT(x, cy);
+ BOOST_CHECK_LE(cx, y);
+ BOOST_CHECK_LE(cx, y - 1);
+ BOOST_CHECK_LT(cx, y);
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(greater_than_or_equal_to, T, all_float_types) {
+ constexpr std::size_t m = 3;
+ constexpr std::size_t n = 4;
+ const T cx = 11.0;
+ const T cy = 10.0;
+ const auto x = make_fvar<T, m>(cx);
+ const auto y = make_fvar<T, 0, n>(cy);
+ BOOST_CHECK_GE(x, y);
+ BOOST_CHECK_GE(x, y + 1);
+ BOOST_CHECK_GT(x, y);
+ BOOST_CHECK_GE(x, cy);
+ BOOST_CHECK_GE(x, cy + 1);
+ BOOST_CHECK_GT(x, cy);
+ BOOST_CHECK_GE(cx, y);
+ BOOST_CHECK_GE(cx, y + 1);
+ BOOST_CHECK_GT(cx, y);
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(fabs_test, T, all_float_types) {
+ using bmp::fabs;
+ using detail::fabs;
+ using std::fabs;
+ constexpr std::size_t m = 3;
+ const T cx = 11.0;
+ const auto x = make_fvar<T, m>(cx);
+ auto a = fabs(x);
+ BOOST_CHECK_EQUAL(a.derivative(0), fabs(cx));
+ BOOST_CHECK_EQUAL(a.derivative(1), 1.0);
+ BOOST_CHECK_EQUAL(a.derivative(2), 0.0);
+ BOOST_CHECK_EQUAL(a.derivative(3), 0.0);
+ a = fabs(-x);
+ BOOST_CHECK_EQUAL(a.derivative(0), fabs(cx));
+ BOOST_CHECK_EQUAL(a.derivative(1), 1.0); // fabs(-x) = fabs(x)
+ BOOST_CHECK_EQUAL(a.derivative(2), 0.0);
+ BOOST_CHECK_EQUAL(a.derivative(3), 0.0);
+ const auto xneg = make_fvar<T, m>(-cx);
+ a = fabs(xneg);
+ BOOST_CHECK_EQUAL(a.derivative(0), fabs(cx));
+ BOOST_CHECK_EQUAL(a.derivative(1), -1.0);
+ BOOST_CHECK_EQUAL(a.derivative(2), 0.0);
+ BOOST_CHECK_EQUAL(a.derivative(3), 0.0);
+ const auto zero = make_fvar<T, m>(0);
+ a = fabs(zero);
+ for (auto i : boost::irange(m + 1)) {
+ BOOST_CHECK_EQUAL(a.derivative(i), 0.0);
+ }
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(ceil_and_floor, T, all_float_types) {
+ using bmp::ceil;
+ using bmp::floor;
+ using std::ceil;
+ using std::floor;
+ constexpr std::size_t m = 3;
+ T tests[]{-1.5, 0.0, 1.5};
+ for (T &test : tests) {
+ const auto x = make_fvar<T, m>(test);
+ auto c = ceil(x);
+ auto f = floor(x);
+ BOOST_CHECK_EQUAL(c.derivative(0), ceil(test));
+ BOOST_CHECK_EQUAL(f.derivative(0), floor(test));
+ for (auto i : boost::irange(std::size_t(1), m + 1)) {
+ BOOST_CHECK_EQUAL(c.derivative(i), 0.0);
+ BOOST_CHECK_EQUAL(f.derivative(i), 0.0);
+ }
+ }
+}
+
+BOOST_AUTO_TEST_SUITE_END()