From: Daniel Arndt Date: Tue, 11 Jul 2023 18:42:53 +0000 (-0400) Subject: Add 1d support for ArborXWrappers X-Git-Tag: relicensing~683^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=74709e1564c57f30b87d4599436f01ae0470c610;p=dealii.git Add 1d support for ArborXWrappers --- diff --git a/source/arborx/access_traits.cc b/source/arborx/access_traits.cc index 1520624d49..3c37bc6982 100644 --- a/source/arborx/access_traits.cc +++ b/source/arborx/access_traits.cc @@ -26,16 +26,15 @@ namespace ArborXWrappers PointPredicate::PointPredicate( const std::vector> &dim_points) { - static_assert(dim != 1, "dim equal to one is not supported."); - const unsigned int size = dim_points.size(); points.reserve(size); for (unsigned int i = 0; i < size; ++i) { points.emplace_back(static_cast(dim_points[i][0]), - static_cast(dim_points[i][1]), - dim == 2 ? 0.f : - static_cast(dim_points[i][2])); + dim < 2 ? 0.f : + static_cast(dim_points[i][1]), + dim < 3 ? 0.f : + static_cast(dim_points[i][2])); } } @@ -153,20 +152,18 @@ namespace ArborXWrappers const std::vector, Number>> &dim_spheres) { - static_assert(dim != 1, "dim equal to one is not supported."); - const unsigned int size = dim_spheres.size(); spheres.reserve(size); for (unsigned int i = 0; i < size; ++i) { // ArborX assumes that the center coordinates and the radius use float // and the sphere is 3d - spheres.emplace_back(std::make_pair( + spheres.emplace_back( dealii::Point<3, float>( static_cast(dim_spheres[i].first[0]), - static_cast(dim_spheres[i].first[1]), - dim == 2 ? 0.f : static_cast(dim_spheres[i].first[2])), - static_cast(dim_spheres[i].second))); + dim < 2 ? 0.f : static_cast(dim_spheres[i].first[1]), + dim < 3 ? 0.f : static_cast(dim_spheres[i].first[2])), + static_cast(dim_spheres[i].second)); } } @@ -237,8 +234,8 @@ namespace ArborX // ArborX assumes that the point coordinates use float and that the point // is 3d return {static_cast(v[i][0]), - static_cast(v[i][1]), - dim == 2 ? 0 : static_cast(v[i][2])}; + dim < 2 ? 0 : static_cast(v[i][1]), + dim < 3 ? 0 : static_cast(v[i][2])}; } @@ -265,11 +262,11 @@ namespace ArborX // ArborX assumes that the bounding box coordinates use float and that the // bounding box is 3d return {{static_cast(min_corner[0]), - static_cast(min_corner[1]), - dim == 2 ? 0.f : static_cast(min_corner[2])}, + dim < 2 ? 0.f : static_cast(min_corner[1]), + dim < 3 ? 0.f : static_cast(min_corner[2])}, {static_cast(max_corner[0]), - static_cast(max_corner[1]), - dim == 2 ? 0.f : static_cast(max_corner[2])}}; + dim < 2 ? 0.f : static_cast(max_corner[1]), + dim < 3 ? 0.f : static_cast(max_corner[2])}}; } @@ -296,8 +293,8 @@ namespace ArborX // ArborX assumes that the center coordinates and the radius use float and // the sphere is 3d return {{static_cast(v[i].first[0]), - static_cast(v[i].first[1]), - dim == 2 ? 0 : static_cast(v[i].first[2])}, + dim < 2 ? 0 : static_cast(v[i].first[1]), + dim < 3 ? 0 : static_cast(v[i].first[2])}, static_cast(v[i].second)}; } } // namespace ArborX diff --git a/source/arborx/access_traits.inst.in b/source/arborx/access_traits.inst.in index f82c483792..45ca6c5b54 100644 --- a/source/arborx/access_traits.inst.in +++ b/source/arborx/access_traits.inst.in @@ -16,7 +16,6 @@ for (DIM : DIMENSIONS; SCALAR : REAL_SCALARS) { -#if DIM > 1 DEAL_II_NAMESPACE_OPEN namespace ArborXWrappers { @@ -55,5 +54,4 @@ for (DIM : DIMENSIONS; SCALAR : REAL_SCALARS) std::vector, SCALAR>>, PrimitivesTag>; \} -#endif } diff --git a/tests/arborx/bounding_box_intersect.cc b/tests/arborx/bounding_box_intersect.cc index 30aec63c20..43f53c2a62 100644 --- a/tests/arborx/bounding_box_intersect.cc +++ b/tests/arborx/bounding_box_intersect.cc @@ -24,6 +24,62 @@ #include "../tests.h" +void +test_1d() +{ + std::vector> bounding_boxes; + std::vector> points; + + unsigned int n_points_1d = 5; + for (unsigned int i = 0; i < n_points_1d; ++i) + points.emplace_back(i); + + for (unsigned int i = 0; i < n_points_1d - 1; ++i) + bounding_boxes.emplace_back(std::make_pair(points[i], points[i + 1])); + + Point<1> point_a(0.5); + Point<1> point_b(1.5); + Point<1> point_c(2.2); + Point<1> point_d(2.6); + + std::vector> query_bounding_boxes; + query_bounding_boxes.emplace_back(std::make_pair(point_a, point_b)); + query_bounding_boxes.emplace_back(std::make_pair(point_c, point_d)); + + ArborXWrappers::BVH bvh(bounding_boxes); + ArborXWrappers::BoundingBoxIntersectPredicate bb_intersect( + query_bounding_boxes); + auto indices_offset = bvh.query(bb_intersect); + std::vector indices = indices_offset.first; + std::vector offset = indices_offset.second; + + std::vector indices_ref = {0, 1, 2}; + std::vector offset_ref = {0, 2, 3}; + + AssertDimension(indices.size(), indices_ref.size()); + AssertDimension(offset.size(), offset_ref.size()); + for (unsigned int i = 0; i < offset.size() - 1; ++i) + { + for (unsigned int j = offset[i]; j < offset[i + 1]; ++j) + { + // The indices associated to each query are not ordered. + bool found = false; + for (unsigned int k = offset[i]; k < offset[i + 1]; ++k) + { + if (indices[j] == indices_ref[k]) + { + found = true; + break; + } + } + AssertThrow(found, ExcInternalError()); + } + } + for (unsigned int i = 0; i < offset.size(); ++i) + AssertThrow(offset[i] == offset_ref[i], ExcInternalError()); + + deallog << "OK" << std::endl; +} void test_2d() @@ -181,6 +237,7 @@ main(int argc, char **argv) initlog(); // tests + test_1d(); test_2d(); test_3d(); diff --git a/tests/arborx/bounding_box_intersect.output b/tests/arborx/bounding_box_intersect.output index 8b3b075900..fb71de2867 100644 --- a/tests/arborx/bounding_box_intersect.output +++ b/tests/arborx/bounding_box_intersect.output @@ -1,3 +1,4 @@ DEAL::OK DEAL::OK +DEAL::OK diff --git a/tests/arborx/distributed_tree_intersect.cc b/tests/arborx/distributed_tree_intersect.cc index b3ec93a610..76cb632e65 100644 --- a/tests/arborx/distributed_tree_intersect.cc +++ b/tests/arborx/distributed_tree_intersect.cc @@ -25,6 +25,70 @@ #include "../tests.h" +void +test_1d() +{ + std::vector> bounding_boxes; + std::vector> points; + + const unsigned int n_points_1d = 5; + const int rank = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + const int n_procs = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + const int offset_bb = 2 * rank * n_points_1d; + const int offset_query = 2 * ((rank + 1) % n_procs) * n_points_1d; + for (unsigned int i = 0; i < n_points_1d; ++i) + points.emplace_back(i + offset_bb); + + for (unsigned int i = 0; i < n_points_1d - 1; ++i) + bounding_boxes.emplace_back(std::make_pair(points[i], points[i + 1])); + + Point<1> point_a(0.5 + offset_query); + Point<1> point_b(1.5 + offset_query); + Point<1> point_c(2.2 + offset_query); + Point<1> point_d(2.6 + offset_query); + + std::vector> query_bounding_boxes; + query_bounding_boxes.emplace_back(std::make_pair(point_a, point_b)); + query_bounding_boxes.emplace_back(std::make_pair(point_c, point_d)); + + ArborXWrappers::DistributedTree distributed_tree(MPI_COMM_WORLD, + bounding_boxes); + ArborXWrappers::BoundingBoxIntersectPredicate bb_intersect( + query_bounding_boxes); + auto indices_ranks_offset = distributed_tree.query(bb_intersect); + auto indices_ranks = indices_ranks_offset.first; + auto offsets = indices_ranks_offset.second; + + std::vector indices_ref = {1, 0, 2}; + std::vector offsets_ref = {0, 2, 3}; + + AssertThrow(indices_ranks.size() == indices_ref.size(), ExcInternalError()); + AssertThrow(offsets.size() == offsets_ref.size(), ExcInternalError()); + for (unsigned int i = 0; i < offsets.size() - 1; ++i) + { + for (int j = offsets[i]; j < offsets[i + 1]; ++j) + { + // The indices associated to each query are not ordered. + bool found = false; + for (int k = offsets[i]; k < offsets[i + 1]; ++k) + { + if ((indices_ranks[j].first == indices_ref[k]) && + (indices_ranks[j].second == (rank + 1) % n_procs)) + { + found = true; + break; + } + } + AssertThrow(found, ExcInternalError()); + } + } + for (unsigned int i = 0; i < offsets.size(); ++i) + AssertThrow(offsets[i] == offsets_ref[i], ExcInternalError()); + + deallog << "OK" << std::endl; +} + + void test_2d() { @@ -192,6 +256,7 @@ main(int argc, char **argv) initlog(); // tests + test_1d(); test_2d(); test_3d(); } diff --git a/tests/arborx/distributed_tree_intersect.with_arborx_with_mpi=on.mpirun=2.output b/tests/arborx/distributed_tree_intersect.with_arborx_with_mpi=on.mpirun=2.output index 9282f706bc..fb71de2867 100644 --- a/tests/arborx/distributed_tree_intersect.with_arborx_with_mpi=on.mpirun=2.output +++ b/tests/arborx/distributed_tree_intersect.with_arborx_with_mpi=on.mpirun=2.output @@ -1,3 +1,4 @@ -DEAL::OK +DEAL::OK +DEAL::OK DEAL::OK diff --git a/tests/arborx/point_intersect.cc b/tests/arborx/point_intersect.cc index e0ff970f26..401d1bac7c 100644 --- a/tests/arborx/point_intersect.cc +++ b/tests/arborx/point_intersect.cc @@ -25,6 +25,61 @@ #include "../tests.h" +void +test_1d() +{ + std::vector> bounding_boxes; + std::vector> points; + + unsigned int n_points_1d = 5; + for (unsigned int i = 0; i < n_points_1d; ++i) + points.emplace_back(i); + + for (unsigned int i = 0; i < n_points_1d - 1; ++i) + bounding_boxes.emplace_back(std::make_pair(points[i], points[i + 1])); + + std::vector> query_points; + query_points.emplace_back(0.5); + query_points.emplace_back(1.5); + query_points.emplace_back(2.2); + query_points.emplace_back(2.6); + + + ArborXWrappers::BVH bvh(bounding_boxes); + ArborXWrappers::PointIntersectPredicate pt_intersect(query_points); + auto indices_offset = bvh.query(pt_intersect); + std::vector indices = indices_offset.first; + std::vector offset = indices_offset.second; + + std::vector indices_ref = {0, 1, 2, 2}; + std::vector offset_ref = {0, 1, 2, 3, 4}; + + AssertThrow(indices.size() == indices_ref.size(), ExcInternalError()); + AssertThrow(offset.size() == offset_ref.size(), ExcInternalError()); + for (unsigned int i = 0; i < offset.size() - 1; ++i) + { + for (unsigned int j = offset[i]; j < offset[i + 1]; ++j) + { + // The indices associated to each query are not ordered. + bool found = false; + for (unsigned int k = offset[i]; k < offset[i + 1]; ++k) + { + if (indices[j] == indices_ref[k]) + { + found = true; + break; + } + } + AssertThrow(found, ExcInternalError()); + } + } + for (unsigned int i = 0; i < offset.size(); ++i) + AssertThrow(offset[i] == offset_ref[i], ExcInternalError()); + + deallog << "OK" << std::endl; +} + + void test_2d() { @@ -176,6 +231,7 @@ main(int argc, char **argv) Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv); // tests + test_1d(); test_2d(); test_3d(); } diff --git a/tests/arborx/point_intersect.output b/tests/arborx/point_intersect.output index 8b3b075900..fb71de2867 100644 --- a/tests/arborx/point_intersect.output +++ b/tests/arborx/point_intersect.output @@ -1,3 +1,4 @@ DEAL::OK DEAL::OK +DEAL::OK diff --git a/tests/arborx/sphere_intersect.cc b/tests/arborx/sphere_intersect.cc index a9a807d18d..b46ecee416 100644 --- a/tests/arborx/sphere_intersect.cc +++ b/tests/arborx/sphere_intersect.cc @@ -25,6 +25,56 @@ #include "../tests.h" +void +test_1d() +{ + std::vector> points; + + unsigned int n_points_1d = 5; + for (unsigned int i = 0; i < n_points_1d; ++i) + points.emplace_back(i); + + std::vector, double>> query_spheres; + query_spheres.emplace_back(Point<1>(0.5), 0.8); + query_spheres.emplace_back(Point<1>(1.5), 0.8); + query_spheres.emplace_back(Point<1>(2.2), 1.2); + query_spheres.emplace_back(Point<1>(2.6), 0.9); + + + ArborXWrappers::BVH bvh(points); + ArborXWrappers::SphereIntersectPredicate sph_intersect(query_spheres); + auto indices_offsets = bvh.query(sph_intersect); + std::vector indices = indices_offsets.first; + std::vector offsets = indices_offsets.second; + + std::vector indices_ref = {0, 1, 1, 2, 1, 2, 3, 2, 3}; + std::vector offsets_ref = {0, 2, 4, 7, 9}; + + AssertThrow(indices.size() == indices_ref.size(), ExcInternalError()); + AssertThrow(offsets.size() == offsets_ref.size(), ExcInternalError()); + for (unsigned int i = 0; i < offsets.size() - 1; ++i) + { + for (int j = offsets[i]; j < offsets[i + 1]; ++j) + { + // The indices associated to each query are not ordered. + bool found = false; + for (int k = offsets[i]; k < offsets[i + 1]; ++k) + { + if (indices[j] == indices_ref[k]) + { + found = true; + break; + } + } + AssertThrow(found, ExcInternalError()); + } + } + for (unsigned int i = 0; i < offsets.size(); ++i) + AssertThrow(offsets[i] == offsets_ref[i], ExcInternalError()); + + deallog << "OK" << std::endl; +} + void test_2d() { @@ -150,6 +200,7 @@ main(int argc, char **argv) Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv); // tests + test_1d(); test_2d(); test_3d(); } diff --git a/tests/arborx/sphere_intersect.output b/tests/arborx/sphere_intersect.output index 8b3b075900..fb71de2867 100644 --- a/tests/arborx/sphere_intersect.output +++ b/tests/arborx/sphere_intersect.output @@ -1,3 +1,4 @@ DEAL::OK DEAL::OK +DEAL::OK diff --git a/tests/arborx/sphere_nearest.cc b/tests/arborx/sphere_nearest.cc index dae5af87c3..8ff19ea7a7 100644 --- a/tests/arborx/sphere_nearest.cc +++ b/tests/arborx/sphere_nearest.cc @@ -26,6 +26,57 @@ +void +test_1d() +{ + std::vector> points; + + unsigned int n_points_1d = 5; + for (unsigned int i = 0; i < n_points_1d; ++i) + points.emplace_back(i); + + std::vector, double>> query_spheres; + query_spheres.emplace_back(Point<1>(0.5), 0.1); + query_spheres.emplace_back(Point<1>(1.5), 0.1); + query_spheres.emplace_back(Point<1>(2.2), 0.1); + query_spheres.emplace_back(Point<1>(2.6), 0.1); + + + ArborXWrappers::BVH bvh(points); + ArborXWrappers::SphereNearestPredicate sph_nearest(query_spheres, 1); + auto indices_offsets = bvh.query(sph_nearest); + std::vector indices = indices_offsets.first; + std::vector offsets = indices_offsets.second; + + std::vector indices_ref = {0, 1, 2, 3}; + std::vector offsets_ref = {0, 1, 2, 3, 4}; + + AssertThrow(indices.size() == indices_ref.size(), ExcInternalError()); + AssertThrow(offsets.size() == offsets_ref.size(), ExcInternalError()); + for (unsigned int i = 0; i < offsets.size() - 1; ++i) + { + for (int j = offsets[i]; j < offsets[i + 1]; ++j) + { + // The indices associated to each query are not ordered. + bool found = false; + for (int k = offsets[i]; k < offsets[i + 1]; ++k) + { + if (indices[j] == indices_ref[k]) + { + found = true; + break; + } + } + AssertThrow(found, ExcInternalError()); + } + } + for (unsigned int i = 0; i < offsets.size(); ++i) + AssertThrow(offsets[i] == offsets_ref[i], ExcInternalError()); + + deallog << "OK" << std::endl; +} + + void test_2d() { @@ -149,6 +200,7 @@ main(int argc, char **argv) Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv); // tests + test_1d(); test_2d(); test_3d(); } diff --git a/tests/arborx/sphere_nearest.output b/tests/arborx/sphere_nearest.output index 8b3b075900..fb71de2867 100644 --- a/tests/arborx/sphere_nearest.output +++ b/tests/arborx/sphere_nearest.output @@ -1,3 +1,4 @@ DEAL::OK DEAL::OK +DEAL::OK