From: Bruno Turcksin Date: Fri, 11 Mar 2022 16:35:54 +0000 (+0000) Subject: Add wrappers for ArborX::Sphere X-Git-Tag: v9.4.0-rc1~335^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=5c52f332750a5675b9460a3b098719aa166db6c0;p=dealii.git Add wrappers for ArborX::Sphere --- diff --git a/include/deal.II/arborx/access_traits.h b/include/deal.II/arborx/access_traits.h index 0f3509ba85..cf9caeadf9 100644 --- a/include/deal.II/arborx/access_traits.h +++ b/include/deal.II/arborx/access_traits.h @@ -23,6 +23,8 @@ # include +# include + DEAL_II_NAMESPACE_OPEN @@ -42,7 +44,7 @@ namespace ArborXWrappers PointPredicate(const std::vector> &points); /** - * Number of points stored in the structure. + * The number of points stored in the structure. */ std::size_t size() const; @@ -133,7 +135,7 @@ namespace ArborXWrappers const std::vector> &bounding_boxes); /** - * Number of bounding boxes stored in the structure. + * The number of bounding boxes stored in the structure. */ std::size_t size() const; @@ -205,6 +207,100 @@ namespace ArborXWrappers private: unsigned int n_nearest_neighbors; }; + + + + /** + * Base class for Sphere-based predicates providing basic functionality for + * derived classes; not supposed to be used on its own. + */ + class SpherePredicate + { + protected: + /** + * Constructor. @p spheres is a list of spheres, i.e. a center and radius + * pair, used by the predicate. + */ + template + SpherePredicate( + const std::vector, Number>> + &spheres); + + /** + * The number of spheres stored in the structure. + */ + std::size_t + size() const; + + /** + * Return the `i`th sphere stored in the object. + */ + const std::pair, float> & + get(unsigned int) const; + + private: + std::vector, float>> spheres; + }; + + + + /** + * This class defines a predicate used by ArborXWrappers::BVH to determine + * for given spheres which of the bounding boxes used to build the + * ArborXWrappers::BVH intersect with them. + * @note The class is not supposed to be used in a polymorphic context. + */ + class SphereIntersectPredicate : private SpherePredicate + { + public: + /** + * Constructor. @p spheres is a list of spheres which we are interested in + * knowing if they intersect ArborXWrappers::BVH bounding boxes. + */ + template + SphereIntersectPredicate( + const std::vector, Number>> + &spheres); + + // We need these since we inherit privately to avoid polymorphic use. + using SpherePredicate::get; + using SpherePredicate::size; + }; + + + + /** + * This class defines a predicate used by ArborXWrappers::BVH to determine, + * for the given spheres, which are the nearest bounding boxes/points among + * the ones used to build the ArborXWrappers::BVH. + * @note The class is not supposed to be used in a polymorphic context. + */ + class SphereNearestPredicate : private SpherePredicate + { + public: + /** + * Constructor. @p spheres is a list of spheres for which we are + * interested in the @p n_nearest_neighbors in the ArborXWrappers::BVH + * bounding boxes/points. + */ + template + SphereNearestPredicate( + const std::vector, Number>> &spheres, + const unsigned int n_nearest_neighbors); + + /** + * Return the number of nearest neighbors we are looking for. + */ + unsigned int + get_n_nearest_neighbors() const; + + // We need these since we inherit privately to avoid polymorphic use. + using SpherePredicate::get; + using SpherePredicate::size; + + private: + unsigned int n_nearest_neighbors; + }; } // namespace ArborXWrappers DEAL_II_NAMESPACE_CLOSE @@ -263,6 +359,34 @@ namespace ArborX + /** + * This struct allows ArborX to use + * std::vector as primitive. + */ + template + struct AccessTraits< + std::vector, Number>>, + PrimitivesTag> + { + using memory_space = Kokkos::HostSpace; + + /** + * Return the size of the vector @p v. + */ + static std::size_t + size(const std::vector, Number>> &v); + + /** + * Return an ArborX::Sphere from the std::pair + * `v[i]`. + */ + static Sphere + get(const std::vector, Number>> &v, + std::size_t i); + }; + + + /** * This struct allows ArborX to use PointIntersectPredicate as a predicate. */ @@ -273,13 +397,13 @@ namespace ArborX using memory_space = Kokkos::HostSpace; /** - * Number of Point stored in @p pt_intersect. + * The number of points stored in @p pt_intersect. */ static std::size_t size(const dealii::ArborXWrappers::PointIntersectPredicate &pt_intersect); /** - * Return an Arbox::intersects(ArborX::Point) object constructed from the + * Return an ArborX::intersects(ArborX::Point) object constructed from the * `i`th dealii::Point stored in @p pt_intersect. */ static auto @@ -288,6 +412,7 @@ namespace ArborX }; + /** * This struct allows ArborX to use PointNearestPredicate as a predicate. */ @@ -298,13 +423,13 @@ namespace ArborX using memory_space = Kokkos::HostSpace; /** - * Number of Point stored in @p pt_nearest. + * The number of points stored in @p pt_nearest. */ static std::size_t size(const dealii::ArborXWrappers::PointNearestPredicate &pt_nearest); /** - * Return an Arbox::nearest(ArborX::Point, + * Return an ArborX::nearest(ArborX::Point, * PointNearestPredicate::get_n_nearest_neighbors) object constructed from * the `i`th dealii::Point stored in @p pt_nearest. */ @@ -314,6 +439,7 @@ namespace ArborX }; + /** * This struct allows ArborX to use BoundingBoxIntersectPredicate as a * predicate. @@ -325,7 +451,7 @@ namespace ArborX using memory_space = Kokkos::HostSpace; /** - * Number of BoundingBox stored in @p bb_intersect. + * The number of bounding boxes stored in @p bb_intersect. */ static std::size_t size(const dealii::ArborXWrappers::BoundingBoxIntersectPredicate @@ -342,6 +468,7 @@ namespace ArborX }; + /** * This struct allows ArborX to use BoundingBoxNearstPredicate as a * predicate. @@ -353,7 +480,7 @@ namespace ArborX using memory_space = Kokkos::HostSpace; /** - * Number of BoundingBox stored in @p bb_nearest. + * The number of bounding boxes stored in @p bb_nearest. */ static std::size_t size(const dealii::ArborXWrappers::BoundingBoxNearestPredicate &bb_nearest); @@ -370,6 +497,59 @@ namespace ArborX std::size_t i); }; + + + /** + * This struct allows ArborX to use SphereIntersectPredicate as a predicate. + */ + template <> + struct AccessTraits + { + using memory_space = Kokkos::HostSpace; + + /** + * The number of points stored in @p sph_intersect. + */ + static std::size_t + size(const dealii::ArborXWrappers::SphereIntersectPredicate &sph_intersect); + + /** + * Return an ArborX::intersects(ArborX::Sphere) object constructed from the + * `i`th sphere stored in @p sph_intersect. + */ + static auto + get(const dealii::ArborXWrappers::SphereIntersectPredicate &sph_intersect, + std::size_t i); + }; + + + + /** + * This struct allows ArborX to use SphereNearestPredicate as a predicate. + */ + template <> + struct AccessTraits + { + using memory_space = Kokkos::HostSpace; + + /** + * The number of spheres stored in @p sph_nearest. + */ + static std::size_t + size(const dealii::ArborXWrappers::SphereNearestPredicate &sph_nearest); + + /** + * Return an ArborX::nearest(ArborX::Sphere, + * SphereNearestPredicate::get_n_nearest_neightbors) object constructed from + * the `i`th sphere stored in @p sph_nearest. + */ + static auto + get(const dealii::ArborXWrappers::SphereNearestPredicate &sph_nearest, + std::size_t i); + }; + // ------------------------------- Inline ----------------------------------// // The implementation of AccessTraits<..., PredicatesTag> needs to be in the @@ -469,6 +649,52 @@ namespace ArborX {max_corner[0], max_corner[1], max_corner[2]}}, bb_nearest.get_n_nearest_neighbors()); } + + + + inline std::size_t + AccessTraits:: + size(const dealii::ArborXWrappers::SphereIntersectPredicate &sph_intersect) + { + return sph_intersect.size(); + } + + + + inline auto + AccessTraits:: + get(const dealii::ArborXWrappers::SphereIntersectPredicate &sph_intersect, + std::size_t i) + { + const auto sphere = sph_intersect.get(i); + return intersects( + Sphere{{sphere.first[0], sphere.first[1], sphere.first[2]}, + sphere.second}); + } + + + + inline std::size_t + AccessTraits:: + size(const dealii::ArborXWrappers::SphereNearestPredicate &sph_nearest) + { + return sph_nearest.size(); + } + + + + inline auto + AccessTraits:: + get(const dealii::ArborXWrappers::SphereNearestPredicate &sph_nearest, + std::size_t i) + { + const auto sphere = sph_nearest.get(i); + return nearest(Sphere{{sphere.first[0], sphere.first[1], sphere.first[2]}, + sphere.second}, + sph_nearest.get_n_nearest_neighbors()); + } } // namespace ArborX #endif diff --git a/source/arborx/access_traits.cc b/source/arborx/access_traits.cc index c8b3d3473d..51e008f38d 100644 --- a/source/arborx/access_traits.cc +++ b/source/arborx/access_traits.cc @@ -146,6 +146,71 @@ namespace ArborXWrappers { return n_nearest_neighbors; } + + // ------------------- SpherePredicate ------------------- // + template + SpherePredicate::SpherePredicate( + 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( + 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))); + } + } + + + + std::size_t + SpherePredicate::size() const + { + return spheres.size(); + } + + + + const std::pair, float> & + SpherePredicate::get(unsigned int i) const + { + return spheres[i]; + } + + + + template + SphereIntersectPredicate::SphereIntersectPredicate( + const std::vector, Number>> &spheres) + : SpherePredicate(spheres) + {} + + + + template + SphereNearestPredicate::SphereNearestPredicate( + const std::vector, Number>> &spheres, + const unsigned int n_nearest_neighbors) + : SpherePredicate(spheres) + , n_nearest_neighbors(n_nearest_neighbors) + {} + + + + unsigned int + SphereNearestPredicate::get_n_nearest_neighbors() const + { + return n_nearest_neighbors; + } } // namespace ArborXWrappers DEAL_II_NAMESPACE_CLOSE @@ -206,6 +271,35 @@ namespace ArborX static_cast(max_corner[1]), dim == 2 ? 0.f : static_cast(max_corner[2])}}; } + + + + // ---------------------- Sphere Primitives AccessTraits ----------------- // + template + std::size_t + AccessTraits, Number>>, + PrimitivesTag>:: + size(const std::vector, Number>> &v) + { + return v.size(); + } + + + + template + Sphere + AccessTraits, Number>>, + PrimitivesTag>:: + get(const std::vector, Number>> &v, + std::size_t i) + { + // 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])}, + static_cast(v[i].second)}; + } } // namespace ArborX // ----------------------- Instantiations --------------------// diff --git a/source/arborx/access_traits.inst.in b/source/arborx/access_traits.inst.in index cc949b0f12..8cd79a4198 100644 --- a/source/arborx/access_traits.inst.in +++ b/source/arborx/access_traits.inst.in @@ -31,6 +31,14 @@ for (DIM : DIMENSIONS; SCALAR : REAL_SCALARS) template BoundingBoxNearestPredicate::BoundingBoxNearestPredicate( const std::vector> &bb, const unsigned int n_nearest_neighbors); + + template SphereIntersectPredicate::SphereIntersectPredicate( + const std::vector, SCALAR>> + &spheres); + template SphereNearestPredicate::SphereNearestPredicate( + const std::vector, SCALAR>> + & spheres, + const unsigned int n_nearest_neighbors); \} DEAL_II_NAMESPACE_CLOSE @@ -38,9 +46,14 @@ for (DIM : DIMENSIONS; SCALAR : REAL_SCALARS) { template struct AccessTraits>, PrimitivesTag>; + template struct AccessTraits< std::vector>, PrimitivesTag>; + + template struct AccessTraits< + std::vector, SCALAR>>, + PrimitivesTag>; \} #endif }