Grid implementation that addresses some issues found in the local

navigation (grid navigation must be possible without reference
surface). This PR provides:

 - 3D grids
 - Grid local coordinates for a global to local transform
 - Memory serialization: standalone grids and grids in grid collections
 - Implement populators, serializers and unittests, except for
   irregular attach populator.

Due to the 'operator==' not being a __device__ function, the
complete and attach populator unittests cannot pass the CI.

The PR also splits the geometric shape information (cylinder,
rectangle etc) from the navigation linking and boundary values, so
that they can be reused for the grid geometry definition. With that
comes the stringent use of local coordinates in the intersections
and some refactoring (renaming etc).

core/CMakeLists.txt

@@ -56,21 +56,21 @@ detray_add_library( detray_core core
    "include/detray/intersection/line_intersector.hpp"
    "include/detray/intersection/plane_intersector.hpp"
    # masks include(s)
-   "include/detray/masks/annulus2.hpp"
-   "include/detray/masks/cylinder3.hpp"
+   "include/detray/masks/annulus2D.hpp"
+   "include/detray/masks/cuboid3D.hpp"
+   "include/detray/masks/cylinder3D.hpp"
    "include/detray/masks/line.hpp"
    "include/detray/masks/masks.hpp"
-   "include/detray/masks/masks_base.hpp"
-   "include/detray/masks/rectangle2.hpp"
-   "include/detray/masks/ring2.hpp"
-   "include/detray/masks/single3.hpp"
-   "include/detray/masks/trapezoid2.hpp"
+   "include/detray/masks/rectangle2D.hpp"
+   "include/detray/masks/ring2D.hpp"
+   "include/detray/masks/single3D.hpp"
+   "include/detray/masks/trapezoid2D.hpp"
    "include/detray/masks/unmasked.hpp"
    # materials include(s)
-   "include/detray/materials/mixture.hpp" 
+   "include/detray/materials/mixture.hpp"
    "include/detray/materials/material.hpp"
-   "include/detray/materials/material_rod.hpp"   
-   "include/detray/materials/material_slab.hpp"   
+   "include/detray/materials/material_rod.hpp"
+   "include/detray/materials/material_slab.hpp"
    "include/detray/materials/predefined_materials.hpp"
    # propagator include(s)
    "include/detray/propagator/aborters.hpp"
@@ -89,6 +89,18 @@ detray_add_library( detray_core core
    "include/detray/surface_finders/brute_force_finder.hpp"
    "include/detray/surface_finders/grid2_finder.hpp"
    "include/detray/surface_finders/neighborhood_kernel.hpp"
+   # surface finder include(s)
+   "include/detray/surface_finders/grid/detail/axis_binning.hpp"
+   "include/detray/surface_finders/grid/detail/axis_helpers.hpp"
+   "include/detray/surface_finders/grid/detail/axis_shape.hpp"
+   "include/detray/surface_finders/grid/detail/grid_helpers.hpp"
+   "include/detray/surface_finders/grid/detail/populator_impl.hpp"
+   "include/detray/surface_finders/grid/detail/serializer_impl.hpp"
+   "include/detray/surface_finders/grid/axis.hpp"
+   "include/detray/surface_finders/grid/grid_builder.hpp"
+   "include/detray/surface_finders/grid/grid.hpp"
+   "include/detray/surface_finders/grid/populator.hpp"
+   "include/detray/surface_finders/grid/serializer.hpp"
    # tools include(s)
    "include/detray/tools/associator.hpp"
    "include/detray/tools/bin_association.hpp"

core/include/detray/coordinates/cylindrical2.hpp

@@ -16,7 +16,7 @@
 
 namespace detray {
 
-/** Local frame projection into a polar coordinate frame
+/** Local frame projection into a 2D cylindrical coordinate frame
  */
 template <typename transform3_t>
 struct cylindrical2 : public coordinate_base<cylindrical2, transform3_t> {
@@ -60,7 +60,7 @@ struct cylindrical2 : public coordinate_base<cylindrical2, transform3_t> {
     DETRAY_HOST_DEVICE inline point3 local_to_global(
         const transform3_t &trf, const mask_t &mask, const point2 &p,
         const vector3 & /*d*/) const {
-        const scalar_type r = mask.radius();
+        const scalar_type r = mask[0];
         const scalar_type phi = p[0] / r;
         const scalar_type x = r * std::cos(phi);
         const scalar_type y = r * std::sin(phi);

core/include/detray/core/type_registry.hpp

@@ -79,6 +79,7 @@ class registry_base<ID, true, registered_types...> {
     template <typename object_t>
     struct get_index {
         static constexpr ID value = get_id<object_t>();
+        DETRAY_HOST_DEVICE
         constexpr bool operator()() noexcept { return is_valid(value); }
     };
 
@@ -96,7 +97,7 @@ class registry_base<ID, true, registered_types...> {
                 is_valid(ref_idx),
                 "Index out of range: Please make sure that indices and type "
                 "enums match the number of types in container.");
-            return static_cast<ID>(index);
+            return static_cast<ID>(ref_idx);
         }
         if constexpr (ref_idx < sizeof...(registered_types) - 1) {
             return to_id<ref_idx + 1>(index);

core/include/detray/definitions/containers.hpp

@@ -1,20 +1,32 @@
 /** Detray library, part of the ACTS project (R&D line)
  *
- * (c) 2020 CERN for the benefit of the ACTS project
+ * (c) 2020-2022 CERN for the benefit of the ACTS project
  *
  * Mozilla Public License Version 2.0
  */
 
 #pragma once
 
+// Project include(s)
+#include "detray/definitions/qualifiers.hpp"
+#include "detray/utils/type_traits.hpp"
+
+// Vecmem include(s)
+#include <vecmem/containers/device_vector.hpp>
+#include <vecmem/containers/jagged_device_vector.hpp>
+#include <vecmem/containers/jagged_vector.hpp>
+#include <vecmem/containers/vector.hpp>
+
+// Thrust include(s)
+#include <thrust/tuple.h>
+
+// System include(s)
 #include <array>
 #include <map>
 #include <tuple>
-#include <vecmem/containers/vector.hpp>
+#include <type_traits>
 #include <vector>
 
-#include "vecmem/containers/jagged_vector.hpp"
-
 namespace detray {
 template <typename value_t, std::size_t kDIM>
 using darray = std::array<value_t, kDIM>;
@@ -31,4 +43,84 @@ using dmap = std::map<key_t, value_t>;
 template <class... types>
 using dtuple = std::tuple<types...>;
 
+/// @brief Bundle container type definitions
+template <template <typename...> class vector_t = dvector,
+          template <typename...> class tuple_t = dtuple,
+          template <typename, std::size_t> class array_t = darray,
+          template <typename...> class jagged_vector_t = djagged_vector,
+          template <typename, typename> class map_t = dmap>
+struct container_types {
+    template <typename T>
+    using vector_type = vector_t<T>;
+
+    template <class... T>
+    using tuple_type = tuple_t<T...>;
+
+    template <typename T, std::size_t kDIM>
+    using array_type = array_t<T, kDIM>;
+
+    template <typename T>
+    using jagged_vector_type = jagged_vector_t<T>;
+
+    template <typename K, typename T>
+    using map_type = map_t<K, T>;
+};
+
+/// Defining some common types
+using host_container_types = container_types<>;
+using device_container_types =
+    container_types<vecmem::device_vector, thrust::tuple, darray,
+                    vecmem::jagged_device_vector>;
+
+// Views for types that aggregate containers
+
+/// Empty view type for inheritance template resolution
+struct dbase_view {};
+
+/// Tag a vecmem view as @c dbase_view , so that it becomes recognizable in
+/// detray as a view type
+template <typename value_t, template <typename> class view_t>
+struct view_wrapper : public dbase_view {
+    /// The vecmem data view
+    view_t<value_t> m_view{};
+
+    /// Default constructor
+    view_wrapper() = default;
+
+    /// Conversion operator from a view of the same value type
+    DETRAY_HOST
+    view_wrapper(view_t<value_t>&& view) : m_view{view} {}
+};
+
+/// Container view helper that aggregates multiple views and performs
+/// compile-time checks.
+template <bool /*check value*/, typename... view_ts>
+class dmulti_view_helper {};
+
+/// In case the checks fail
+template <typename... view_ts>
+class dmulti_view_helper<false, view_ts...> {};
+
+/// General view type that aggregates vecmem based view implementations
+template <typename... view_ts>
+struct dmulti_view_helper<true, view_ts...> : public dbase_view {
+    dtuple<view_ts...> m_views;
+
+    dmulti_view_helper() = default;
+
+    /// Tie multiple views together
+    DETRAY_HOST
+    dmulti_view_helper(view_ts&&... views) { m_views = std::tie(views...); }
+};
+
+/// The detray container view exists, if all contained view types also derive
+/// from @c dbase_view.
+template <typename... view_ts>
+using dmulti_view = dmulti_view_helper<
+    std::conjunction_v<std::is_base_of<dbase_view, view_ts>...>, view_ts...>;
+
+/// Detray view for vector containers
+template <typename value_t>
+using dvector_view = view_wrapper<value_t, vecmem::data::vector_view>;
+
 }  // namespace detray

core/include/detray/definitions/detail/accessor.hpp

@@ -7,6 +7,7 @@
 #pragma once
 
 #if defined(__CUDACC__)
+#include <thrust/binary_search.h>
 #include <thrust/execution_policy.h>
 #include <thrust/find.h>
 #include <thrust/sort.h>
@@ -36,6 +37,24 @@ namespace vtuple = std;
 
 namespace detail {
 
+/// Extract the first type from a parameter pack.
+template <typename first_t = void, typename... other_t>
+struct first_type {
+    using type = first_t;
+};
+
+template <typename... TYPES>
+using first_t = typename first_type<TYPES...>::type;
+
+/// Extract the first value from an index pack.
+template <std::size_t first_t, std::size_t... other_t>
+struct first_idx {
+    static constexpr std::size_t value = first_t;
+};
+
+template <std::size_t... TYPES>
+inline constexpr std::size_t first_idx_v = first_idx<TYPES...>::value;
+
 namespace {
 struct empty_type {};
 
@@ -254,5 +273,27 @@ DETRAY_HOST_DEVICE auto find_if(RandomIt first, RandomIt last, Predicate comp) {
 #endif
 }
 
+/// @brief lower_bound implementation for device
+template <class ForwardIt, typename Value>
+DETRAY_HOST_DEVICE auto lower_bound(ForwardIt first, ForwardIt last,
+                                    Value value) {
+#if defined(__CUDACC__)
+    return thrust::lower_bound(thrust::seq, first, last, value);
+#elif !defined(__CUDACC__)
+    return std::lower_bound(first, last, value);
+#endif
+}
+
+/// @brief upper_bound implementation for device
+template <class ForwardIt, typename Value>
+DETRAY_HOST_DEVICE auto upper_bound(ForwardIt first, ForwardIt last,
+                                    Value value) {
+#if defined(__CUDACC__)
+    return thrust::upper_bound(thrust::seq, first, last, value);
+#elif !defined(__CUDACC__)
+    return std::upper_bound(first, last, value);
+#endif
+}
+
 }  // namespace detail
 }  // namespace detray

core/include/detray/definitions/indexing.hpp

@@ -1,6 +1,6 @@
 /** Detray library, part of the ACTS project (R&D line)
  *
- * (c) 2020 CERN for the benefit of the ACTS project
+ * (c) 2020-2022 CERN for the benefit of the ACTS project
  *
  * Mozilla Public License Version 2.0
  */
@@ -18,6 +18,32 @@ dindex constexpr dindex_invalid = std::numeric_limits<dindex>::max();
 using dindex_range = darray<dindex, 2>;
 using dindex_sequence = dvector<dindex>;
 
+/// @brief Simple multi-index structure
+///
+/// @tparam DIM number of indices that are held by this type
+/// @tparam value_t type of indices
+template <typename value_t = dindex, std::size_t DIM = 3>
+struct dmulti_index {
+    std::array<value_t, DIM> indices{};
+
+    /// Elementwise access.
+    DETRAY_HOST_DEVICE
+    inline auto operator[](const std::size_t i) -> value_t& {
+        return indices[i];
+    }
+    DETRAY_HOST_DEVICE
+    inline auto operator[](const std::size_t i) const -> const value_t& {
+        return indices[i];
+    }
+
+    /// Equality operator @returns true if all bin indices match.
+    DETRAY_HOST_DEVICE
+    inline auto operator==(const dmulti_index<value_t, DIM>& rhs) const
+        -> bool {
+        return (indices == rhs.indices);
+    }
+};
+
 /// @brief Ties an object type and an index into a container together.
 ///
 /// @tparam id_type Represents the type of object that is being indexed

core/include/detray/grids/axis.hpp

@@ -43,14 +43,14 @@ struct regular {
     /** Defualt constructor for dummy axis **/
     DETRAY_HOST_DEVICE
     regular()
-        : n_bins(invalid_value<dindex>()),
+        : n_bins(detail::invalid_value<dindex>()),
           min(0.),
           max(static_cast<scalar>(n_bins)) {}
 
     /** Constructor with vecmem memory resource **/
     DETRAY_HOST
     regular(vecmem::memory_resource & /*resource*/)
-        : n_bins(invalid_value<dindex>()),
+        : n_bins(detail::invalid_value<dindex>()),
           min(0.),
           max(static_cast<scalar>(n_bins)) {}
 
@@ -225,14 +225,14 @@ struct circular {
     /** Defualt constructor for dummy axis **/
     DETRAY_HOST_DEVICE
     circular()
-        : n_bins(invalid_value<dindex>()),
+        : n_bins(detail::invalid_value<dindex>()),
           min(0.),
           max(static_cast<scalar>(n_bins)) {}
 
     /** Constructor with vecmem memory resource **/
     DETRAY_HOST
     circular(vecmem::memory_resource & /*resource*/)
-        : n_bins(invalid_value<dindex>()),
+        : n_bins(detail::invalid_value<dindex>()),
           min(0.),
           max(static_cast<scalar>(n_bins)) {}
 
@@ -432,15 +432,15 @@ struct irregular {
     /** Defualt constructor for dummy axis **/
     DETRAY_HOST_DEVICE
     irregular()
-        : n_bins(invalid_value<dindex>()),
+        : n_bins(detail::invalid_value<dindex>()),
           min(0.),
           max(static_cast<scalar>(n_bins)),
           boundaries({}) {}
 
     /** Constructor with vecmem memory resource **/
     DETRAY_HOST
     irregular(vecmem::memory_resource &resource)
-        : n_bins(invalid_value<dindex>()),
+        : n_bins(detail::invalid_value<dindex>()),
           min(0.),
           max(static_cast<scalar>(n_bins)),
           boundaries(&resource) {}
@@ -628,7 +628,7 @@ struct axis_data<axis_t, scalar_t,
     /// Construct a const data object from a non-const one
     template <
         typename other_scalar_t,
-        std::enable_if_t<details::is_same_nc<scalar_t, other_scalar_t>::value,
+        std::enable_if_t<detail::is_same_nc<scalar_t, other_scalar_t>::value,
                          bool> = true>
     DETRAY_HOST_DEVICE axis_data(
         const axis_data<axis_t, other_scalar_t, void> &parent)
@@ -654,7 +654,7 @@ struct axis_data<axis_t, scalar_t,
     /// Construct a const data object from a non-const one
     template <
         typename other_scalar_t,
-        std::enable_if_t<details::is_same_nc<scalar_t, other_scalar_t>::value,
+        std::enable_if_t<detail::is_same_nc<scalar_t, other_scalar_t>::value,
                          bool> = true>
     DETRAY_HOST_DEVICE axis_data(
         const axis_data<axis_t, other_scalar_t, void> &parent)