add detail::indices_{of,array} to obtain indices of a mask

Change-Id: Id94ba16a00352a56bca915b0793b31bb08412bef

rigid_geometric_algebra/BUILD.bazel

@@ -31,6 +31,8 @@ cc_library(
         "detail/even.hpp",
         "detail/has_type.hpp",
         "detail/has_value.hpp",
+        "detail/indices_array.hpp",
+        "detail/indices_of.hpp",
         "detail/invoke_prioritized_overload.hpp",
         "detail/is_complete.hpp",
         "detail/is_defined.hpp",

rigid_geometric_algebra/blade.hpp

@@ -74,6 +74,7 @@ class blade
           detail::get_coefficient>,
       detail::derive_subtraction<blade<A, Is...>>
 {
+  // TODO replace
   static constexpr auto sorted_dimensions = [] {
     auto dims = std::array<std::size_t, sizeof...(Is)>{Is...};
     std::ranges::sort(dims);

rigid_geometric_algebra/blade_sum.hpp

@@ -2,15 +2,14 @@
 
 #include "rigid_geometric_algebra/canonical_type.hpp"
 #include "rigid_geometric_algebra/common_algebra_type.hpp"
+#include "rigid_geometric_algebra/detail/indices_array.hpp"
 #include "rigid_geometric_algebra/detail/type_list.hpp"
 #include "rigid_geometric_algebra/is_blade.hpp"
 #include "rigid_geometric_algebra/multivector.hpp"
 #include "rigid_geometric_algebra/sorted_canonical_blades.hpp"
 
-#include <algorithm>
 #include <array>
 #include <cstddef>
-#include <ranges>
 #include <tuple>
 #include <type_traits>
 #include <utility>
@@ -26,26 +25,14 @@ class blade_sum_fn
     template <detail::blade... Bs>
     static constexpr auto operator()(const std::tuple<Bs...>& bs) -> Canonical
     {
-      static constexpr auto mask =
-          std::array{std::is_same_v<canonical_type_t<Bs>, Canonical>...};
-      static constexpr auto count = std::ranges::count(mask, true);
-      static_assert(count >= 1);
+      static constexpr auto indices = detail::indices_array<std::array{
+          std::is_same_v<canonical_type_t<Bs>, Canonical>...}>;
 
-      static constexpr auto indices = []() consteval {
-        auto indices = std::array<std::size_t, count>{};
-        auto it = indices.begin();
-        for (const auto [i, convertible] :
-             std::views::zip(std::views::iota(0UZ), mask)) {
-          if (convertible) {
-            *it++ = i;
-          }
-        }
-        return indices;
-      }();
+      static_assert(indices.size() >= 1);
 
       return [&bs]<std::size_t... Is>(std::index_sequence<Is...>) {
         return (std::get<indices[Is]>(bs).canonical() + ...);
-      }(std::make_index_sequence<count>{});
+      }(std::make_index_sequence<indices.size()>{});
     }
   };
 

rigid_geometric_algebra/blade_type_from.hpp

@@ -3,9 +3,8 @@
 #include "rigid_geometric_algebra/algebra_type.hpp"
 #include "rigid_geometric_algebra/blade.hpp"
 #include "rigid_geometric_algebra/blade_ordering.hpp"
-#include "rigid_geometric_algebra/detail/contract.hpp"
+#include "rigid_geometric_algebra/detail/indices_array.hpp"
 
-#include <array>
 #include <cstddef>
 #include <utility>
 
@@ -21,25 +20,9 @@ struct blade_type_from
   template <std::size_t... Is>
   static constexpr auto impl(std::index_sequence<Is...>)
   {
-    using A = algebra_type_t<decltype(ord)>;
-
-    static constexpr auto dims = [] {
-      const auto& mask = ord.mask;
-      auto dims = std::array<std::size_t, sizeof...(Is)>{};
-      auto it = dims.begin();
-
-      for (auto index = 0UZ; index != mask.size(); ++index) {
-        if (mask.test(index)) {
-          detail::invariant(it != dims.end());
-          *it++ = index;
-        }
-      }
-
-      detail::postcondition(it == dims.end());
-      return dims;
-    }();
-
-    return blade<A, std::get<Is>(dims)...>{};
+    return blade<
+        algebra_type_t<decltype(ord)>,
+        std::get<Is>(detail::indices_array<ord.mask>)...>{};
   }
 
   using type = decltype(impl(std::make_index_sequence<ord.mask.count()>{}));

rigid_geometric_algebra/detail/indices_array.hpp

@@ -0,0 +1,25 @@
+#pragma once
+
+#include "rigid_geometric_algebra/detail/contract.hpp"
+#include "rigid_geometric_algebra/detail/indices_of.hpp"
+
+#include <algorithm>
+#include <array>
+#include <cstddef>
+#include <functional>
+
+namespace rigid_geometric_algebra::detail {
+
+/// create an array containing the indices of a mask
+///
+template <auto mask, auto pred = std::identity{}>
+inline constexpr auto indices_array = [] {
+  auto indices = std::array<std::size_t, std::ranges::count_if(mask, pred)>{};
+
+  auto it = detail::indices_of(mask, indices.begin());
+  detail::postcondition(it == indices.end());
+
+  return indices;
+}();
+
+}  // namespace rigid_geometric_algebra::detail

rigid_geometric_algebra/detail/indices_of.hpp

@@ -0,0 +1,38 @@
+#pragma once
+
+#include <concepts>
+#include <cstddef>
+#include <functional>
+#include <iterator>
+#include <ranges>
+
+namespace rigid_geometric_algebra::detail {
+
+/// determine indices of a range that satisfy a unary predicate
+///
+inline constexpr class
+{
+public:
+  template <
+      std::ranges::forward_range R,
+      std::weakly_incrementable O,
+      class UnaryPredicate = std::identity>
+    requires std::unsigned_integral<std::iter_value_t<O>> and
+             std::indirect_unary_predicate<
+                 UnaryPredicate,
+                 std::ranges::iterator_t<R>>
+  static constexpr auto
+  operator()(R&& rng, O dest, UnaryPredicate pred = {}) -> O
+  {
+    for (auto index = std::iter_value_t<O>{};
+         const auto& element : std::forward<R>(rng)) {
+      if (pred(element)) {
+        *dest++ = index;
+      }
+      ++index;
+    }
+    return dest;
+  }
+} indices_of{};
+
+}  // namespace rigid_geometric_algebra::detail

rigid_geometric_algebra/detail/structural_bitset.hpp

@@ -5,6 +5,7 @@
 #include <bit>
 #include <cstddef>
 #include <cstdint>
+#include <iterator>
 #include <type_traits>
 
 namespace rigid_geometric_algebra::detail {
@@ -18,6 +19,50 @@ struct structural_bitset
   static_assert(N <= 8);
   using value_type = std::uint8_t;
 
+  class const_iterator
+  {
+    const structural_bitset* parent_{};
+    structural_bitset::value_type pos_{};
+
+  public:
+    using iterator_category = std::forward_iterator_tag;
+    using value_type = bool;
+    using reference = bool;
+    using difference_type = std::ptrdiff_t;
+    using pointer = void;
+
+    constexpr const_iterator()
+    {
+      // https://stackoverflow.com/questions/28832492/why-do-iterators-need-to-be-default-constructible
+      detail::precondition(true);
+    }
+
+    constexpr const_iterator(
+        const structural_bitset& parent, structural_bitset::value_type pos)
+        : parent_{&parent}, pos_{pos}
+    {
+      detail::precondition(pos_ <= parent_->size());
+    }
+
+    constexpr auto operator*() const noexcept -> reference
+    {
+      return parent_->test(pos_);
+    }
+    constexpr auto operator++() & -> const_iterator&
+    {
+      ++pos_;
+      return *this;
+    }
+    constexpr auto operator++(int) & -> const_iterator
+    {
+      auto tmp = *this;
+      ++(*this);
+      return tmp;
+    }
+
+    friend auto operator==(const_iterator, const_iterator) -> bool = default;
+  };
+
   /// the number of bits the `structural_bitset` holds
   ///
   static constexpr auto size = std::integral_constant<std::size_t, N>{};
@@ -87,10 +132,27 @@ struct structural_bitset
     return value;
   }
 
+  /// iterators for range-based access
+  ///
+  /// @{
+
+  [[nodiscard]]
+  constexpr auto begin() const noexcept -> const_iterator
+  {
+    return {*this, 0};
+  }
+  [[nodiscard]]
+  constexpr auto end() const noexcept -> const_iterator
+  {
+    return {*this, size};
+  }
+
+  /// @}
+
   /// equality comparison
   ///
-  friend constexpr auto
-  operator==(structural_bitset, structural_bitset) noexcept -> bool = default;
+  friend auto
+  operator==(structural_bitset, structural_bitset) -> bool = default;
 };
 
 }  // namespace rigid_geometric_algebra::detail

