From ef2420b2c48c50b13a6f329ef820d15c138ef798 Mon Sep 17 00:00:00 2001 From: Sietze van Buuren Date: Tue, 29 Apr 2025 14:29:43 +0200 Subject: [PATCH] test: Generalize assertion helper for N dimensions Signed-off-by: Sietze van Buuren --- include/linear_interp.hpp | 10 +++ tests/assertion_helpers.hpp | 158 ++++++++++++++++++++++++++-------- tests/test_linear_interp.cpp | 162 +++++++++++++++++++++++++++++++++++ 3 files changed, 293 insertions(+), 37 deletions(-) diff --git a/include/linear_interp.hpp b/include/linear_interp.hpp index 5c0dec2..b8cab8a 100644 --- a/include/linear_interp.hpp +++ b/include/linear_interp.hpp @@ -280,6 +280,16 @@ class LinearInterp3D : public LinearInterpND { ~LinearInterp3D() { } }; +template +class LinearInterp4D : public LinearInterpND { + using Vector = std::vector; + using Vector4 = cip::VectorN; +public: + explicit LinearInterp4D(const Vector &x, const Vector &y, const Vector &z, const Vector &w, const Vector4 &f) + : LinearInterpND(f, x, y, z, w) + {} + ~LinearInterp4D() { } +}; } // namespace cip diff --git a/tests/assertion_helpers.hpp b/tests/assertion_helpers.hpp index d93f804..cf7282d 100644 --- a/tests/assertion_helpers.hpp +++ b/tests/assertion_helpers.hpp @@ -1,63 +1,147 @@ #include #include - +#include +#include +#include +#include +#include +#include namespace cip { - using Vector = std::vector; using Vector2 = std::vector; using Vector3 = std::vector; +using Vector4 = std::vector; constexpr double TOLERANCE_DEFAULT = 5.0e-12; -template -testing::AssertionResult Interp1DAssertions(Vector x, Vector f, Vector x_fine, Vector f_fine, double tol=TOLERANCE_DEFAULT) { - T interp(x, f); - for ( auto i = 0; i < x_fine.size(); i++ ) { - auto val = interp.eval(x_fine[i]); - if (!testing::internal::DoubleNearPredFormat("expected", "actual", "tolerance", val, f_fine[i], tol) ) { - return testing::AssertionFailure() - << "for x = " << x_fine[i] << " expected " << f_fine[i] << " but got " << val; +// Helper class for N-dimensional interpolation assertions +template +class InterpNDAssertions { +private: + // Helper to access nested vector value at specific indices + template + static double get_value(const VecT& val) { + return val; + } + + template + static double get_value(const VecT& vec, size_t idx, Indices... rest) { + return get_value(vec[idx], rest...); + } + + // Helper for dimension name (x, y, z, w, etc.) + static constexpr char dim_name(size_t i) { + return i < 3 ? 'x' + i : 'w' + (i - 3); + } + + // Recursive template to build nested loops over N dimensions + template + static testing::AssertionResult check_points( + const T& interp, + const std::array& fine_grid_vectors, + const FuncValT& f_fine, + std::array& indices, + std::array& coords, + double tol) { + + // Base case: we've set up all N dimension indices + if constexpr (Dim == N) { + // Call eval with the coordinates using apply + double val = std::apply([&interp](auto... args) { + return interp.eval(args...); + }, coords); + + // Get expected value from nested vector using indices + double expected = std::apply([&f_fine](auto... args) { + return get_value(f_fine, args...); + }, indices); + + // Check if the values match within tolerance + if (!testing::internal::DoubleNearPredFormat("expected", "actual", "tolerance", val, expected, tol)) { + testing::AssertionResult failure = testing::AssertionFailure(); + + // Format dimensions for error message (x=1.23, y=4.56, etc.) + for (size_t i = 0; i < N; ++i) { + failure << (i == 0 ? "for " : ", ") << dim_name(i) << " = " << coords[i]; + } + + failure << " expected " << expected << " but got " << val; + return failure; + } + + return testing::AssertionSuccess(); + } + // Recursive case: loop over the current dimension + else { + for (size_t i = 0; i < fine_grid_vectors[Dim].size(); ++i) { + indices[Dim] = i; + coords[Dim] = fine_grid_vectors[Dim][i]; + + // Recursively process the next dimension + auto result = check_points(interp, fine_grid_vectors, f_fine, indices, coords, tol); + if (!result) { + return result; + } + } + return testing::AssertionSuccess(); } } - return testing::AssertionSuccess(); -} +public: + // Main assertion method + static testing::AssertionResult test( + const std::array& grid_vectors, + const FuncValT& f, + const std::array& fine_grid_vectors, + const FuncValT& f_fine, + double tol = TOLERANCE_DEFAULT) { + + // Create the interpolator with a simple lambda + std::apply + T interp = std::apply([&f](const auto&... grid_args) { + return T(grid_args..., f); + }, grid_vectors); + + // Set up indices and coordinates arrays + std::array indices{}; + std::array coords{}; + + // Start the recursive dimension traversal + return check_points<0>(interp, fine_grid_vectors, f_fine, indices, coords, tol); + } +}; +// Convenience wrapper functions to maintain backward compatibility +template +testing::AssertionResult Interp1DAssertions(Vector x, Vector f, Vector x_fine, Vector f_fine, double tol=TOLERANCE_DEFAULT) { + std::array grid_vectors = {x}; + std::array fine_grid_vectors = {x_fine}; + return InterpNDAssertions::test(grid_vectors, f, fine_grid_vectors, f_fine, tol); +} template testing::AssertionResult Interp2DAssertions(const Vector &x, const Vector &y, const Vector2 &f, const Vector &x_fine, const Vector &y_fine, const Vector2 &f_fine, double tol=TOLERANCE_DEFAULT) { - T interp2(x, y, f); - for ( auto i = 0; i < x_fine.size(); ++i ) { - for ( auto j = 0; j < y_fine.size(); ++j ) { - auto val = interp2.eval(x_fine[i], y_fine[j]); - if (!testing::internal::DoubleNearPredFormat("expected", "actual", "tolerance", val, f_fine[i][j], tol) ) { - return testing::AssertionFailure() - << "for x = " << x_fine[i] << ", y = " << y_fine[j] << " expected " << f_fine[i][j] << " but got " << val; - } - } - } - return testing::AssertionSuccess(); + std::array grid_vectors = {x, y}; + std::array fine_grid_vectors = {x_fine, y_fine}; + return InterpNDAssertions::test(grid_vectors, f, fine_grid_vectors, f_fine, tol); } - template testing::AssertionResult Interp3DAssertions(const Vector &x, const Vector &y, const Vector &z, const Vector3 &f, const Vector &x_fine, const Vector &y_fine, const Vector &z_fine, const Vector3 &f_fine, double tol=TOLERANCE_DEFAULT) { - T interp3(x, y, z, f); - for ( auto i = 0; i < x_fine.size(); ++i ) { - for ( auto j = 0; j < y_fine.size(); ++j ) { - for ( auto k = 0; k < y_fine.size(); ++k ) { - auto val = interp3.eval(x_fine[i], y_fine[j], z_fine[k]); - if (!testing::internal::DoubleNearPredFormat("expected", "actual", "tolerance", val, f_fine[i][j][k], tol) ) { - return testing::AssertionFailure() - << "for x = " << x_fine[i] << ", y = " << y_fine[j] << ", z = " << z_fine[j] << " expected " << f_fine[i][j][k] << " but got " << val; - } - } - } - } - return testing::AssertionSuccess(); + std::array grid_vectors = {x, y, z}; + std::array fine_grid_vectors = {x_fine, y_fine, z_fine}; + return InterpNDAssertions::test(grid_vectors, f, fine_grid_vectors, f_fine, tol); } +template +testing::AssertionResult Interp4DAssertions(const Vector &x, const Vector &y, const Vector &z, const Vector &w, + const Vector4 &f, const Vector &x_fine, const Vector &y_fine, + const Vector &z_fine, const Vector &w_fine, const Vector4 &f_fine, + double tol=TOLERANCE_DEFAULT) { + std::array grid_vectors = {x, y, z, w}; + std::array fine_grid_vectors = {x_fine, y_fine, z_fine, w_fine}; + return InterpNDAssertions::test(grid_vectors, f, fine_grid_vectors, f_fine, tol); +} } // namespace cip diff --git a/tests/test_linear_interp.cpp b/tests/test_linear_interp.cpp index 1c16fb6..c6e86a9 100644 --- a/tests/test_linear_interp.cpp +++ b/tests/test_linear_interp.cpp @@ -7,6 +7,7 @@ using VectorN1 = cip::VectorN; using VectorN2 = cip::VectorN; using VectorN3 = cip::VectorN; +using VectorN4 = cip::VectorN; using Span = std::span; using Pr = std::pair; @@ -197,3 +198,164 @@ TEST(TestInterp3D, test_linear_interp_3d) { ASSERT_TRUE(cip::Interp3DAssertions>(x, y, z, f, x_fine, y_fine, z_fine, f_fine)); } + +/** + * Test for 4D linear interpolation + * + * This test demonstrates using the InterpNDAssertions class with N=4 dimensions. + * We use a simple 4D linear function f(x,y,z,w) = a*x + b*y + c*z + d*w + e + * to generate the test data, which should be perfectly interpolated by a linear interpolator. + */ +TEST(TestInterp4D, test_linear_interp_4d) { + // Create grid vectors for each dimension + cip::Vector x = { 0.0, 1.0, 2.0 }; + cip::Vector y = { 0.0, 1.0, 2.0 }; + cip::Vector z = { 0.0, 1.0, 2.0 }; + cip::Vector w = { 0.0, 1.0, 2.0 }; + + // Define a linear function coefficients: f(x,y,z,w) = 2*x + 3*y - z + 0.5*w + 1 + const double a = 2.0; + const double b = 3.0; + const double c = -1.0; + const double d = 0.5; + const double e = 1.0; + + // Create 4D function values array + cip::Vector4 f; + f.resize(x.size()); + for (size_t i = 0; i < x.size(); ++i) { + f[i].resize(y.size()); + for (size_t j = 0; j < y.size(); ++j) { + f[i][j].resize(z.size()); + for (size_t k = 0; k < z.size(); ++k) { + f[i][j][k].resize(w.size()); + for (size_t l = 0; l < w.size(); ++l) { + // Function: f(x,y,z,w) = 2*x + 3*y - z + 0.5*w + 1 + f[i][j][k][l] = a*x[i] + b*y[j] + c*z[k] + d*w[l] + e; + } + } + } + } + + // Fine grid for testing interpolation + cip::Vector x_fine = { 0.0, 0.5, 1.0, 1.5, 2.0 }; + cip::Vector y_fine = { 0.0, 0.5, 1.0, 1.5, 2.0 }; + cip::Vector z_fine = { 0.0, 0.5, 1.0, 1.5, 2.0 }; + cip::Vector w_fine = { 0.0, 0.5, 1.0, 1.5, 2.0 }; + + // Expected values at fine grid points + cip::Vector4 f_fine; + f_fine.resize(x_fine.size()); + for (size_t i = 0; i < x_fine.size(); ++i) { + f_fine[i].resize(y_fine.size()); + for (size_t j = 0; j < y_fine.size(); ++j) { + f_fine[i][j].resize(z_fine.size()); + for (size_t k = 0; k < z_fine.size(); ++k) { + f_fine[i][j][k].resize(w_fine.size()); + for (size_t l = 0; l < w_fine.size(); ++l) { + // Function: f(x,y,z,w) = 2*x + 3*y - z + 0.5*w + 1 + f_fine[i][j][k][l] = a*x_fine[i] + b*y_fine[j] + c*z_fine[k] + d*w_fine[l] + e; + } + } + } + } + + // Test the 4D interpolation using our InterpNDAssertions class + ASSERT_TRUE(cip::Interp4DAssertions>( + x, y, z, w, f, x_fine, y_fine, z_fine, w_fine, f_fine)); +} + +// Direct test using the InterpNDAssertions template with N=4 +TEST(TestInterp4D, test_direct_nd_assertions) { + // Create grid vectors for each dimension + cip::Vector x = { 0.0, 1.0, 2.0 }; + cip::Vector y = { 0.0, 1.0, 2.0 }; + cip::Vector z = { 0.0, 1.0, 2.0 }; + cip::Vector w = { 0.0, 1.0, 2.0 }; + + // Create 4D function values array + // Function is f(x,y,z,w) = x + y + z + w + cip::Vector4 f; + f.resize(x.size()); + for (size_t i = 0; i < x.size(); ++i) { + f[i].resize(y.size()); + for (size_t j = 0; j < y.size(); ++j) { + f[i][j].resize(z.size()); + for (size_t k = 0; k < z.size(); ++k) { + f[i][j][k].resize(w.size()); + for (size_t l = 0; l < w.size(); ++l) { + // Function: f(x,y,z,w) = x + y + z + w + f[i][j][k][l] = x[i] + y[j] + z[k] + w[l]; + } + } + } + } + + // Fine grid for testing interpolation + cip::Vector x_fine = { 0.0, 0.5, 1.0, 1.5, 2.0 }; + cip::Vector y_fine = { 0.0, 0.5, 1.0, 1.5, 2.0 }; + cip::Vector z_fine = { 0.0, 0.5, 1.0, 1.5, 2.0 }; + cip::Vector w_fine = { 0.0, 0.5, 1.0, 1.5, 2.0 }; + + // Expected values at fine grid points + cip::Vector4 f_fine; + f_fine.resize(x_fine.size()); + for (size_t i = 0; i < x_fine.size(); ++i) { + f_fine[i].resize(y_fine.size()); + for (size_t j = 0; j < y_fine.size(); ++j) { + f_fine[i][j].resize(z_fine.size()); + for (size_t k = 0; k < z_fine.size(); ++k) { + f_fine[i][j][k].resize(w_fine.size()); + for (size_t l = 0; l < w_fine.size(); ++l) { + // Function: f(x,y,z,w) = x + y + z + w + f_fine[i][j][k][l] = x_fine[i] + y_fine[j] + z_fine[k] + w_fine[l]; + } + } + } + } + + // Pack grid vectors and fine grid vectors into arrays + std::array grid_vectors = {x, y, z, w}; + std::array fine_grid_vectors = {x_fine, y_fine, z_fine, w_fine}; + + // Test the 4D interpolation directly using the InterpNDAssertions class + auto result = cip::InterpNDAssertions, 4, cip::Vector4>::test( + grid_vectors, f, fine_grid_vectors, f_fine); + ASSERT_TRUE(result); +} + +// Test a non-linear function with 4D interpolation +TEST(TestInterp4D, test_nonlinear_4d) { + // Create grid vectors for each dimension + cip::Vector x = { 0.0, 1.0, 2.0 }; + cip::Vector y = { 0.0, 1.0, 2.0 }; + cip::Vector z = { 0.0, 1.0, 2.0 }; + cip::Vector w = { 0.0, 1.0, 2.0 }; + + // Create 4D function values array for a quadratic function + // Function is f(x,y,z,w) = x^2 + y^2 + z^2 + w^2 + cip::Vector4 f; + f.resize(x.size()); + for (size_t i = 0; i < x.size(); ++i) { + f[i].resize(y.size()); + for (size_t j = 0; j < y.size(); ++j) { + f[i][j].resize(z.size()); + for (size_t k = 0; k < z.size(); ++k) { + f[i][j][k].resize(w.size()); + for (size_t l = 0; l < w.size(); ++l) { + // Quadratic function + f[i][j][k][l] = x[i]*x[i] + y[j]*y[j] + z[k]*z[k] + w[l]*w[l]; + } + } + } + } + + // Create interpolator + cip::LinearInterp4D interp(x, y, z, w, f); + + // Test a few specific points without using fine grid + // For interpolation at grid points, result should be exact + EXPECT_DOUBLE_EQ(interp.eval(0.0, 0.0, 0.0, 0.0), 0.0); + EXPECT_DOUBLE_EQ(interp.eval(1.0, 1.0, 1.0, 1.0), 4.0); + EXPECT_DOUBLE_EQ(interp.eval(2.0, 2.0, 2.0, 2.0), 16.0); +}