From 6d907d9bd3704cc0b9685ff5fdafacb951f98951 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Thu, 21 Dec 2023 08:20:43 -0500 Subject: [PATCH 01/15] Remove unnecessary `KOKKOS_FUNCTION` and `KOKKOS_INLINE_FUNCTION` --- src/geometry/ArborX_GeometryTraits.hpp | 2 +- src/interpolation/ArborX_InterpMovingLeastSquares.hpp | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/geometry/ArborX_GeometryTraits.hpp b/src/geometry/ArborX_GeometryTraits.hpp index 6315e5e3a..7c13b9bb7 100644 --- a/src/geometry/ArborX_GeometryTraits.hpp +++ b/src/geometry/ArborX_GeometryTraits.hpp @@ -90,7 +90,7 @@ struct is_triangle : std::is_same::type, TriangleTag> {}; template -KOKKOS_FUNCTION constexpr void check_valid_geometry_traits(Geometry const &) +void check_valid_geometry_traits(Geometry const &) { static_assert( !Kokkos::is_detected{}, diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index d540d675a..65f819e93 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -48,13 +48,13 @@ struct AccessTraits< Interpolation::Details::MLSTargetPointsPredicateWrapper, PredicatesTag> { - KOKKOS_INLINE_FUNCTION static auto size( + KOKKOS_FUNCTION static auto size( Interpolation::Details::MLSTargetPointsPredicateWrapper const &tp) { return AccessTraits::size(tp.target_points); } - KOKKOS_INLINE_FUNCTION static auto + KOKKOS_FUNCTION static auto get(Interpolation::Details::MLSTargetPointsPredicateWrapper const &tp, int const i) { From 26ba073d70d78c2f25841010f062e34bc3625ece Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Thu, 21 Dec 2023 09:07:50 -0500 Subject: [PATCH 02/15] Using AccessValues instead of AccessTraits --- .../ArborX_InterpMovingLeastSquares.hpp | 71 ++++++++++--------- ...pDetailsMovingLeastSquaresCoefficients.hpp | 23 +++--- 2 files changed, 50 insertions(+), 44 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index 65f819e93..795f9cb79 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -31,10 +31,10 @@ namespace ArborX::Interpolation::Details { // This is done to avoid a clash with another predicates access trait -template +template struct MLSTargetPointsPredicateWrapper { - Points target_points; + ArborX::Details::AccessValues target_access; int num_neighbors; }; @@ -51,20 +51,18 @@ struct AccessTraits< KOKKOS_FUNCTION static auto size( Interpolation::Details::MLSTargetPointsPredicateWrapper const &tp) { - return AccessTraits::size(tp.target_points); + return tp.target_access.size(); } KOKKOS_FUNCTION static auto get(Interpolation::Details::MLSTargetPointsPredicateWrapper const &tp, int const i) { - return nearest( - AccessTraits::get(tp.target_points, i), - tp.num_neighbors); + return nearest(tp.target_access(i), tp.num_neighbors); } using memory_space = - typename AccessTraits::memory_space; + typename Details::AccessValues::memory_space; }; } // namespace ArborX @@ -95,29 +93,31 @@ class MovingLeastSquares // SourcePoints is an access trait of points ArborX::Details::check_valid_access_traits(PrimitivesTag{}, source_points); - using src_acc = AccessTraits; - static_assert(KokkosExt::is_accessible_from::value, - "Source points must be accessible from the execution space"); - using src_point = - typename ArborX::Details::AccessTraitsHelper::type; - GeometryTraits::check_valid_geometry_traits(src_point{}); - static_assert(GeometryTraits::is_point::value, + using SourceAccess = + ArborX::Details::AccessValues; + static_assert( + KokkosExt::is_accessible_from::value, + "Source points must be accessible from the execution space"); + using SourcePoint = typename SourceAccess::value_type; + GeometryTraits::check_valid_geometry_traits(SourcePoint{}); + static_assert(GeometryTraits::is_point::value, "Source points elements must be points"); - static constexpr int dimension = GeometryTraits::dimension_v; + static constexpr int dimension = GeometryTraits::dimension_v; // TargetPoints is an access trait of points ArborX::Details::check_valid_access_traits(PrimitivesTag{}, target_points); - using tgt_acc = AccessTraits; - static_assert(KokkosExt::is_accessible_from::value, - "Target points must be accessible from the execution space"); - using tgt_point = - typename ArborX::Details::AccessTraitsHelper::type; - GeometryTraits::check_valid_geometry_traits(tgt_point{}); - static_assert(GeometryTraits::is_point::value, + using TargetAccess = + ArborX::Details::AccessValues; + static_assert( + KokkosExt::is_accessible_from::value, + "Target points must be accessible from the execution space"); + using TargetPoint = typename TargetAccess::value_type; + GeometryTraits::check_valid_geometry_traits(TargetPoint{}); + static_assert(GeometryTraits::is_point::value, "Target points elements must be points"); - static_assert(dimension == GeometryTraits::dimension_v, + static_assert(dimension == GeometryTraits::dimension_v, "Target and source points must have the same dimension"); int num_neighbors_val = @@ -125,18 +125,20 @@ class MovingLeastSquares : Details::polynomialBasisSize(); - int const num_targets = tgt_acc::size(target_points); - _source_size = src_acc::size(source_points); + TargetAccess target_access{target_points}; + SourceAccess source_access{source_points}; + int const num_targets = target_access.size(); + _source_size = source_access.size(); // There must be enough source points KOKKOS_ASSERT(0 < num_neighbors_val && num_neighbors_val <= _source_size); // Organize the source points as a tree - BoundingVolumeHierarchy> + BoundingVolumeHierarchy> source_tree(space, ArborX::Experimental::attach_indices(source_points)); // Create the predicates Details::MLSTargetPointsPredicateWrapper predicates{ - target_points, num_neighbors_val}; + {target_points}, num_neighbors_val}; // Query the source Kokkos::View indices( @@ -171,18 +173,19 @@ class MovingLeastSquares auto guard = Kokkos::Profiling::ScopedRegion( "ArborX::MovingLeastSquares::fillValuesIndicesAndGetSourceView"); - using src_acc = AccessTraits; - using src_point = - typename ArborX::Details::AccessTraitsHelper::type; + using SourceAccess = + ArborX::Details::AccessValues; + using SourcePoint = typename SourceAccess::value_type; _values_indices = Kokkos::View( Kokkos::view_alloc(Kokkos::WithoutInitializing, "ArborX::MovingLeastSquares::values_indices"), num_targets, num_neighbors); - Kokkos::View source_view( + Kokkos::View source_view( Kokkos::view_alloc(Kokkos::WithoutInitializing, "ArborX::MovingLeastSquares::source_view"), num_targets, num_neighbors); + SourceAccess source_access{source_points}; Kokkos::parallel_for( "ArborX::MovingLeastSquares::values_indices_and_source_view_fill", Kokkos::MDRangePolicy>( @@ -190,7 +193,7 @@ class MovingLeastSquares KOKKOS_CLASS_LAMBDA(int const i, int const j) { auto index = indices(offsets(i) + j); _values_indices(i, j) = index; - source_view(i, j) = src_acc::get(source_points, index); + source_view(i, j) = source_access(index); }); return source_view; diff --git a/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp b/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp index 35fd048fa..5fe216712 100644 --- a/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp +++ b/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp @@ -57,18 +57,21 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, // TargetPoints is an access trait of points ArborX::Details::check_valid_access_traits(PrimitivesTag{}, target_points); - using tgt_acc = AccessTraits; - static_assert(KokkosExt::is_accessible_from::value, - "target points must be accessible from the execution space"); - using tgt_point = typename ArborX::Details::AccessTraitsHelper::type; - GeometryTraits::check_valid_geometry_traits(tgt_point{}); - static_assert(GeometryTraits::is_point::value, + using TargetAccess = + ArborX::Details::AccessValues; + static_assert( + KokkosExt::is_accessible_from::value, + "target points must be accessible from the execution space"); + using TargetPoint = typename TargetAccess::value_type; + GeometryTraits::check_valid_geometry_traits(TargetPoint{}); + static_assert(GeometryTraits::is_point::value, "target points elements must be points"); - static_assert(dimension == GeometryTraits::dimension_v, + static_assert(dimension == GeometryTraits::dimension_v, "target and source points must have the same dimension"); - int const num_targets = tgt_acc::size(target_points); + TargetAccess target_access{target_points}; + int const num_targets = target_access.size(); int const num_neighbors = source_points.extent(1); // There must be a set of neighbors for each target @@ -104,7 +107,7 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, space, {0, 0}, {num_targets, num_neighbors}), KOKKOS_LAMBDA(int const i, int const j) { auto src = source_points(i, j); - auto tgt = tgt_acc::get(target_points, i); + auto tgt = target_access(i); point_t t{}; for (int k = 0; k < dimension; k++) From 612bb43132a48cce5fceb861b1c63868f8a4fc1f Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Thu, 21 Dec 2023 09:31:19 -0500 Subject: [PATCH 03/15] Slice lengths fix --- .../details/ArborX_InterpDetailsPolynomialBasis.hpp | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/src/interpolation/details/ArborX_InterpDetailsPolynomialBasis.hpp b/src/interpolation/details/ArborX_InterpDetailsPolynomialBasis.hpp index 84776973f..42dbc995b 100644 --- a/src/interpolation/details/ArborX_InterpDetailsPolynomialBasis.hpp +++ b/src/interpolation/details/ArborX_InterpDetailsPolynomialBasis.hpp @@ -116,9 +116,8 @@ KOKKOS_FUNCTION auto evaluatePolynomialBasis(Point const &p) if constexpr (Degree > 0) { // Cannot use structured binding with constexpr - static constexpr auto slice_lengths_struct = + static constexpr auto slice_lengths = polynomialBasisSliceLengths(); - auto &slice_lengths = slice_lengths_struct.arr; std::size_t prev_col = 0; std::size_t curr_col = 1; @@ -129,10 +128,10 @@ KOKKOS_FUNCTION auto evaluatePolynomialBasis(Point const &p) for (std::size_t dim = 0; dim < DIM; dim++) { // copy the previous column and multply by p[dim] - for (std::size_t i = 0; i < slice_lengths[deg][dim]; i++) + for (std::size_t i = 0; i < slice_lengths.arr[deg][dim]; i++) arr[loc_offset + i] = arr[prev_col + i] * p[dim]; - loc_offset += slice_lengths[deg][dim]; + loc_offset += slice_lengths.arr[deg][dim]; } prev_col = curr_col; From ad84b1e47c6179c2173a9d69809779a10bbaa0c4 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Thu, 21 Dec 2023 09:42:56 -0500 Subject: [PATCH 04/15] PascalCase on internal types --- .../ArborX_InterpMovingLeastSquares.hpp | 4 +-- ...pDetailsMovingLeastSquaresCoefficients.hpp | 18 +++++----- .../ArborX_InterpDetailsPolynomialBasis.hpp | 6 ++-- ...InterpDetailsSymmetricPseudoInverseSVD.hpp | 24 ++++++------- test/tstInterpDetailsMLSCoefficients.cpp | 34 +++++++++---------- test/tstInterpDetailsPolyBasis.cpp | 22 ++++++------ test/tstInterpDetailsSVD.cpp | 12 +++---- test/tstInterpMovingLeastSquares.cpp | 34 +++++++++---------- 8 files changed, 77 insertions(+), 77 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index 795f9cb79..1d384a1c3 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -236,7 +236,7 @@ class MovingLeastSquares // original input KOKKOS_ASSERT(_source_size == source_values.extent_int(0)); - using value_t = typename ApproxValues::non_const_value_type; + using Value = typename ApproxValues::non_const_value_type; int const num_targets = _values_indices.extent(0); int const num_neighbors = _values_indices.extent(1); @@ -247,7 +247,7 @@ class MovingLeastSquares "ArborX::MovingLeastSquares::target_interpolation", Kokkos::RangePolicy(space, 0, num_targets), KOKKOS_CLASS_LAMBDA(int const i) { - value_t tmp = 0; + Value tmp = 0; for (int j = 0; j < num_neighbors; j++) tmp += _coeffs(i, j) * source_values(_values_indices(i, j)); approx_values(i) = tmp; diff --git a/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp b/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp index 5fe216712..6dfc1276d 100644 --- a/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp +++ b/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp @@ -49,11 +49,11 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, KokkosExt::is_accessible_from::value, "source points must be accessible from the execution space"); - using src_point = typename SourcePoints::non_const_value_type; - GeometryTraits::check_valid_geometry_traits(src_point{}); - static_assert(GeometryTraits::is_point::value, + using SourcePoint = typename SourcePoints::non_const_value_type; + GeometryTraits::check_valid_geometry_traits(SourcePoint{}); + static_assert(GeometryTraits::is_point::value, "source points elements must be points"); - static constexpr int dimension = GeometryTraits::dimension_v; + static constexpr int dimension = GeometryTraits::dimension_v; // TargetPoints is an access trait of points ArborX::Details::check_valid_access_traits(PrimitivesTag{}, target_points); @@ -77,7 +77,7 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, // There must be a set of neighbors for each target KOKKOS_ASSERT(num_targets == source_points.extent_int(0)); - using point_t = ExperimentalHyperGeometry::Point; + using Point = ExperimentalHyperGeometry::Point; static constexpr auto epsilon = Kokkos::Experimental::epsilon_v; static constexpr int degree = PolynomialDegree::value; @@ -96,7 +96,7 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, // We first change the origin of the evaluation to be at the target point. // This lets us use p(0) which is [1 0 ... 0]. - Kokkos::View source_ref_target( + Kokkos::View source_ref_target( Kokkos::view_alloc( space, Kokkos::WithoutInitializing, "ArborX::MovingLeastSquaresCoefficients::source_ref_target"), @@ -108,7 +108,7 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, KOKKOS_LAMBDA(int const i, int const j) { auto src = source_points(i, j); auto tgt = target_access(i); - point_t t{}; + Point t{}; for (int k = 0; k < dimension; k++) t[k] = src[k] - tgt[k]; @@ -135,7 +135,7 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, for (int j = 0; j < num_neighbors; j++) { CoefficientsType norm = - ArborX::Details::distance(source_ref_target(i, j), point_t{}); + ArborX::Details::distance(source_ref_target(i, j), Point{}); radius = Kokkos::max(radius, norm); } @@ -158,7 +158,7 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, space, {0, 0}, {num_targets, num_neighbors}), KOKKOS_LAMBDA(int const i, int const j) { CoefficientsType norm = - ArborX::Details::distance(source_ref_target(i, j), point_t{}); + ArborX::Details::distance(source_ref_target(i, j), Point{}); phi(i, j) = CRBF::evaluate(norm / radii(i)); }); diff --git a/src/interpolation/details/ArborX_InterpDetailsPolynomialBasis.hpp b/src/interpolation/details/ArborX_InterpDetailsPolynomialBasis.hpp index 42dbc995b..2efed72f2 100644 --- a/src/interpolation/details/ArborX_InterpDetailsPolynomialBasis.hpp +++ b/src/interpolation/details/ArborX_InterpDetailsPolynomialBasis.hpp @@ -107,11 +107,11 @@ KOKKOS_FUNCTION auto evaluatePolynomialBasis(Point const &p) static_assert(GeometryTraits::is_point::value, "point must be a point"); static constexpr std::size_t DIM = GeometryTraits::dimension_v; - using value_t = typename GeometryTraits::coordinate_type::type; + using Value = typename GeometryTraits::coordinate_type::type; static_assert(DIM > 0, "Polynomial basis with no dimension is invalid"); - Kokkos::Array()> arr{}; - arr[0] = value_t(1); + Kokkos::Array()> arr{}; + arr[0] = Value(1); if constexpr (Degree > 0) { diff --git a/src/interpolation/details/ArborX_InterpDetailsSymmetricPseudoInverseSVD.hpp b/src/interpolation/details/ArborX_InterpDetailsSymmetricPseudoInverseSVD.hpp index 47a1b0cd0..029408001 100644 --- a/src/interpolation/details/ArborX_InterpDetailsSymmetricPseudoInverseSVD.hpp +++ b/src/interpolation/details/ArborX_InterpDetailsSymmetricPseudoInverseSVD.hpp @@ -52,11 +52,11 @@ template KOKKOS_FUNCTION auto argmaxUpperTriangle(Matrix const &mat) { ensureIsSquareMatrix(mat); - using value_t = typename Matrix::non_const_value_type; + using Value = typename Matrix::non_const_value_type; struct { - value_t max = 0; + Value max = 0; int row = 0; int col = 0; } result; @@ -65,7 +65,7 @@ KOKKOS_FUNCTION auto argmaxUpperTriangle(Matrix const &mat) for (int i = 0; i < size; i++) for (int j = i + 1; j < size; j++) { - value_t val = Kokkos::abs(mat(i, j)); + Value val = Kokkos::abs(mat(i, j)); if (result.max < val) { result.max = val; @@ -101,18 +101,18 @@ symmetricPseudoInverseSVDSerialKernel(AMatrix &A, ESMatrix &ES, UMatrix &U) typename UMatrix::value_type>, "All input matrices must have the same value type"); KOKKOS_ASSERT(A.extent(0) == ES.extent(0) && ES.extent(0) == U.extent(0)); - using value_t = typename AMatrix::non_const_value_type; + using Value = typename AMatrix::non_const_value_type; int const size = A.extent(0); // We first initialize U as the identity matrix and copy A to ES for (int i = 0; i < size; i++) for (int j = 0; j < size; j++) { - U(i, j) = value_t(i == j); + U(i, j) = Value(i == j); ES(i, j) = A(i, j); } - static constexpr value_t epsilon = Kokkos::Experimental::epsilon_v; + static constexpr Value epsilon = Kokkos::Experimental::epsilon_v; while (true) { // We have a guarantee that p < q @@ -138,13 +138,13 @@ symmetricPseudoInverseSVDSerialKernel(AMatrix &A, ESMatrix &ES, UMatrix &U) // | b | c | | 0 | y | // +---+---+ +---+---+ - value_t cos_theta; - value_t sin_theta; - value_t x; - value_t y; + Value cos_theta; + Value sin_theta; + Value x; + Value y; if (a == c) { - cos_theta = Kokkos::sqrt(value_t(2)) / 2; + cos_theta = Kokkos::sqrt(Value(2)) / 2; sin_theta = cos_theta; x = a + b; y = a - b; @@ -212,7 +212,7 @@ symmetricPseudoInverseSVDSerialKernel(AMatrix &A, ESMatrix &ES, UMatrix &U) for (int i = 0; i < size; i++) for (int j = 0; j < size; j++) { - value_t tmp = 0; + Value tmp = 0; for (int k = 0; k < size; k++) tmp += ES(k, k) * U(i, k) * U(j, k); A(i, j) = tmp; diff --git a/test/tstInterpDetailsMLSCoefficients.cpp b/test/tstInterpDetailsMLSCoefficients.cpp index ee2ff5473..015b5bac2 100644 --- a/test/tstInterpDetailsMLSCoefficients.cpp +++ b/test/tstInterpDetailsMLSCoefficients.cpp @@ -51,9 +51,9 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients, DeviceType, ARBORX_DEVICE_TYPES) // -------0---------------> // SRC: 0 2 4 6 // TGT: 1 3 5 - using point0 = ArborX::ExperimentalHyperGeometry::Point<1, double>; - Kokkos::View srcp0("Testing::srcp0", 3, 2); - Kokkos::View tgtp0("Testing::tgtp0", 3); + using Point0 = ArborX::ExperimentalHyperGeometry::Point<1, double>; + Kokkos::View srcp0("Testing::srcp0", 3, 2); + Kokkos::View tgtp0("Testing::tgtp0", 3); Kokkos::View srcv0("Testing::srcv0", 3, 2); Kokkos::View tgtv0("Testing::tgtv0", 3); Kokkos::parallel_for( @@ -64,7 +64,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients, DeviceType, ARBORX_DEVICE_TYPES) srcp0(i, 1) = {{2. * i + 2}}; tgtp0(i) = {{2. * i + 1}}; - auto f = [](const point0 &) { return 3.; }; + auto f = [](const Point0 &) { return 3.; }; srcv0(i, 0) = f(srcp0(i, 0)); srcv0(i, 1) = f(srcp0(i, 1)); @@ -86,9 +86,9 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients, DeviceType, ARBORX_DEVICE_TYPES) // T | T // S S S // | - using point1 = ArborX::ExperimentalHyperGeometry::Point<2, double>; - Kokkos::View srcp1("Testing::srcp1", 4, 8); - Kokkos::View tgtp1("Testing::tgtp1", 4); + using Point1 = ArborX::ExperimentalHyperGeometry::Point<2, double>; + Kokkos::View srcp1("Testing::srcp1", 4, 8); + Kokkos::View tgtp1("Testing::tgtp1", 4); Kokkos::View srcv1("Testing::srcv1", 4, 8); Kokkos::View tgtv1("Testing::tgtv1", 4); Kokkos::parallel_for( @@ -109,7 +109,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients, DeviceType, ARBORX_DEVICE_TYPES) } tgtp1(i) = {{double(u), double(v)}}; - auto f = [](const point1 &p) { return p[0] * p[1] + 4 * p[0]; }; + auto f = [](const Point1 &p) { return p[0] * p[1] + 4 * p[0]; }; for (int j = 0; j < 8; j++) srcv1(i, j) = f(srcp1(i, j)); @@ -131,9 +131,9 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients_edge_cases, DeviceType, ExecutionSpace space{}; // Case 1: Same as previous case 1, but points are 2D and locked on y=0 - using point0 = ArborX::ExperimentalHyperGeometry::Point<2, double>; - Kokkos::View srcp0("Testing::srcp0", 3, 2); - Kokkos::View tgtp0("Testing::tgtp0", 3); + using Point0 = ArborX::ExperimentalHyperGeometry::Point<2, double>; + Kokkos::View srcp0("Testing::srcp0", 3, 2); + Kokkos::View tgtp0("Testing::tgtp0", 3); Kokkos::View srcv0("Testing::srcv0", 3, 2); Kokkos::View tgtv0("Testing::tgtv0", 3); Kokkos::parallel_for( @@ -144,7 +144,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients_edge_cases, DeviceType, srcp0(i, 1) = {{2. * i + 2, 0.}}; tgtp0(i) = {{2. * i + 1, 0.}}; - auto f = [](const point0 &) { return 3.; }; + auto f = [](const Point0 &) { return 3.; }; srcv0(i, 0) = f(srcp0(i, 0)); srcv0(i, 1) = f(srcp0(i, 1)); @@ -158,9 +158,9 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients_edge_cases, DeviceType, ARBORX_MDVIEW_TEST_TOL(eval0, tgtv0, Kokkos::Experimental::epsilon_v); // Case 2: Same but corner source points are also targets - using point1 = ArborX::ExperimentalHyperGeometry::Point<2, double>; - Kokkos::View srcp1("Testing::srcp1", 4, 8); - Kokkos::View tgtp1("Testing::tgtp1", 4); + using Point1 = ArborX::ExperimentalHyperGeometry::Point<2, double>; + Kokkos::View srcp1("Testing::srcp1", 4, 8); + Kokkos::View tgtp1("Testing::tgtp1", 4); Kokkos::View srcv1("Testing::srcv1", 4, 8); Kokkos::View tgtv1("Testing::tgtv1", 4); Kokkos::parallel_for( @@ -181,7 +181,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients_edge_cases, DeviceType, } tgtp1(i) = {{u * 2., v * 2.}}; - auto f = [](const point1 &p) { return p[0] * p[1] + 4 * p[0]; }; + auto f = [](const Point1 &p) { return p[0] * p[1] + 4 * p[0]; }; for (int j = 0; j < 8; j++) srcv1(i, j) = f(srcp1(i, j)); @@ -193,4 +193,4 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients_edge_cases, DeviceType, space, srcp1, tgtp1); auto eval1 = interpolate(space, srcv1, coeffs1); ARBORX_MDVIEW_TEST_TOL(eval1, tgtv1, Kokkos::Experimental::epsilon_v); -} \ No newline at end of file +} diff --git a/test/tstInterpDetailsPolyBasis.cpp b/test/tstInterpDetailsPolyBasis.cpp index 4aba902c1..22f9499f6 100644 --- a/test/tstInterpDetailsPolyBasis.cpp +++ b/test/tstInterpDetailsPolyBasis.cpp @@ -17,21 +17,21 @@ BOOST_AUTO_TEST_CASE(polynomial_basis_slice_lengths) { - using view = Kokkos::View; + using View = Kokkos::View; auto [arr0] = ArborX::Interpolation::Details::polynomialBasisSliceLengths<5, 3>(); std::size_t ref0[3][5] = { {1, 1, 1, 1, 1}, {1, 2, 3, 4, 5}, {1, 3, 6, 10, 15}}; - view arr0_view(&arr0[0][0], 3, 5); - view ref0_view(&ref0[0][0], 3, 5); + View arr0_view(&arr0[0][0], 3, 5); + View ref0_view(&ref0[0][0], 3, 5); ARBORX_MDVIEW_TEST(arr0_view, ref0_view); auto [arr1] = ArborX::Interpolation::Details::polynomialBasisSliceLengths<2, 3>(); std::size_t ref1[3][2] = {{1, 1}, {1, 2}, {1, 3}}; - view arr1_view(&arr1[0][0], 3, 2); - view ref1_view(&ref1[0][0], 3, 2); + View arr1_view(&arr1[0][0], 3, 2); + View ref1_view(&ref1[0][0], 3, 2); ARBORX_MDVIEW_TEST(arr1_view, ref1_view); } @@ -47,7 +47,7 @@ BOOST_AUTO_TEST_CASE(polynomial_basis_size) BOOST_AUTO_TEST_CASE(polynomial_basis) { - using view = Kokkos::View; + using View = Kokkos::View; ArborX::ExperimentalHyperGeometry::Point<5, double> point0 = {1, 2, 3, 4, 5}; auto arr0 = @@ -56,15 +56,15 @@ BOOST_AUTO_TEST_CASE(polynomial_basis) 12, 16, 5, 10, 15, 20, 25, 1, 2, 4, 8, 3, 6, 12, 9, 18, 27, 4, 8, 16, 12, 24, 36, 16, 32, 48, 64, 5, 10, 20, 15, 30, 45, 20, 40, 60, 80, 25, 50, 75, 100, 125}; - view arr0_view(arr0.data(), 56); - view ref0_view(&ref0[0], 56); + View arr0_view(arr0.data(), 56); + View ref0_view(&ref0[0], 56); ARBORX_MDVIEW_TEST(arr0_view, ref0_view); ArborX::ExperimentalHyperGeometry::Point<2, double> point1 = {-2, 0}; auto arr1 = ArborX::Interpolation::Details::evaluatePolynomialBasis<3>(point1); double ref1[10] = {1, -2, 0, 4, 0, 0, -8, 0, 0, 0}; - view arr1_view(arr1.data(), 10); - view ref1_view(&ref1[0], 10); + View arr1_view(arr1.data(), 10); + View ref1_view(&ref1[0], 10); ARBORX_MDVIEW_TEST(arr1_view, ref1_view); -} \ No newline at end of file +} diff --git a/test/tstInterpDetailsSVD.cpp b/test/tstInterpDetailsSVD.cpp index f0f56af1f..df8b26033 100644 --- a/test/tstInterpDetailsSVD.cpp +++ b/test/tstInterpDetailsSVD.cpp @@ -20,12 +20,12 @@ template void makeCase(ES const &es, V const (&src_arr)[M][N][N], V const (&ref_arr)[M][N][N]) { - using device_view = Kokkos::View; - using host_view = typename device_view::HostMirror; + using DeviceView = Kokkos::View; + using HostView = typename DeviceView::HostMirror; - host_view src("Testing::src"); - host_view ref("Testing::ref"); - device_view inv("Testing::inv"); + HostView src("Testing::src"); + HostView ref("Testing::ref"); + DeviceView inv("Testing::inv"); for (int i = 0; i < M; i++) for (int j = 0; j < N; j++) @@ -125,4 +125,4 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(pseudo_inv_empty, DeviceType, ARBORX_DEVICE_TYPES) Kokkos::View mat("mat", 0, 0, 0); ArborX::Interpolation::Details::symmetricPseudoInverseSVD(space, mat); BOOST_TEST(mat.size() == 0); -} \ No newline at end of file +} diff --git a/test/tstInterpMovingLeastSquares.cpp b/test/tstInterpMovingLeastSquares.cpp index ff9e4a767..9f6a786a4 100644 --- a/test/tstInterpMovingLeastSquares.cpp +++ b/test/tstInterpMovingLeastSquares.cpp @@ -29,9 +29,9 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares, DeviceType, // -------0---------------> // SRC: 0 2 4 6 // TGT: 1 3 5 - using point0 = ArborX::ExperimentalHyperGeometry::Point<1, double>; - Kokkos::View srcp0("Testing::srcp0", 4); - Kokkos::View tgtp0("Testing::tgtp0", 3); + using Point0 = ArborX::ExperimentalHyperGeometry::Point<1, double>; + Kokkos::View srcp0("Testing::srcp0", 4); + Kokkos::View tgtp0("Testing::tgtp0", 3); Kokkos::View srcv0("Testing::srcv0", 4); Kokkos::View tgtv0("Testing::tgtv0", 3); Kokkos::View eval0("Testing::eval0", 0); @@ -39,7 +39,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares, DeviceType, "Testing::moving_least_squares::for0", Kokkos::RangePolicy(space, 0, 4), KOKKOS_LAMBDA(int const i) { - auto f = [](const point0 &) { return 3.; }; + auto f = [](const Point0 &) { return 3.; }; srcp0(i) = {{2. * i}}; srcv0(i) = f(srcp0(i)); @@ -65,9 +65,9 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares, DeviceType, // T | T // S S S // | - using point1 = ArborX::ExperimentalHyperGeometry::Point<2, double>; - Kokkos::View srcp1("Testing::srcp1", 9); - Kokkos::View tgtp1("Testing::tgtp1", 4); + using Point1 = ArborX::ExperimentalHyperGeometry::Point<2, double>; + Kokkos::View srcp1("Testing::srcp1", 9); + Kokkos::View tgtp1("Testing::tgtp1", 4); Kokkos::View srcv1("Testing::srcv1", 9); Kokkos::View tgtv1("Testing::tgtv1", 4); Kokkos::View eval1("Testing::eval1", 0); @@ -79,7 +79,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares, DeviceType, int v = (i % 2) * 2 - 1; int x = (i / 3) - 1; int y = (i % 3) - 1; - auto f = [](const point1 &p) { return p[0] * p[1] + 4 * p[0]; }; + auto f = [](const Point1 &p) { return p[0] * p[1] + 4 * p[0]; }; srcp1(i) = {{x * 2., y * 2.}}; srcv1(i) = f(srcp1(i)); @@ -104,9 +104,9 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares_edge_cases, DeviceType, ExecutionSpace space{}; // Case 1: Same as previous case 1, but points are 2D and locked on y=0 - using point0 = ArborX::ExperimentalHyperGeometry::Point<2, double>; - Kokkos::View srcp0("Testing::srcp0", 4); - Kokkos::View tgtp0("Testing::tgtp0", 3); + using Point0 = ArborX::ExperimentalHyperGeometry::Point<2, double>; + Kokkos::View srcp0("Testing::srcp0", 4); + Kokkos::View tgtp0("Testing::tgtp0", 3); Kokkos::View srcv0("Testing::srcv0", 4); Kokkos::View tgtv0("Testing::tgtv0", 3); Kokkos::View eval0("Testing::eval0", 0); @@ -114,7 +114,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares_edge_cases, DeviceType, "Testing::moving_least_squares_edge_cases::for0", Kokkos::RangePolicy(space, 0, 4), KOKKOS_LAMBDA(int const i) { - auto f = [](const point0 &) { return 3.; }; + auto f = [](const Point0 &) { return 3.; }; srcp0(i) = {{2. * i, 0.}}; srcv0(i) = f(srcp0(i)); @@ -132,9 +132,9 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares_edge_cases, DeviceType, ARBORX_MDVIEW_TEST_TOL(eval0, tgtv0, Kokkos::Experimental::epsilon_v); // Case 2: Same but corner source points are also targets - using point1 = ArborX::ExperimentalHyperGeometry::Point<2, double>; - Kokkos::View srcp1("Testing::srcp1", 9); - Kokkos::View tgtp1("Testing::tgtp1", 4); + using Point1 = ArborX::ExperimentalHyperGeometry::Point<2, double>; + Kokkos::View srcp1("Testing::srcp1", 9); + Kokkos::View tgtp1("Testing::tgtp1", 4); Kokkos::View srcv1("Testing::srcv1", 9); Kokkos::View tgtv1("Testing::tgtv1", 4); Kokkos::View eval1("Testing::eval1", 0); @@ -146,7 +146,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares_edge_cases, DeviceType, int v = (i % 2) * 2 - 1; int x = (i / 3) - 1; int y = (i % 3) - 1; - auto f = [](const point1 &p) { return p[0] * p[1] + 4 * p[0]; }; + auto f = [](const Point1 &p) { return p[0] * p[1] + 4 * p[0]; }; srcp1(i) = {{x * 2., y * 2.}}; srcv1(i) = f(srcp1(i)); @@ -161,4 +161,4 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares_edge_cases, DeviceType, ArborX::Interpolation::PolynomialDegree<2>{}, 8); mls1.interpolate(space, srcv1, eval1); ARBORX_MDVIEW_TEST_TOL(eval1, tgtv1, Kokkos::Experimental::epsilon_v); -} \ No newline at end of file +} From 5b86976e8add5b037adf5420de7a7744cf838715 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Thu, 21 Dec 2023 10:07:21 -0500 Subject: [PATCH 05/15] Template CRBFs on input (point) type --- .../ArborX_InterpMovingLeastSquares.hpp | 6 ++--- ...nterpDetailsCompactRadialBasisFunction.hpp | 25 ++++++++++++++++--- ...pDetailsMovingLeastSquaresCoefficients.hpp | 11 ++++---- test/tstInterpDetailsCRBF.cpp | 8 +++++- 4 files changed, 36 insertions(+), 14 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index 1d384a1c3..f37c40a7b 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -75,11 +75,11 @@ class MovingLeastSquares { public: template , + typename TargetPoints, typename CRBFunc = CRBF::Wendland<0>, typename PolynomialDegree = PolynomialDegree<2>> MovingLeastSquares(ExecutionSpace const &space, SourcePoints const &source_points, - TargetPoints const &target_points, CRBF = {}, + TargetPoints const &target_points, CRBFunc = {}, PolynomialDegree = {}, std::optional num_neighbors = std::nullopt) { @@ -156,7 +156,7 @@ class MovingLeastSquares // Compute the moving least squares coefficients _coeffs = Details::movingLeastSquaresCoefficients< - CRBF, PolynomialDegree, FloatingCalculationType, MemorySpace>( + CRBFunc, PolynomialDegree, FloatingCalculationType, MemorySpace>( space, source_view, target_points); } diff --git a/src/interpolation/details/ArborX_InterpDetailsCompactRadialBasisFunction.hpp b/src/interpolation/details/ArborX_InterpDetailsCompactRadialBasisFunction.hpp index 23d12bd14..019cb29e3 100644 --- a/src/interpolation/details/ArborX_InterpDetailsCompactRadialBasisFunction.hpp +++ b/src/interpolation/details/ArborX_InterpDetailsCompactRadialBasisFunction.hpp @@ -12,6 +12,9 @@ #ifndef ARBORX_INTERP_DETAILS_COMPACT_RADIAL_BASIS_FUNCTION_HPP #define ARBORX_INTERP_DETAILS_COMPACT_RADIAL_BASIS_FUNCTION_HPP +#include +#include + #include #include @@ -40,20 +43,22 @@ namespace CRBF { #define CRBF_DECL(NAME) \ - template \ + template \ struct NAME; #define CRBF_DEF(NAME, N, FUNC) \ template <> \ struct NAME \ { \ - template \ + template \ KOKKOS_INLINE_FUNCTION static constexpr T evaluate(T const y) \ { \ /* We force the input to be between 0 and 1. \ Because CRBF(-a) = CRBF(a) = CRBF(|a|), we take the absolute value \ and clamp the range to [0, 1] before entering in the definition of \ - the CRBF. */ \ + the CRBF. \ + We also template the internal function on the dimension as CRBFs \ + depend on the point's dimensionality. */ \ T const x = Kokkos::min(Kokkos::abs(y), T(1)); \ return Kokkos::abs(FUNC); \ } \ @@ -90,8 +95,20 @@ CRBF_DEF(Buhmann, 4, #undef CRBF_DEF #undef CRBF_DECL +template +KOKKOS_INLINE_FUNCTION constexpr auto +evaluate(Point const &point, + typename GeometryTraits::coordinate_type::type const radius) +{ + static_assert(GeometryTraits::is_point::value, + "point must be a point"); + constexpr std::size_t dim = GeometryTraits::dimension_v; + return CRBFunc::template evaluate( + ArborX::Details::distance(point, Point{}) / radius); +} + } // namespace CRBF } // namespace ArborX::Interpolation -#endif \ No newline at end of file +#endif diff --git a/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp b/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp index 6dfc1276d..823b30e62 100644 --- a/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp +++ b/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp @@ -16,6 +16,7 @@ #include #include #include +#include #include #include @@ -25,9 +26,9 @@ namespace ArborX::Interpolation::Details { -template +template Kokkos::View movingLeastSquaresCoefficients(ExecutionSpace const &space, SourcePoints const &source_points, @@ -157,9 +158,7 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, Kokkos::MDRangePolicy>( space, {0, 0}, {num_targets, num_neighbors}), KOKKOS_LAMBDA(int const i, int const j) { - CoefficientsType norm = - ArborX::Details::distance(source_ref_target(i, j), Point{}); - phi(i, j) = CRBF::evaluate(norm / radii(i)); + phi(i, j) = CRBF::evaluate(source_ref_target(i, j), radii(i)); }); Kokkos::Profiling::popRegion(); diff --git a/test/tstInterpDetailsCRBF.cpp b/test/tstInterpDetailsCRBF.cpp index 39d1ed457..54b9fb4ee 100644 --- a/test/tstInterpDetailsCRBF.cpp +++ b/test/tstInterpDetailsCRBF.cpp @@ -11,6 +11,8 @@ #include "ArborX_EnableDeviceTypes.hpp" #include "ArborX_EnableViewComparison.hpp" +#include +#include #include #include "BoostTest_CUDA_clang_workarounds.hpp" @@ -20,6 +22,7 @@ template void makeCase(ES const &es, std::function const &tf, T tol = 1e-5) { + using Point = ArborX::ExperimentalHyperGeometry::Point<1, T>; using View = Kokkos::View; using HostView = typename View::HostMirror; static constexpr int range = 15; @@ -37,7 +40,10 @@ void makeCase(ES const &es, std::function const &tf, T tol = 1e-5) Kokkos::deep_copy(es, eval, input); Kokkos::parallel_for( "Testing::eval_crbf", Kokkos::RangePolicy(es, 0, 4 * range), - KOKKOS_LAMBDA(int const i) { eval(i) = CRBF::evaluate(eval(i)); }); + KOKKOS_LAMBDA(int const i) { + eval(i) = + ArborX::Interpolation::CRBF::evaluate(Point{eval(i)}, 1); + }); if (bool(tf)) { From 27dd128213d75dcda878841c8707e1bc0aadde1f Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Thu, 21 Dec 2023 10:14:21 -0500 Subject: [PATCH 06/15] Removing reallocation for `approx_values` --- src/interpolation/ArborX_InterpMovingLeastSquares.hpp | 3 ++- test/tstInterpMovingLeastSquares.cpp | 8 ++++---- 2 files changed, 6 insertions(+), 5 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index f37c40a7b..2c4340f07 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -241,7 +241,8 @@ class MovingLeastSquares int const num_targets = _values_indices.extent(0); int const num_neighbors = _values_indices.extent(1); - KokkosExt::reallocWithoutInitializing(space, approx_values, num_targets); + // Approx values must have the correct size + KOKKOS_ASSERT(approx_values.extent_int(0) == num_targets); Kokkos::parallel_for( "ArborX::MovingLeastSquares::target_interpolation", diff --git a/test/tstInterpMovingLeastSquares.cpp b/test/tstInterpMovingLeastSquares.cpp index 9f6a786a4..ceee3e7b4 100644 --- a/test/tstInterpMovingLeastSquares.cpp +++ b/test/tstInterpMovingLeastSquares.cpp @@ -34,7 +34,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares, DeviceType, Kokkos::View tgtp0("Testing::tgtp0", 3); Kokkos::View srcv0("Testing::srcv0", 4); Kokkos::View tgtv0("Testing::tgtv0", 3); - Kokkos::View eval0("Testing::eval0", 0); + Kokkos::View eval0("Testing::eval0", 3); Kokkos::parallel_for( "Testing::moving_least_squares::for0", Kokkos::RangePolicy(space, 0, 4), @@ -70,7 +70,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares, DeviceType, Kokkos::View tgtp1("Testing::tgtp1", 4); Kokkos::View srcv1("Testing::srcv1", 9); Kokkos::View tgtv1("Testing::tgtv1", 4); - Kokkos::View eval1("Testing::eval1", 0); + Kokkos::View eval1("Testing::eval1", 4); Kokkos::parallel_for( "Testing::moving_least_squares::for1", Kokkos::RangePolicy(space, 0, 9), @@ -109,7 +109,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares_edge_cases, DeviceType, Kokkos::View tgtp0("Testing::tgtp0", 3); Kokkos::View srcv0("Testing::srcv0", 4); Kokkos::View tgtv0("Testing::tgtv0", 3); - Kokkos::View eval0("Testing::eval0", 0); + Kokkos::View eval0("Testing::eval0", 3); Kokkos::parallel_for( "Testing::moving_least_squares_edge_cases::for0", Kokkos::RangePolicy(space, 0, 4), @@ -137,7 +137,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(moving_least_squares_edge_cases, DeviceType, Kokkos::View tgtp1("Testing::tgtp1", 4); Kokkos::View srcv1("Testing::srcv1", 9); Kokkos::View tgtv1("Testing::tgtv1", 4); - Kokkos::View eval1("Testing::eval1", 0); + Kokkos::View eval1("Testing::eval1", 4); Kokkos::parallel_for( "Testing::moving_least_squares_edge_cases::for0", Kokkos::RangePolicy(space, 0, 9), From fc7a2ed083955f4d9e6427af9c374b0bb954f660 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Thu, 21 Dec 2023 10:26:31 -0500 Subject: [PATCH 07/15] Using raw indices and avoiding remapping --- .../ArborX_InterpMovingLeastSquares.hpp | 55 ++++++++----------- 1 file changed, 22 insertions(+), 33 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index 2c4340f07..8b40d6352 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -120,17 +120,17 @@ class MovingLeastSquares static_assert(dimension == GeometryTraits::dimension_v, "Target and source points must have the same dimension"); - int num_neighbors_val = + _num_neighbors = num_neighbors ? *num_neighbors : Details::polynomialBasisSize(); TargetAccess target_access{target_points}; SourceAccess source_access{source_points}; - int const num_targets = target_access.size(); + _num_targets = target_access.size(); _source_size = source_access.size(); // There must be enough source points - KOKKOS_ASSERT(0 < num_neighbors_val && num_neighbors_val <= _source_size); + KOKKOS_ASSERT(0 < _num_neighbors && _num_neighbors <= _source_size); // Organize the source points as a tree BoundingVolumeHierarchy> @@ -138,21 +138,20 @@ class MovingLeastSquares // Create the predicates Details::MLSTargetPointsPredicateWrapper predicates{ - {target_points}, num_neighbors_val}; + {target_points}, _num_neighbors}; // Query the source - Kokkos::View indices( + _indices = Kokkos::View( "ArborX::MovingLeastSquares::indices", 0); Kokkos::View offsets( "ArborX::MovingLeastSquares::offsets", 0); source_tree.query(space, predicates, - ArborX::Details::LegacyDefaultCallback{}, indices, + ArborX::Details::LegacyDefaultCallback{}, _indices, offsets); // Fill in the value indices object so values can be transferred from a 1D // source data to a properly distributed 2D array for each target. - auto const source_view = fillValuesIndicesAndGetSourceView( - space, indices, offsets, num_targets, num_neighbors_val, source_points); + auto const source_view = getSourceView(space, source_points); // Compute the moving least squares coefficients _coeffs = Details::movingLeastSquaresCoefficients< @@ -164,36 +163,26 @@ class MovingLeastSquares Kokkos::View>::type **, MemorySpace> - fillValuesIndicesAndGetSourceView( - ExecutionSpace const &space, - Kokkos::View const &indices, - Kokkos::View const &offsets, int const num_targets, - int const num_neighbors, SourcePoints const &source_points) + getSourceView(ExecutionSpace const &space, SourcePoints const &source_points) { auto guard = Kokkos::Profiling::ScopedRegion( - "ArborX::MovingLeastSquares::fillValuesIndicesAndGetSourceView"); + "ArborX::MovingLeastSquares::getSourceView"); using SourceAccess = ArborX::Details::AccessValues; using SourcePoint = typename SourceAccess::value_type; - _values_indices = Kokkos::View( - Kokkos::view_alloc(Kokkos::WithoutInitializing, - "ArborX::MovingLeastSquares::values_indices"), - num_targets, num_neighbors); Kokkos::View source_view( Kokkos::view_alloc(Kokkos::WithoutInitializing, "ArborX::MovingLeastSquares::source_view"), - num_targets, num_neighbors); + _num_targets, _num_neighbors); SourceAccess source_access{source_points}; Kokkos::parallel_for( "ArborX::MovingLeastSquares::values_indices_and_source_view_fill", Kokkos::MDRangePolicy>( - space, {0, 0}, {num_targets, num_neighbors}), + space, {0, 0}, {_num_targets, _num_neighbors}), KOKKOS_CLASS_LAMBDA(int const i, int const j) { - auto index = indices(offsets(i) + j); - _values_indices(i, j) = index; - source_view(i, j) = source_access(index); + source_view(i, j) = source_access(_indices(i * _num_neighbors + j)); }); return source_view; @@ -236,28 +225,28 @@ class MovingLeastSquares // original input KOKKOS_ASSERT(_source_size == source_values.extent_int(0)); - using Value = typename ApproxValues::non_const_value_type; - - int const num_targets = _values_indices.extent(0); - int const num_neighbors = _values_indices.extent(1); - // Approx values must have the correct size - KOKKOS_ASSERT(approx_values.extent_int(0) == num_targets); + KOKKOS_ASSERT(approx_values.extent_int(0) == _num_targets); + + using Value = typename ApproxValues::non_const_value_type; Kokkos::parallel_for( "ArborX::MovingLeastSquares::target_interpolation", - Kokkos::RangePolicy(space, 0, num_targets), + Kokkos::RangePolicy(space, 0, _num_targets), KOKKOS_CLASS_LAMBDA(int const i) { Value tmp = 0; - for (int j = 0; j < num_neighbors; j++) - tmp += _coeffs(i, j) * source_values(_values_indices(i, j)); + for (int j = 0; j < _num_neighbors; j++) + tmp += + _coeffs(i, j) * source_values(_indices(i * _num_neighbors + j)); approx_values(i) = tmp; }); } private: Kokkos::View _coeffs; - Kokkos::View _values_indices; + Kokkos::View _indices; + int _num_targets; + int _num_neighbors; int _source_size; }; From 2f814b5569cdd5c6201ea7638d37f270ed1c6926 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Thu, 21 Dec 2023 10:55:42 -0500 Subject: [PATCH 08/15] Private get source view --- .../ArborX_InterpMovingLeastSquares.hpp | 73 +++++++++++-------- 1 file changed, 44 insertions(+), 29 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index 8b40d6352..c6d9ff1af 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -38,6 +38,21 @@ struct MLSTargetPointsPredicateWrapper int num_neighbors; }; +// Private kernel to compute the 2D source view +template +struct GetSourceViewKernel +{ + SourceView source_view; + ArborX::Details::AccessValues source_access; + Indices indices; + int num_neighbors; + + KOKKOS_FUNCTION void operator()(int const i, int const j) const + { + source_view(i, j) = source_access(indices(i * num_neighbors + j)); + } +}; + } // namespace ArborX::Interpolation::Details namespace ArborX @@ -159,35 +174,6 @@ class MovingLeastSquares space, source_view, target_points); } - template - Kokkos::View>::type **, - MemorySpace> - getSourceView(ExecutionSpace const &space, SourcePoints const &source_points) - { - auto guard = Kokkos::Profiling::ScopedRegion( - "ArborX::MovingLeastSquares::getSourceView"); - - using SourceAccess = - ArborX::Details::AccessValues; - using SourcePoint = typename SourceAccess::value_type; - - Kokkos::View source_view( - Kokkos::view_alloc(Kokkos::WithoutInitializing, - "ArborX::MovingLeastSquares::source_view"), - _num_targets, _num_neighbors); - SourceAccess source_access{source_points}; - Kokkos::parallel_for( - "ArborX::MovingLeastSquares::values_indices_and_source_view_fill", - Kokkos::MDRangePolicy>( - space, {0, 0}, {_num_targets, _num_neighbors}), - KOKKOS_CLASS_LAMBDA(int const i, int const j) { - source_view(i, j) = source_access(_indices(i * _num_neighbors + j)); - }); - - return source_view; - } - template void interpolate(ExecutionSpace const &space, @@ -243,6 +229,35 @@ class MovingLeastSquares } private: + template + auto getSourceView(ExecutionSpace const &space, + SourcePoints const &source_points) + { + auto guard = Kokkos::Profiling::ScopedRegion( + "ArborX::MovingLeastSquares::getSourceView"); + + using SourceAccess = + ArborX::Details::AccessValues; + using SourcePoint = typename SourceAccess::value_type; + + Kokkos::View source_view( + Kokkos::view_alloc(Kokkos::WithoutInitializing, + "ArborX::MovingLeastSquares::source_view"), + _num_targets, _num_neighbors); + + Details::GetSourceViewKernel + kernel{source_view, {source_points}, _indices, _num_neighbors}; + + Kokkos::parallel_for( + "ArborX::MovingLeastSquares::values_indices_and_source_view_fill", + Kokkos::MDRangePolicy>( + space, {0, 0}, {_num_targets, _num_neighbors}), + kernel); + + return source_view; + } + Kokkos::View _coeffs; Kokkos::View _indices; int _num_targets; From fcaefb0bc3fcf8791a710cd64912b0a88c2abbf8 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Thu, 21 Dec 2023 11:41:27 -0500 Subject: [PATCH 09/15] Handling of tree destruction and unused views --- .../ArborX_InterpMovingLeastSquares.hpp | 48 ++++++++++++------- 1 file changed, 32 insertions(+), 16 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index c6d9ff1af..6ab19429c 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -147,22 +147,8 @@ class MovingLeastSquares // There must be enough source points KOKKOS_ASSERT(0 < _num_neighbors && _num_neighbors <= _source_size); - // Organize the source points as a tree - BoundingVolumeHierarchy> - source_tree(space, ArborX::Experimental::attach_indices(source_points)); - - // Create the predicates - Details::MLSTargetPointsPredicateWrapper predicates{ - {target_points}, _num_neighbors}; - - // Query the source - _indices = Kokkos::View( - "ArborX::MovingLeastSquares::indices", 0); - Kokkos::View offsets( - "ArborX::MovingLeastSquares::offsets", 0); - source_tree.query(space, predicates, - ArborX::Details::LegacyDefaultCallback{}, _indices, - offsets); + // Search for neighbors + searchNeighbors(space, source_points, target_points); // Fill in the value indices object so values can be transferred from a 1D // source data to a properly distributed 2D array for each target. @@ -258,6 +244,36 @@ class MovingLeastSquares return source_view; } + template + void searchNeighbors(ExecutionSpace const &space, + SourcePoints const &source_points, + TargetPoints const &target_points) + { + auto guard = Kokkos::Profiling::ScopedRegion( + "ArborX::MovingLeastSquares::searchNeighbors"); + + // Organize the source points as a tree + using SourcePoint = + typename ArborX::Details::AccessValues::value_type; + BoundingVolumeHierarchy> + source_tree(space, ArborX::Experimental::attach_indices(source_points)); + + // Create the predicates + Details::MLSTargetPointsPredicateWrapper predicates{ + {target_points}, _num_neighbors}; + + // Query the source + _indices = Kokkos::View( + "ArborX::MovingLeastSquares::indices", 0); + Kokkos::View offsets( + "ArborX::MovingLeastSquares::offsets", 0); + source_tree.query(space, predicates, + ArborX::Details::LegacyDefaultCallback{}, _indices, + offsets); + } + Kokkos::View _coeffs; Kokkos::View _indices; int _num_targets; From c38f97fdd3031dabae36e0f797adf0fc81c2f051 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Thu, 21 Dec 2023 13:47:41 -0500 Subject: [PATCH 10/15] Deriving MemorySpace from source points --- .../ArborX_InterpMovingLeastSquares.hpp | 4 ++-- ...X_InterpDetailsMovingLeastSquaresCoefficients.hpp | 7 ++++--- test/tstInterpDetailsMLSCoefficients.cpp | 12 ++++-------- 3 files changed, 10 insertions(+), 13 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index 6ab19429c..12f90e1cb 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -155,8 +155,8 @@ class MovingLeastSquares auto const source_view = getSourceView(space, source_points); // Compute the moving least squares coefficients - _coeffs = Details::movingLeastSquaresCoefficients< - CRBFunc, PolynomialDegree, FloatingCalculationType, MemorySpace>( + _coeffs = Details::movingLeastSquaresCoefficients( space, source_view, target_points); } diff --git a/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp b/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp index 823b30e62..47152418e 100644 --- a/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp +++ b/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp @@ -27,9 +27,9 @@ namespace ArborX::Interpolation::Details { template -Kokkos::View + typename CoefficientsType, typename ExecutionSpace, + typename SourcePoints, typename TargetPoints> +Kokkos::View movingLeastSquaresCoefficients(ExecutionSpace const &space, SourcePoints const &source_points, TargetPoints const &target_points) @@ -39,6 +39,7 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, namespace KokkosExt = ::ArborX::Details::KokkosExt; + using MemorySpace = typename SourcePoints::memory_space; static_assert( KokkosExt::is_accessible_from::value, "Memory space must be accessible from the execution space"); diff --git a/test/tstInterpDetailsMLSCoefficients.cpp b/test/tstInterpDetailsMLSCoefficients.cpp index 015b5bac2..c94ec70f6 100644 --- a/test/tstInterpDetailsMLSCoefficients.cpp +++ b/test/tstInterpDetailsMLSCoefficients.cpp @@ -72,8 +72,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients, DeviceType, ARBORX_DEVICE_TYPES) }); auto coeffs0 = ArborX::Interpolation::Details::movingLeastSquaresCoefficients< ArborX::Interpolation::CRBF::Wendland<0>, - ArborX::Interpolation::PolynomialDegree<1>, double, MemorySpace>( - space, srcp0, tgtp0); + ArborX::Interpolation::PolynomialDegree<1>, double>(space, srcp0, tgtp0); auto eval0 = interpolate(space, srcv0, coeffs0); ARBORX_MDVIEW_TEST_TOL(eval0, tgtv0, Kokkos::Experimental::epsilon_v); @@ -117,8 +116,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients, DeviceType, ARBORX_DEVICE_TYPES) }); auto coeffs1 = ArborX::Interpolation::Details::movingLeastSquaresCoefficients< ArborX::Interpolation::CRBF::Wendland<2>, - ArborX::Interpolation::PolynomialDegree<2>, double, MemorySpace>( - space, srcp1, tgtp1); + ArborX::Interpolation::PolynomialDegree<2>, double>(space, srcp1, tgtp1); auto eval1 = interpolate(space, srcv1, coeffs1); ARBORX_MDVIEW_TEST_TOL(eval1, tgtv1, Kokkos::Experimental::epsilon_v); } @@ -152,8 +150,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients_edge_cases, DeviceType, }); auto coeffs0 = ArborX::Interpolation::Details::movingLeastSquaresCoefficients< ArborX::Interpolation::CRBF::Wendland<0>, - ArborX::Interpolation::PolynomialDegree<1>, double, MemorySpace>( - space, srcp0, tgtp0); + ArborX::Interpolation::PolynomialDegree<1>, double>(space, srcp0, tgtp0); auto eval0 = interpolate(space, srcv0, coeffs0); ARBORX_MDVIEW_TEST_TOL(eval0, tgtv0, Kokkos::Experimental::epsilon_v); @@ -189,8 +186,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(mls_coefficients_edge_cases, DeviceType, }); auto coeffs1 = ArborX::Interpolation::Details::movingLeastSquaresCoefficients< ArborX::Interpolation::CRBF::Wendland<2>, - ArborX::Interpolation::PolynomialDegree<2>, double, MemorySpace>( - space, srcp1, tgtp1); + ArborX::Interpolation::PolynomialDegree<2>, double>(space, srcp1, tgtp1); auto eval1 = interpolate(space, srcv1, coeffs1); ARBORX_MDVIEW_TEST_TOL(eval1, tgtv1, Kokkos::Experimental::epsilon_v); } From af5193681fdb80d3492bfeef2cb14b88fde04908 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Fri, 22 Dec 2023 09:42:31 -0500 Subject: [PATCH 11/15] Renaming CRBF --- test/CMakeLists.txt | 2 +- ...RBF.cpp => tstInterpDetailsCompactRadialBasisFunction.cpp} | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) rename test/{tstInterpDetailsCRBF.cpp => tstInterpDetailsCompactRadialBasisFunction.cpp} (96%) diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 4e6cf2016..b1cfa4c18 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -253,7 +253,7 @@ add_test(NAME ArborX_Test_BoostAdapters COMMAND ArborX_Test_BoostAdapters.exe) add_executable(ArborX_Test_InterpMovingLeastSquares.exe tstInterpDetailsSVD.cpp - tstInterpDetailsCRBF.cpp + tstInterpDetailsCompactRadialBasisFunction.cpp tstInterpDetailsPolyBasis.cpp tstInterpDetailsMLSCoefficients.cpp tstInterpMovingLeastSquares.cpp diff --git a/test/tstInterpDetailsCRBF.cpp b/test/tstInterpDetailsCompactRadialBasisFunction.cpp similarity index 96% rename from test/tstInterpDetailsCRBF.cpp rename to test/tstInterpDetailsCompactRadialBasisFunction.cpp index 54b9fb4ee..b3167cea3 100644 --- a/test/tstInterpDetailsCRBF.cpp +++ b/test/tstInterpDetailsCompactRadialBasisFunction.cpp @@ -65,8 +65,8 @@ void makeCase(ES const &es, std::function const &tf, T tol = 1e-5) } #define MAKE_TEST(F, I, TF) \ - BOOST_AUTO_TEST_CASE_TEMPLATE(crbf_##F##_##I, DeviceType, \ - ARBORX_DEVICE_TYPES) \ + BOOST_AUTO_TEST_CASE_TEMPLATE(compact_radial_basis_function_##F##_##I, \ + DeviceType, ARBORX_DEVICE_TYPES) \ { \ makeCase>( \ typename DeviceType::execution_space{}, TF); \ From 8d9abeb851fe42ef82286f243676013773362c26 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Fri, 22 Dec 2023 09:51:27 -0500 Subject: [PATCH 12/15] Adding new lines at the end of files --- examples/moving_least_squares/CMakeLists.txt | 2 +- examples/moving_least_squares/moving_least_squares.cpp | 2 +- test/tstInterpDetailsCompactRadialBasisFunction.cpp | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/examples/moving_least_squares/CMakeLists.txt b/examples/moving_least_squares/CMakeLists.txt index 9f3ec8d3f..a22bdbe30 100644 --- a/examples/moving_least_squares/CMakeLists.txt +++ b/examples/moving_least_squares/CMakeLists.txt @@ -1,3 +1,3 @@ add_executable(ArborX_Example_MovingLeastSquares.exe moving_least_squares.cpp) target_link_libraries(ArborX_Example_MovingLeastSquares.exe ArborX::ArborX) -add_test(NAME ArborX_Example_MovingLeastSquares COMMAND ArborX_Example_MovingLeastSquares.exe) \ No newline at end of file +add_test(NAME ArborX_Example_MovingLeastSquares COMMAND ArborX_Example_MovingLeastSquares.exe) diff --git a/examples/moving_least_squares/moving_least_squares.cpp b/examples/moving_least_squares/moving_least_squares.cpp index ef7cf161e..dcc15cc12 100644 --- a/examples/moving_least_squares/moving_least_squares.cpp +++ b/examples/moving_least_squares/moving_least_squares.cpp @@ -110,4 +110,4 @@ int main(int argc, char *argv[]) << diff(2) << '\n'; return 0; -} \ No newline at end of file +} diff --git a/test/tstInterpDetailsCompactRadialBasisFunction.cpp b/test/tstInterpDetailsCompactRadialBasisFunction.cpp index b3167cea3..d17419bdb 100644 --- a/test/tstInterpDetailsCompactRadialBasisFunction.cpp +++ b/test/tstInterpDetailsCompactRadialBasisFunction.cpp @@ -100,4 +100,4 @@ MAKE_TEST_NONE(Buhmann, 4) #undef MAKE_TEST_NONE #undef MAKE_TEST_POLY -#undef MAKE_TEST \ No newline at end of file +#undef MAKE_TEST From 9933846e52aec34f3da873b3778fcc0db1495e0a Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Fri, 22 Dec 2023 13:59:19 -0500 Subject: [PATCH 13/15] Moving source view computation in private callback Using primitives instead of access values in callback --- .../ArborX_InterpMovingLeastSquares.hpp | 97 ++++++++----------- 1 file changed, 41 insertions(+), 56 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index 12f90e1cb..373cf3949 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -38,18 +38,26 @@ struct MLSTargetPointsPredicateWrapper int num_neighbors; }; -// Private kernel to compute the 2D source view -template -struct GetSourceViewKernel +// Functor used in the tree query to create the 2D source view and indices +template +struct SearchNeighborsCallback { SourceView source_view; - ArborX::Details::AccessValues source_access; - Indices indices; - int num_neighbors; + IndicesView indices; + CounterView counter; + + using SourcePoint = typename SourceView::non_const_value_type; - KOKKOS_FUNCTION void operator()(int const i, int const j) const + template + KOKKOS_FUNCTION void + operator()(Predicate const &predicate, + ArborX::PairValueIndex const &primitive) const { - source_view(i, j) = source_access(indices(i * num_neighbors + j)); + int const target = getData(predicate); + int const source = primitive.index; + auto count = Kokkos::atomic_fetch_add(&counter(target), 1); + indices(target, count) = source; + source_view(target, count) = primitive.value; } }; @@ -73,7 +81,7 @@ struct AccessTraits< get(Interpolation::Details::MLSTargetPointsPredicateWrapper const &tp, int const i) { - return nearest(tp.target_access(i), tp.num_neighbors); + return attach(nearest(tp.target_access(i), tp.num_neighbors), i); } using memory_space = @@ -147,12 +155,8 @@ class MovingLeastSquares // There must be enough source points KOKKOS_ASSERT(0 < _num_neighbors && _num_neighbors <= _source_size); - // Search for neighbors - searchNeighbors(space, source_points, target_points); - - // Fill in the value indices object so values can be transferred from a 1D - // source data to a properly distributed 2D array for each target. - auto const source_view = getSourceView(space, source_points); + // Search for neighbors and get the arranged source points + auto source_view = searchNeighbors(space, source_points, target_points); // Compute the moving least squares coefficients _coeffs = Details::movingLeastSquaresCoefficients - auto getSourceView(ExecutionSpace const &space, - SourcePoints const &source_points) - { - auto guard = Kokkos::Profiling::ScopedRegion( - "ArborX::MovingLeastSquares::getSourceView"); - - using SourceAccess = - ArborX::Details::AccessValues; - using SourcePoint = typename SourceAccess::value_type; - - Kokkos::View source_view( - Kokkos::view_alloc(Kokkos::WithoutInitializing, - "ArborX::MovingLeastSquares::source_view"), - _num_targets, _num_neighbors); - - Details::GetSourceViewKernel - kernel{source_view, {source_points}, _indices, _num_neighbors}; - - Kokkos::parallel_for( - "ArborX::MovingLeastSquares::values_indices_and_source_view_fill", - Kokkos::MDRangePolicy>( - space, {0, 0}, {_num_targets, _num_neighbors}), - kernel); - - return source_view; - } - template - void searchNeighbors(ExecutionSpace const &space, + auto searchNeighbors(ExecutionSpace const &space, SourcePoints const &source_points, TargetPoints const &target_points) { @@ -264,18 +238,29 @@ class MovingLeastSquares Details::MLSTargetPointsPredicateWrapper predicates{ {target_points}, _num_neighbors}; - // Query the source - _indices = Kokkos::View( - "ArborX::MovingLeastSquares::indices", 0); - Kokkos::View offsets( - "ArborX::MovingLeastSquares::offsets", 0); - source_tree.query(space, predicates, - ArborX::Details::LegacyDefaultCallback{}, _indices, - offsets); + // Create the callback + Kokkos::View source_view( + Kokkos::view_alloc(space, Kokkos::WithoutInitializing, + "ArborX::MovingLeastSquares::source_view"), + _num_targets, _num_neighbors); + _indices = Kokkos::View( + Kokkos::view_alloc(space, Kokkos::WithoutInitializing, + "ArborX::MovingLeastSquares::indices"), + _num_targets, _num_neighbors); + Kokkos::View counter( + "ArborX::MovingLeastSquares::counter", _num_targets); + Details::SearchNeighborsCallback + callback{source_view, _indices, counter}; + + // Query the source tree + source_tree.query(space, predicates, callback); + + return source_view; } Kokkos::View _coeffs; - Kokkos::View _indices; + Kokkos::View _indices; int _num_targets; int _num_neighbors; int _source_size; From a1f89a50ade4e2efd005a5b8f13e56ff609f6575 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Sun, 24 Dec 2023 16:57:32 -0500 Subject: [PATCH 14/15] Passing around AccessValues instead of raw user inputs --- .../ArborX_InterpMovingLeastSquares.hpp | 44 +++++++++---------- ...pDetailsMovingLeastSquaresCoefficients.hpp | 14 +++--- 2 files changed, 26 insertions(+), 32 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index 373cf3949..7919d0203 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -31,10 +31,10 @@ namespace ArborX::Interpolation::Details { // This is done to avoid a clash with another predicates access trait -template +template struct MLSTargetPointsPredicateWrapper { - ArborX::Details::AccessValues target_access; + TargetAccess target_access; int num_neighbors; }; @@ -66,26 +66,25 @@ struct SearchNeighborsCallback namespace ArborX { -template +template struct AccessTraits< - Interpolation::Details::MLSTargetPointsPredicateWrapper, + Interpolation::Details::MLSTargetPointsPredicateWrapper, PredicatesTag> { - KOKKOS_FUNCTION static auto size( - Interpolation::Details::MLSTargetPointsPredicateWrapper const &tp) + using Self = + Interpolation::Details::MLSTargetPointsPredicateWrapper; + + KOKKOS_FUNCTION static auto size(Self const &tp) { return tp.target_access.size(); } - KOKKOS_FUNCTION static auto - get(Interpolation::Details::MLSTargetPointsPredicateWrapper const &tp, - int const i) + KOKKOS_FUNCTION static auto get(Self const &tp, int const i) { return attach(nearest(tp.target_access(i), tp.num_neighbors), i); } - using memory_space = - typename Details::AccessValues::memory_space; + using memory_space = typename TargetAccess::memory_space; }; } // namespace ArborX @@ -150,18 +149,19 @@ class MovingLeastSquares TargetAccess target_access{target_points}; SourceAccess source_access{source_points}; + _num_targets = target_access.size(); _source_size = source_access.size(); // There must be enough source points KOKKOS_ASSERT(0 < _num_neighbors && _num_neighbors <= _source_size); // Search for neighbors and get the arranged source points - auto source_view = searchNeighbors(space, source_points, target_points); + auto source_view = searchNeighbors(space, source_access, target_access); // Compute the moving least squares coefficients _coeffs = Details::movingLeastSquaresCoefficients( - space, source_view, target_points); + space, source_view, target_access._values); } template + template auto searchNeighbors(ExecutionSpace const &space, - SourcePoints const &source_points, - TargetPoints const &target_points) + SourceAccess const &source_access, + TargetAccess const &target_access) { auto guard = Kokkos::Profiling::ScopedRegion( "ArborX::MovingLeastSquares::searchNeighbors"); // Organize the source points as a tree - using SourcePoint = - typename ArborX::Details::AccessValues::value_type; + using SourcePoint = typename SourceAccess::value_type; BoundingVolumeHierarchy> - source_tree(space, ArborX::Experimental::attach_indices(source_points)); + source_tree(space, ArborX::Experimental::attach_indices(source_access)); // Create the predicates - Details::MLSTargetPointsPredicateWrapper predicates{ - {target_points}, _num_neighbors}; + Details::MLSTargetPointsPredicateWrapper predicates{ + target_access, _num_neighbors}; // Create the callback Kokkos::View source_view( diff --git a/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp b/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp index 47152418e..8261520e7 100644 --- a/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp +++ b/src/interpolation/details/ArborX_InterpDetailsMovingLeastSquaresCoefficients.hpp @@ -28,11 +28,11 @@ namespace ArborX::Interpolation::Details template + typename SourcePoints, typename TargetAccess> Kokkos::View movingLeastSquaresCoefficients(ExecutionSpace const &space, SourcePoints const &source_points, - TargetPoints const &target_points) + TargetAccess const &target_access) { auto guard = Kokkos::Profiling::ScopedRegion("ArborX::MovingLeastSquaresCoefficients"); @@ -57,22 +57,18 @@ movingLeastSquaresCoefficients(ExecutionSpace const &space, "source points elements must be points"); static constexpr int dimension = GeometryTraits::dimension_v; - // TargetPoints is an access trait of points - ArborX::Details::check_valid_access_traits(PrimitivesTag{}, target_points); - using TargetAccess = - ArborX::Details::AccessValues; + // TargetAccess is an access values of points static_assert( KokkosExt::is_accessible_from::value, - "target points must be accessible from the execution space"); + "target access must be accessible from the execution space"); using TargetPoint = typename TargetAccess::value_type; GeometryTraits::check_valid_geometry_traits(TargetPoint{}); static_assert(GeometryTraits::is_point::value, - "target points elements must be points"); + "target access elements must be points"); static_assert(dimension == GeometryTraits::dimension_v, "target and source points must have the same dimension"); - TargetAccess target_access{target_points}; int const num_targets = target_access.size(); int const num_neighbors = source_points.extent(1); From 9f99b90beaa8a9bdfc9627ce8268371dfd355f60 Mon Sep 17 00:00:00 2001 From: Yohann Bosqued Date: Tue, 26 Dec 2023 21:30:43 -0500 Subject: [PATCH 15/15] Renaming external predicates / callback --- .../ArborX_InterpMovingLeastSquares.hpp | 20 +++++++++---------- 1 file changed, 9 insertions(+), 11 deletions(-) diff --git a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp index 7919d0203..c6f418301 100644 --- a/src/interpolation/ArborX_InterpMovingLeastSquares.hpp +++ b/src/interpolation/ArborX_InterpMovingLeastSquares.hpp @@ -32,7 +32,7 @@ namespace ArborX::Interpolation::Details // This is done to avoid a clash with another predicates access trait template -struct MLSTargetPointsPredicateWrapper +struct MLSPredicateWrapper { TargetAccess target_access; int num_neighbors; @@ -40,7 +40,7 @@ struct MLSTargetPointsPredicateWrapper // Functor used in the tree query to create the 2D source view and indices template -struct SearchNeighborsCallback +struct MLSSearchNeighborsCallback { SourceView source_view; IndicesView indices; @@ -67,12 +67,10 @@ namespace ArborX { template -struct AccessTraits< - Interpolation::Details::MLSTargetPointsPredicateWrapper, - PredicatesTag> +struct AccessTraits, + PredicatesTag> { - using Self = - Interpolation::Details::MLSTargetPointsPredicateWrapper; + using Self = Interpolation::Details::MLSPredicateWrapper; KOKKOS_FUNCTION static auto size(Self const &tp) { @@ -233,8 +231,8 @@ class MovingLeastSquares source_tree(space, ArborX::Experimental::attach_indices(source_access)); // Create the predicates - Details::MLSTargetPointsPredicateWrapper predicates{ - target_access, _num_neighbors}; + Details::MLSPredicateWrapper predicates{target_access, + _num_neighbors}; // Create the callback Kokkos::View source_view( @@ -247,8 +245,8 @@ class MovingLeastSquares _num_targets, _num_neighbors); Kokkos::View counter( "ArborX::MovingLeastSquares::counter", _num_targets); - Details::SearchNeighborsCallback + Details::MLSSearchNeighborsCallback callback{source_view, _indices, counter}; // Query the source tree