]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add wrappers for ArborX::Sphere
authorBruno Turcksin <bruno.turcksin@gmail.com>
Fri, 11 Mar 2022 16:35:54 +0000 (16:35 +0000)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Fri, 25 Mar 2022 15:47:20 +0000 (15:47 +0000)
include/deal.II/arborx/access_traits.h
source/arborx/access_traits.cc
source/arborx/access_traits.inst.in

index 0f3509ba85b2ec835678462f3bc4a1a3593cf7ff..cf9caeadf92e10e8d3f004ec86ccb9170ea761e9 100644 (file)
@@ -23,6 +23,8 @@
 
 #  include <ArborX.hpp>
 
+#  include <utility>
+
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -42,7 +44,7 @@ namespace ArborXWrappers
     PointPredicate(const std::vector<dealii::Point<dim, Number>> &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<dealii::BoundingBox<dim, Number>> &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 <int dim, typename Number>
+    SpherePredicate(
+      const std::vector<std::pair<dealii::Point<dim, Number>, 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<dealii::Point<3, float>, float> &
+    get(unsigned int) const;
+
+  private:
+    std::vector<std::pair<dealii::Point<3, float>, 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 <int dim, typename Number>
+    SphereIntersectPredicate(
+      const std::vector<std::pair<dealii::Point<dim, Number>, 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 <int dim, typename Number>
+    SphereNearestPredicate(
+      const std::vector<std::pair<dealii::Point<dim, Number>, 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<std::pair<dealii::Point, Number> as primitive.
+   */
+  template <int dim, typename Number>
+  struct AccessTraits<
+    std::vector<std::pair<dealii::Point<dim, Number>, Number>>,
+    PrimitivesTag>
+  {
+    using memory_space = Kokkos::HostSpace;
+
+    /**
+     * Return the size of the vector @p v.
+     */
+    static std::size_t
+    size(const std::vector<std::pair<dealii::Point<dim, Number>, Number>> &v);
+
+    /**
+     * Return an ArborX::Sphere from the std::pair<dealii::Point, Number>
+     * `v[i]`.
+     */
+    static Sphere
+    get(const std::vector<std::pair<dealii::Point<dim, Number>, 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<dealii::ArborXWrappers::SphereIntersectPredicate,
+                      PredicatesTag>
+  {
+    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<dealii::ArborXWrappers::SphereNearestPredicate,
+                      PredicatesTag>
+  {
+    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<dealii::ArborXWrappers::SphereIntersectPredicate,
+               PredicatesTag>::
+    size(const dealii::ArborXWrappers::SphereIntersectPredicate &sph_intersect)
+  {
+    return sph_intersect.size();
+  }
+
+
+
+  inline auto
+  AccessTraits<dealii::ArborXWrappers::SphereIntersectPredicate,
+               PredicatesTag>::
+    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<dealii::ArborXWrappers::SphereNearestPredicate, PredicatesTag>::
+    size(const dealii::ArborXWrappers::SphereNearestPredicate &sph_nearest)
+  {
+    return sph_nearest.size();
+  }
+
+
+
+  inline auto
+  AccessTraits<dealii::ArborXWrappers::SphereNearestPredicate, PredicatesTag>::
+    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
index c8b3d3473d3fc94e886e13b1734fb292d768440c..51e008f38d691a0e458cbe9db73ebcbfdfe7ad21 100644 (file)
@@ -146,6 +146,71 @@ namespace ArborXWrappers
   {
     return n_nearest_neighbors;
   }
+
+  // ------------------- SpherePredicate ------------------- //
+  template <int dim, typename Number>
+  SpherePredicate::SpherePredicate(
+    const std::vector<std::pair<dealii::Point<dim, Number>, 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<float>(dim_spheres[i].first[0]),
+            static_cast<float>(dim_spheres[i].first[1]),
+            dim == 2 ? 0.f : static_cast<float>(dim_spheres[i].first[2])),
+          static_cast<float>(dim_spheres[i].second)));
+      }
+  }
+
+
+
+  std::size_t
+  SpherePredicate::size() const
+  {
+    return spheres.size();
+  }
+
+
+
+  const std::pair<dealii::Point<3, float>, float> &
+  SpherePredicate::get(unsigned int i) const
+  {
+    return spheres[i];
+  }
+
+
+
+  template <int dim, typename Number>
+  SphereIntersectPredicate::SphereIntersectPredicate(
+    const std::vector<std::pair<dealii::Point<dim, Number>, Number>> &spheres)
+    : SpherePredicate(spheres)
+  {}
+
+
+
+  template <int dim, typename Number>
+  SphereNearestPredicate::SphereNearestPredicate(
+    const std::vector<std::pair<dealii::Point<dim, Number>, 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<float>(max_corner[1]),
              dim == 2 ? 0.f : static_cast<float>(max_corner[2])}};
   }
+
+
+
+  // ---------------------- Sphere Primitives AccessTraits ----------------- //
+  template <int dim, typename Number>
+  std::size_t
+  AccessTraits<std::vector<std::pair<dealii::Point<dim, Number>, Number>>,
+               PrimitivesTag>::
+    size(const std::vector<std::pair<dealii::Point<dim, Number>, Number>> &v)
+  {
+    return v.size();
+  }
+
+
+
+  template <int dim, typename Number>
+  Sphere
+  AccessTraits<std::vector<std::pair<dealii::Point<dim, Number>, Number>>,
+               PrimitivesTag>::
+    get(const std::vector<std::pair<dealii::Point<dim, Number>, 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<float>(v[i].first[0]),
+             static_cast<float>(v[i].first[1]),
+             dim == 2 ? 0 : static_cast<float>(v[i].first[2])},
+            static_cast<float>(v[i].second)};
+  }
 } // namespace ArborX
 
 // ----------------------- Instantiations --------------------//
index cc949b0f1285d471215316549154f2927d26c1e3..8cd79a4198b565a95d0d343f7138fbd1db0bb5d3 100644 (file)
@@ -31,6 +31,14 @@ for (DIM : DIMENSIONS; SCALAR : REAL_SCALARS)
       template BoundingBoxNearestPredicate::BoundingBoxNearestPredicate(
         const std::vector<dealii::BoundingBox<DIM, SCALAR>> &bb,
         const unsigned int n_nearest_neighbors);
+
+      template SphereIntersectPredicate::SphereIntersectPredicate(
+        const std::vector<std::pair<dealii::Point<DIM, SCALAR>, SCALAR>>
+          &spheres);
+      template SphereNearestPredicate::SphereNearestPredicate(
+        const std::vector<std::pair<dealii::Point<DIM, SCALAR>, 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<std::vector<dealii::Point<DIM, SCALAR>>,
                                    PrimitivesTag>;
+
       template struct AccessTraits<
         std::vector<dealii::BoundingBox<DIM, SCALAR>>,
         PrimitivesTag>;
+
+      template struct AccessTraits<
+        std::vector<std::pair<dealii::Point<DIM, SCALAR>, SCALAR>>,
+        PrimitivesTag>;
     \}
 #endif
   }

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.