rigid_geometric_algebra/detail/swaps_to_sorted_dimensions.hpp

@@ -46,7 +46,7 @@ inline constexpr class
 
 public:
   template <std::same_as<std::size_t>... Ts>
-  constexpr static auto operator()(Ts... is) -> std::size_t
+  static constexpr auto operator()(Ts... is) -> std::size_t
   {
     auto unsorted = std::array<std::size_t, sizeof...(is)>{is...};
 
@@ -62,7 +62,6 @@ inline constexpr class
         std::ranges::is_sorted(sorted, std::ranges::less_equal{}));
     return impl(std::ranges::subrange(unsorted), std::ranges::subrange(sorted));
   }
-
 } swaps_to_sorted_dimensions{};
 
 }  // namespace rigid_geometric_algebra::detail

rigid_geometric_algebra/detail/type_filter.hpp

@@ -1,8 +1,11 @@
 #pragma once
 
-#include "rigid_geometric_algebra/detail/type_concat.hpp"
+#include "rigid_geometric_algebra/detail/indices_array.hpp"
 
-#include <type_traits>
+#include <array>
+#include <cstddef>
+#include <tuple>
+#include <utility>
 namespace rigid_geometric_algebra::detail {
 
 /// filter types in a type list
@@ -16,22 +19,23 @@ template <class L, auto pred>
 struct type_filter
 {};
 
-template <class L, auto pred>
-using type_filter_t = typename type_filter<L, pred>::type;
-
-template <template <class...> class list, auto pred>
-struct type_filter<list<>, pred>
+template <template <class...> class list, class... Ts, auto pred>
+struct type_filter<list<Ts...>, pred>
 {
-  using type = list<>;
+  static constexpr auto indices = detail::indices_array<
+      std::array<bool, sizeof...(Ts)>{pred.template operator()<Ts>()...}>;
+
+  using type =
+      decltype([]<std::size_t... Is>(std::index_sequence<Is...>)
+                   -> list<std::tuple_element_t<
+                       indices[Is],
+                       std::tuple<Ts...>>...> {
+        return {};
+      }(std::make_index_sequence<indices.size()>{}));
 };
 
-template <template <class...> class list, class T0, class... Ts, auto pred>
-struct type_filter<list<T0, Ts...>, pred>
-{
-  using type = detail::type_concat_t<
-      std::conditional_t<pred.template operator()<T0>(), list<T0>, list<>>,
-      type_filter_t<list<Ts...>, pred>>;
-};
+template <class L, auto pred>
+using type_filter_t = typename type_filter<L, pred>::type;
 
 /// @}