use canonical dimension order in blades

Change-Id: Id8730a7b484435b72188eaa479b443c795a55121

BUILD.bazel

@@ -18,6 +18,8 @@ cc_library(
     name = "geometry",
     srcs = [
         "geometry/src/algebra.hpp",
+        "geometry/src/detail/sort_dimensions.hpp",
+        "geometry/src/detail/unique_sorted.hpp",
     ],
     hdrs = ["geometry/geometry.hpp"],
     visibility = ["//visibility:public"],

geometry/src/algebra.hpp

@@ -1,38 +1,78 @@
 #pragma once
 
+#include "geometry/src/detail/sort_dimensions.hpp"
+#include "geometry/src/detail/unique_sorted.hpp"
+
 #include <cstddef>
 #include <type_traits>
 
 namespace geometry {
 
+/// represents a projective geometric algebra
+/// @tparam S scalar type
+/// @tparam N number of euclidean geometry dimensions
+///
+/// Models an N-dimensional projective geometric algebra denoted as
+/// @f[
+/// P(S^*_{N,0,1})
+/// @f]
+///
 template <class S, std::size_t N>
 struct algebra
 {
+  /// scalar type
+  ///
   using scalar_type = S;
 
+  /// blade
+  /// @tparam Is dimensions
+  ///
+  /// Defines a blade for the given algebra.
+  ///
   template <std::size_t... Is>
   struct blade
   {
     static_assert(sizeof...(Is) <= N + 1, "number of `Is` exceeds rank");
-    static_assert(((Is <= N) and ...), "`Is` exceeds rank");
-    // TODO check Is.. values are unique
+    static_assert(((Is <= N) and ...), "`Is` value exceeds rank");
+    static_assert(
+        detail::unique_sorted(Is...), "`Is` values must be unique and sorted");
 
+    /// algebra type
+    ///
     using algebra_type = algebra;
+
+    /// scalar type
+    ///
     using scalar_type = typename algebra::scalar_type;
 
+    /// grade
+    ///
     static constexpr auto grade = sizeof...(Is);
 
+    /// k-vector coefficient
+    ///
     scalar_type coefficient{};
 
+    /// construct a zero blade
+    ///
     blade() = default;
 
+    /// construct a blade
+    /// @param a coefficient
+    /// @note this constructor is explicit if and only if grade != 0
+    ///
+    /// @{
     template <std::size_t M = grade, std::enable_if_t<M == 0, bool> = true>
     constexpr blade(scalar_type a) : coefficient{a}
     {}
     template <std::size_t M = grade, std::enable_if_t<M != 0, bool> = true>
     constexpr explicit blade(scalar_type a) : coefficient{a}
     {}
+    /// @}
 
+    /// comparison operators
+    ///
+    /// @{
     [[nodiscard]]
     friend constexpr auto
     operator==(blade x, blade y) noexcept -> bool
@@ -69,14 +109,20 @@ struct algebra
     {
       return x.coefficient >= y.coefficient;
     }
+    /// @}
 
+    /// unary negation
+    ///
     [[nodiscard]]
     friend constexpr auto
     operator-(blade x) -> blade
     {
       return blade{-x.coefficient};
     }
 
+    /// scalar multiplication
+    ///
+    /// @{
     [[nodiscard]]
     friend constexpr auto
     operator*(scalar_type a, blade x) -> blade
@@ -89,10 +135,48 @@ struct algebra
     {
       return a * x;
     }
+    /// @}
   };
 
+private:
+  template <std::size_t... Is>
+  static auto rebind(std::index_sequence<Is...>) -> blade<Is...>;
+
+  template <class Seq>
+  using rebind_t = decltype(rebind(Seq{}));
+
+  template <std::size_t... Is>
+  using reified_blade_t =
+      rebind_t<detail::sort_dimensions_t<std::index_sequence<Is...>>>;
+
+  template <std::size_t... Is>
+  static constexpr auto reified_blade_coefficient_v =
+      detail::sort_dimensions_swap_count_v<std::index_sequence<Is...>> % 2 == 0
+          ? scalar_type{1}
+          : -scalar_type{1};
+
+public:
+  /// unit blade
+  /// @tparam Is dimensions
+  ///
+  /// Obtains the unit blade corresponding to dimensions `Is...` for the given
+  /// algebra. `Is...` will be reified if not provided in canonical form (i.e.
+  /// in non-decreasing order).
+  ///
+  /// ~~~{.cpp}
+  /// static_assert(e<1, 2> == -e<2, 1>);
+  /// ~~~
+  ///
+  /// A unit blade specified with zero dimensions (i.e. `e<>`) is semantically
+  /// equivalent to a scalar.
+  ///
+  /// ~~~{.cpp}
+  /// static_assert(2 == 2*e<>);
+  /// ~~~
+  ///
   template <std::size_t... Is>
-  static constexpr auto e = blade<Is...>{scalar_type{1}};
+  static constexpr auto e =
+      reified_blade_t<Is...>{reified_blade_coefficient_v<Is...>};
 };
 
 }  // namespace geometry

geometry/src/detail/sort_dimensions.hpp

@@ -0,0 +1,68 @@
+#pragma once
+
+#include <algorithm>
+#include <cstddef>
+#include <type_traits>
+#include <utility>
+
+namespace geometry::detail {
+
+template <class InSeq, class OutSeq = std::index_sequence<>>
+struct sort_swap_first
+{
+  using type = OutSeq;
+};
+template <std::size_t I0, std::size_t... Js>
+struct sort_swap_first<std::index_sequence<I0>, std::index_sequence<Js...>>
+{
+  using type = std::index_sequence<Js..., I0>;
+};
+template <std::size_t I0, std::size_t I1, std::size_t... Is, std::size_t... Js>
+struct sort_swap_first<
+    std::index_sequence<I0, I1, Is...>,
+    std::index_sequence<Js...>>
+    : std::conditional_t<
+          (I0 > I1),
+          sort_swap_first<std::index_sequence<>,
+                          std::index_sequence<Js..., I1, I0, Is...>>,
+          sort_swap_first<std::index_sequence<I1, Is...>,
+                          std::index_sequence<Js..., I0>>>
+{};
+
+template <class Seq>
+using sort_swap_first_t = typename sort_swap_first<Seq>::type;
+
+template <class Seq, std::size_t Swaps = 0, class = void>
+struct sort_swap : sort_swap<sort_swap_first_t<Seq>, Swaps + 1>
+{};
+template <std::size_t... Is, std::size_t Swaps>
+struct sort_swap<
+    std::index_sequence<Is...>,
+    Swaps,
+    std::enable_if_t<std::is_same_v<
+        std::index_sequence<Is...>,
+        sort_swap_first_t<std::index_sequence<Is...>>>>>
+{
+  using type = std::index_sequence<Is...>;
+  static constexpr auto swap_count = Swaps;
+};
+
+template <class Seq>
+struct sort_dimensions
+{
+  using type = typename sort_swap<Seq>::type;
+};
+
+template <class Seq>
+using sort_dimensions_t = typename sort_dimensions<Seq>::type;
+
+template <class Seq>
+struct sort_dimensions_swap_count
+    : std::integral_constant<std::size_t, sort_swap<Seq>::swap_count>
+{};
+
+template <class Seq>
+inline constexpr auto sort_dimensions_swap_count_v =
+    sort_dimensions_swap_count<Seq>::value;
+
+}  // namespace geometry::detail

geometry/src/detail/unique_sorted.hpp

@@ -0,0 +1,39 @@
+#pragma once
+
+#include <array>
+#include <cstddef>
+#include <utility>
+
+namespace geometry::detail {
+
+inline constexpr class
+{
+  template <class T>
+  static constexpr auto impl(std::index_sequence<>, const std::array<T, 1>&)
+  {
+    return true;
+  }
+  template <
+      std::size_t... Is,
+      class T,
+      class = std::enable_if_t<(sizeof...(Is) != 0)>>
+  static constexpr auto
+  impl(std::index_sequence<Is...>,
+       const std::array<T, sizeof...(Is) + 1>& values)
+  {
+    return ((std::get<Is>(values) < std::get<Is + 1>(values)) and ...);
+  }
+
+public:
+  // 0-args
+  constexpr auto operator()() const { return true; }
+  // 1+-args
+  template <class... Ts>
+  constexpr auto operator()(std::size_t t0, Ts... ts) const
+  {
+    static_assert((std::is_same_v<std::size_t, Ts> and ...));
+    return impl(std::index_sequence_for<Ts...>{}, std::array{t0, ts...});
+  }
+} unique_sorted{};
+
+}  // namespace geometry::detail

test/algebra_test.cpp

@@ -52,4 +52,16 @@ auto main() -> int
             gt(ei, 0 * ei) and  //
             ge(ei, 0 * ei));
       };
+
+  "blade uses canonical dimensions"_ctest = [] {
+    return expect(
+        eq(-e<1, 0>, e<0, 1>) and        //
+        eq(-e<2, 0>, e<0, 2>) and        //
+        eq(-e<2, 1>, e<1, 2>) and        //
+        eq(-e<0, 2, 1>, e<0, 1, 2>) and  //
+        eq(e<2, 0, 1>, e<0, 1, 2>) and   //
+        eq(-e<1, 0, 2>, e<0, 1, 2>) and  //
+        eq(e<1, 2, 0>, e<0, 1, 2>) and   //
+        eq(-e<2, 1, 0>, e<0, 1, 2>));
+  };
 }