diff --git a/geometry/BUILD.bazel b/geometry/BUILD.bazel
index 729589a..c756d5b 100644
--- a/geometry/BUILD.bazel
+++ b/geometry/BUILD.bazel
@@ -41,6 +41,7 @@ cc_library(
     name = "geometry",
     srcs = [
         "src/algebra.hpp",
+        "src/algebra_core.hpp",
         "src/algebra_type.hpp",
         "src/common_algebra_type.hpp",
         "src/detail/blade_list.hpp",
@@ -58,6 +59,7 @@ cc_library(
         "src/exterior_product.hpp",
         "src/geometric_product.hpp",
         "src/get.hpp",
+        "src/interior_product.hpp",
         "src/multivector_for.hpp",
         "src/regressive_product.hpp",
         "src/sum.hpp",
diff --git a/geometry/geometry.hpp b/geometry/geometry.hpp
index 874b779..cb6d38f 100644
--- a/geometry/geometry.hpp
+++ b/geometry/geometry.hpp
@@ -8,6 +8,7 @@
 #include "geometry/src/exterior_product.hpp"
 #include "geometry/src/geometric_product.hpp"
 #include "geometry/src/get.hpp"
+#include "geometry/src/interior_product.hpp"
 #include "geometry/src/multivector_for.hpp"
 #include "geometry/src/regressive_product.hpp"
 #include "geometry/src/sum.hpp"
diff --git a/geometry/src/algebra.hpp b/geometry/src/algebra.hpp
index b1eb0e7..3a82572 100644
--- a/geometry/src/algebra.hpp
+++ b/geometry/src/algebra.hpp
@@ -1,17 +1,10 @@
 #pragma once
 
+#include "geometry/src/algebra_core.hpp"
 #include "geometry/src/detail/blade_list.hpp"
-#include "geometry/src/detail/contract_dimensions.hpp"
-#include "geometry/src/detail/ordered.hpp"
-#include "geometry/src/detail/strictly_increasing.hpp"
-#include "geometry/src/get.hpp"
 #include "geometry/type_metaprogramming.hpp"
 
 #include <cstddef>
-#include <ostream>
-#include <tuple>
-#include <type_traits>
-#include <utility>
 
 namespace geometry {
 
@@ -25,105 +18,28 @@ namespace geometry {
 /// P(S^*_{N,0,1})
 /// @f]
 ///
+/// @note Defines a project geometric algebra, including common constants.
+///
 template <class S, std::size_t N, bool Dual = true>
-struct algebra
+struct algebra : algebra_core<S, N, Dual>
 {
-  /// scalar type
+  /// algebra core type
   ///
-  using scalar_type = S;
+  using algebra_core_type = typename algebra::algebra_core;
 
-  /// algebra dimension
+  /// registration of algebra type for associated core type
+  /// see declaration in algebra_core
   ///
-  static constexpr auto dimension = N + 1;
+  constexpr friend auto algebra_for(algebra_core_type) { return algebra{}; }
 
   /// dual algebra type
   ///
   using dual_algebra_type = algebra<S, N, not Dual>;
 
-  /// blade
-  /// @tparam Is dimensions
-  ///
-  /// Defines a blade for the given algebra.
-  ///
-  template <std::size_t... Is>
-  struct blade;
-
-  /// multivector
-  /// @tparam Bs blade types
-  ///
-  /// A linear combination of blades.
-  ///
-  /// @note a multivector may be composed of only blades with the same grade
-  ///
-  template <class... Bs>
-  struct multivector;
-
-  /// determines if a type is a blade
-  /// @tparam T type
-  ///
-  /// @{
-  template <class T>
-  struct is_blade : std::false_type
-  {};
-  template <std::size_t... Is>
-  struct is_blade<blade<Is...>> : std::true_type
-  {};
-  template <class T>
-  static constexpr auto is_blade_v = is_blade<T>::value;
-  /// @}
-
-  /// determines if a type is a multivector
-  /// @tparam T type
-  ///
-  /// @{
-  template <class T>
-  struct is_multivector : std::false_type
-  {};
-  template <class... Bs>
-  struct is_multivector<multivector<Bs...>> : std::true_type
-  {};
-  template <class T>
-  static constexpr auto is_multivector_v = is_multivector<T>::value;
-  /// @}
-
   /// list of blades in the algebra
   ///
-  using blade_list_type = detail::blade_list_t<blade>;
-
-  /// blade product
-  ///
-  /// @{
-  template <std::size_t... Is>
-  struct reified_blade
-  {
-    using type = tmp::convert_to_sequence_t<
-        detail::contract_dimensions_t<
-            detail::contraction_map::projective,
-            tmp::sort_t<tmp::index_constant_list<Is...>>>,
-        blade>;
-  };
-  template <std::size_t... Is>
-  using reified_blade_t = typename reified_blade<Is...>::type;
-  /// @}
-
-  /// coefficient for blade product
-  ///
-  /// @{
-  template <std::size_t... Is>
-  struct reified_blade_coefficient
-      : std::integral_constant<
-            int,
-            (tmp::sort_swap_count_v<tmp::index_constant_list<Is...>> % 2 == 0
-                 ? 1
-                 : -1) *
-                detail::contract_dimensions_coefficient_v<
-                    detail::contraction_map::projective,
-                    tmp::sort_t<tmp::index_constant_list<Is...>>>>
-  {};
-  template <std::size_t... Is>
-  static constexpr auto reified_blade_coefficient_v =
-      reified_blade_coefficient<Is...>::value;
-  /// @}
+  using blade_list_type =
+      detail::blade_list_t<algebra_core_type::template blade>;
 
   /// unit blade
   /// @tparam Is dimensions
@@ -147,282 +63,14 @@ struct algebra
   /// for `i == 0` and `1` for `i != 0`.
   ///
   template <std::size_t... Is>
-  static constexpr auto e =
-      reified_blade_t<Is...>{scalar_type{reified_blade_coefficient_v<Is...>}};
+  static constexpr auto e = typename algebra_core_type::
+      template reified_blade_t<Is...>{typename algebra_core_type::scalar_type{
+          algebra_core_type::template reified_blade_coefficient_v<Is...>}};
 
   /// pseudoscalar
   ///
-  static constexpr auto i = tmp::last_t<blade_list_type>{scalar_type{1}};
-
-  template <std::size_t... Is>
-  struct blade
-  {
-    static_assert(sizeof...(Is) <= N + 1, "number of `Is` exceeds rank");
-    static_assert(((Is <= N) and ...), "`Is` value exceeds rank");
-    static_assert(
-        detail::strictly_increasing(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
-    {
-      return x.coefficient == y.coefficient;
-    }
-    [[nodiscard]]
-    friend constexpr auto
-    operator!=(blade x, blade y) noexcept -> bool
-    {
-      return x.coefficient != y.coefficient;
-    }
-    [[nodiscard]]
-    friend constexpr auto
-    operator<(blade x, blade y) noexcept -> bool
-    {
-      return x.coefficient < y.coefficient;
-    }
-    [[nodiscard]]
-    friend constexpr auto
-    operator>(blade x, blade y) noexcept -> bool
-    {
-      return x.coefficient > y.coefficient;
-    }
-    [[nodiscard]]
-    friend constexpr auto
-    operator<=(blade x, blade y) noexcept -> bool
-    {
-      return x.coefficient <= y.coefficient;
-    }
-    [[nodiscard]]
-    friend constexpr auto
-    operator>=(blade x, blade y) noexcept -> bool
-    {
-      return x.coefficient >= y.coefficient;
-    }
-    /// @}
-
-    /// unary negation
-    ///
-    [[nodiscard]]
-    friend constexpr auto
-    operator-(blade x) -> blade
-    {
-      return blade{-x.coefficient};
-    }
-
-    /// addition
-    ///
-    /// @{
-    [[nodiscard]]
-    friend constexpr auto
-    operator+(blade x, blade y) -> blade
-    {
-      return blade{x.coefficient + y.coefficient};
-    }
-    /// @}
-
-    /// compound assignment operators
-    ///
-    /// @{
-    friend constexpr auto& operator+=(blade& x, blade y) { return x = x + y; }
-    friend constexpr auto& operator-=(blade& x, blade y) { return x = x - y; }
-    /// @}
-
-    /// geometric product
-    ///
-    template <std::size_t... Js>
-    [[nodiscard]]
-    friend constexpr auto
-    operator*(blade x, blade<Js...> y)
-    {
-      constexpr auto unit_blade_coeff =
-          reified_blade_coefficient_v<Is..., Js...>;
-      using B = reified_blade_t<Is..., Js...>;
-
-      if constexpr (unit_blade_coeff == 0) {
-        return B{};
-      } else {
-        const auto coeff = x.coefficient * y.coefficient;
-        return B{unit_blade_coeff == 1 ? coeff : -coeff};
-      }
-    }
-
-    /// stream insertion
-    ///
-    friend auto operator<<(std::ostream& os, blade x) -> std::enable_if_t<
-        (N < 10),
-        decltype(os << std::declval<const scalar_type&>(), os)>
-    {
-      if (sizeof...(Is) == 0) {
-        os << x.coefficient;
-        return os;
-      }
-
-      if (x.coefficient != scalar_type{1}) {
-        os << x.coefficient << "*";
-      }
-
-      os << "e";
-      std::ignore = ((os << Is, true) and ...);
-
-      return os;
-    }
-  };
-
-  template <class... Bs>
-  struct multivector : Bs...
-  {
-    static_assert(
-        (is_blade_v<Bs> and ...), "`Bs` must be a specialization of blade");
-    static_assert(
-        tmp::same_v<typename Bs::scalar_type...>,
-        "`Bs` must have the same scalar type");
-    static_assert(
-        detail::strictly_increasing(detail::ordered<Bs>{}...),
-        "`Bs` must be lexicographically sorted");
-
-    /// algebra type
-    ///
-    using algebra_type = algebra;
-
-    /// determine if a blade is contained in this multivector
-    ///
-    template <class B>
-    static constexpr auto contains = (std::is_same_v<B, Bs> or ...);
-
-    /// construct a zero multivector
-    ///
-    multivector() = default;
-
-    /// construct a multivector from blades
-    ///
-    constexpr explicit multivector(Bs... bs) : Bs{bs}... {}
-
-    /// construct a multivector from a scalar
-    ///
-    template <
-        class T = multivector<blade<>>,
-        class = std::enable_if_t<std::is_same_v<T, multivector>>>
-    constexpr explicit multivector(scalar_type s) : multivector(blade<>{s})
-    {}
-
-    /// access a specific blade
-    ///
-    /// @{
-    template <class B>
-    friend constexpr auto
-    get(multivector& x) -> std::enable_if_t<contains<B>, B&>
-    {
-      return static_cast<B&>(x);
-    }
-    template <class B>
-    friend constexpr auto
-    get(const multivector& x) -> std::enable_if_t<contains<B>, const B&>
-    {
-      return static_cast<const B&>(x);
-    }
-    template <class B>
-    friend constexpr auto
-    // NOLINTNEXTLINE(cppcoreguidelines-rvalue-reference-param-not-moved)
-    get(multivector&& x) -> std::enable_if_t<contains<B>, B&&>
-    {
-      return static_cast<B&&>(x);
-    }
-    template <class B>
-    friend constexpr auto
-    get(const multivector&& x) -> std::enable_if_t<contains<B>, const B&&>
-    {
-      return static_cast<const B&&>(x);
-    }
-    /// @}
-
-    /// access a specific blade or return a default value
-    ///
-    /// @{
-    template <class B>
-    friend constexpr auto
-    get_or(const multivector& x, const B& value) -> const B&
-    {
-      if constexpr (contains<B>) {
-        return get<B>(x);
-      } else {
-        return value;
-      }
-    }
-    /// @}
-
-    /// equality comparison
-    ///
-    template <class... B2s>
-    [[nodiscard]]
-    friend constexpr auto
-    operator==(const multivector& x, const multivector<B2s...>& y) -> bool
-    {
-      return ((get_or(x, Bs{}) == get_or(y, Bs{})) and ...) and
-             ((get_or(x, B2s{}) == get_or(y, B2s{})) and ...);
-    }
-
-    /// unary negation
-    ///
-    [[nodiscard]]
-    friend constexpr auto
-    operator-(const multivector& x) -> multivector
-    {
-      return multivector{-get<Bs>(x)...};
-    }
-
-    /// stream insertion
-    ///
-    friend auto operator<<(std::ostream& os, const multivector& x)
-        -> decltype(((os << std::declval<const Bs&>(), true) and ...), os)
-    {
-      auto first = true;
-
-      std::ignore =
-          ((os << (std::exchange(first, false) ? "" : " + "),
-            os << get<Bs>(x),
-            true) and
-           ...);
-
-      return os;
-    }
-  };
+  static constexpr auto i =
+      tmp::last_t<blade_list_type>{typename algebra_core_type::scalar_type{1}};
 };
 
 namespace detail {
diff --git a/geometry/src/algebra_core.hpp b/geometry/src/algebra_core.hpp
new file mode 100644
index 0000000..f481d02
--- /dev/null
+++ b/geometry/src/algebra_core.hpp
@@ -0,0 +1,414 @@
+#pragma once
+
+#include "geometry/src/detail/contract_dimensions.hpp"
+#include "geometry/src/detail/ordered.hpp"
+#include "geometry/src/detail/strictly_increasing.hpp"
+#include "geometry/src/get.hpp"
+#include "geometry/type_metaprogramming.hpp"
+
+#include <cstddef>
+#include <ostream>
+#include <tuple>
+#include <type_traits>
+#include <utility>
+
+namespace geometry {
+
+/// represents a projective geometric algebra
+/// @tparam S scalar type
+/// @tparam N number of euclidean geometry dimensions
+/// @tparam Dual true if the algebra is the dual (e.g. G* vs G)
+///
+/// Models an N-dimensional projective geometric algebra denoted as
+/// @f[
+/// P(S^*_{N,0,1})
+/// @f]
+///
+/// @note This type defines the `blade` and `multivector` types of an algebra
+///     does not define common constants.
+///
+template <class S, std::size_t N, bool Dual = true>
+struct algebra_core
+{
+#pragma GCC diagnostic push
+#ifndef __clang__
+#pragma GCC diagnostic ignored "-Wnon-template-friend"
+#endif  // __clang__
+  /// registration of algebra type for associated core type
+  /// see declaration in algebra
+  ///
+  constexpr friend auto algebra_for(algebra_core);
+#pragma GCC diagnostic pop
+
+  /// scalar type
+  ///
+  using scalar_type = S;
+
+  /// algebra dimension
+  ///
+  static constexpr auto dimension = N + 1;
+
+  /// dual algebra type
+  ///
+  using dual_algebra_core_type = algebra_core<S, N, not Dual>;
+
+  /// blade
+  /// @tparam Is dimensions
+  ///
+  /// Defines a blade for the given algebra.
+  ///
+  template <std::size_t... Is>
+  struct blade;
+
+  /// multivector
+  /// @tparam Bs blade types
+  ///
+  /// A linear combination of blades.
+  ///
+  /// @note a multivector may be composed of only blades with the same grade
+  ///
+  template <class... Bs>
+  struct multivector;
+
+  /// reified blade product
+  /// @tparam Is dimensions
+  ///
+  /// Returns the blade type after sorting `Is` and removing duplicates.
+  ///
+  /// @{
+  template <std::size_t... Is>
+  struct reified_blade
+  {
+    using type = tmp::convert_to_sequence_t<
+        detail::contract_dimensions_t<
+            detail::contraction_map::projective,
+            tmp::sort_t<tmp::index_constant_list<Is...>>>,
+        blade>;
+  };
+  template <std::size_t... Is>
+  using reified_blade_t = typename reified_blade<Is...>::type;
+  /// @}
+
+  /// coefficient for blade product
+  /// @tparam Is dimensions
+  ///
+  /// Returns the scalar coefficient for the reified blade type after sorting
+  /// `Is` and removing duplicates.
+  ///
+  /// @{
+  template <std::size_t... Is>
+  struct reified_blade_coefficient
+      : std::integral_constant<
+            int,
+            (tmp::sort_swap_count_v<tmp::index_constant_list<Is...>> % 2 == 0
+                 ? 1
+                 : -1) *
+                detail::contract_dimensions_coefficient_v<
+                    detail::contraction_map::projective,
+                    tmp::sort_t<tmp::index_constant_list<Is...>>>>
+  {};
+  template <std::size_t... Is>
+  static constexpr auto reified_blade_coefficient_v =
+      reified_blade_coefficient<Is...>::value;
+  /// @}
+
+  /// determines if a type is a blade
+  /// @tparam T type
+  ///
+  /// @{
+  template <class T>
+  struct is_blade : std::false_type
+  {};
+  template <std::size_t... Is>
+  struct is_blade<blade<Is...>> : std::true_type
+  {};
+  template <class T>
+  static constexpr auto is_blade_v = is_blade<T>::value;
+  /// @}
+
+  /// determines if a type is a multivector
+  /// @tparam T type
+  ///
+  /// @{
+  template <class T>
+  struct is_multivector : std::false_type
+  {};
+  template <class... Bs>
+  struct is_multivector<multivector<Bs...>> : std::true_type
+  {};
+  template <class T>
+  static constexpr auto is_multivector_v = is_multivector<T>::value;
+  /// @}
+
+  template <std::size_t... Is>
+  struct blade
+  {
+    static_assert(sizeof...(Is) <= N + 1, "number of `Is` exceeds rank");
+    static_assert(((Is <= N) and ...), "`Is` value exceeds rank");
+    static_assert(
+        detail::strictly_increasing(Is...),
+        "`Is` values must be unique and sorted");
+
+    /// algebra core type
+    ///
+    using algebra_core_type = algebra_core;
+
+    /// scalar type
+    ///
+    using scalar_type = typename algebra_core_type::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
+    {
+      return x.coefficient == y.coefficient;
+    }
+    [[nodiscard]]
+    friend constexpr auto
+    operator!=(blade x, blade y) noexcept -> bool
+    {
+      return x.coefficient != y.coefficient;
+    }
+    [[nodiscard]]
+    friend constexpr auto
+    operator<(blade x, blade y) noexcept -> bool
+    {
+      return x.coefficient < y.coefficient;
+    }
+    [[nodiscard]]
+    friend constexpr auto
+    operator>(blade x, blade y) noexcept -> bool
+    {
+      return x.coefficient > y.coefficient;
+    }
+    [[nodiscard]]
+    friend constexpr auto
+    operator<=(blade x, blade y) noexcept -> bool
+    {
+      return x.coefficient <= y.coefficient;
+    }
+    [[nodiscard]]
+    friend constexpr auto
+    operator>=(blade x, blade y) noexcept -> bool
+    {
+      return x.coefficient >= y.coefficient;
+    }
+    /// @}
+
+    /// unary negation
+    ///
+    [[nodiscard]]
+    friend constexpr auto
+    operator-(blade x) -> blade
+    {
+      return blade{-x.coefficient};
+    }
+
+    /// addition
+    ///
+    /// @{
+    [[nodiscard]]
+    friend constexpr auto
+    operator+(blade x, blade y) -> blade
+    {
+      return blade{x.coefficient + y.coefficient};
+    }
+    /// @}
+
+    /// compound assignment operators
+    ///
+    /// @{
+    friend constexpr auto& operator+=(blade& x, blade y) { return x = x + y; }
+    friend constexpr auto& operator-=(blade& x, blade y) { return x = x - y; }
+    /// @}
+
+    /// geometric product
+    ///
+    template <std::size_t... Js>
+    [[nodiscard]]
+    friend constexpr auto
+    operator*(blade x, blade<Js...> y)
+    {
+      constexpr auto unit_blade_coeff =
+          reified_blade_coefficient_v<Is..., Js...>;
+      using B = reified_blade_t<Is..., Js...>;
+
+      if constexpr (unit_blade_coeff == 0) {
+        return B{};
+      } else {
+        const auto coeff = x.coefficient * y.coefficient;
+        return B{unit_blade_coeff == 1 ? coeff : -coeff};
+      }
+    }
+
+    /// stream insertion
+    ///
+    friend auto operator<<(std::ostream& os, blade x) -> std::enable_if_t<
+        (N < 10),
+        decltype(os << std::declval<const scalar_type&>(), os)>
+    {
+      if (sizeof...(Is) == 0) {
+        os << x.coefficient;
+        return os;
+      }
+
+      if (x.coefficient != scalar_type{1}) {
+        os << x.coefficient << "*";
+      }
+
+      os << "e";
+      std::ignore = ((os << Is, true) and ...);
+
+      return os;
+    }
+  };
+
+  template <class... Bs>
+  struct multivector : Bs...
+  {
+    static_assert(
+        (is_blade_v<Bs> and ...), "`Bs` must be a specialization of blade");
+    static_assert(
+        tmp::same_v<typename Bs::scalar_type...>,
+        "`Bs` must have the same scalar type");
+    static_assert(
+        detail::strictly_increasing(detail::ordered<Bs>{}...),
+        "`Bs` must be lexicographically sorted");
+
+    /// algebra core type
+    ///
+    using algebra_core_type = algebra_core;
+
+    /// determine if a blade is contained in this multivector
+    ///
+    template <class B>
+    static constexpr auto contains = (std::is_same_v<B, Bs> or ...);
+
+    /// construct a zero multivector
+    ///
+    multivector() = default;
+
+    /// construct a multivector from blades
+    ///
+    constexpr explicit multivector(Bs... bs) : Bs{bs}... {}
+
+    /// construct a multivector from a scalar
+    ///
+    template <
+        class T = multivector<blade<>>,
+        class = std::enable_if_t<std::is_same_v<T, multivector>>>
+    constexpr explicit multivector(scalar_type s) : multivector(blade<>{s})
+    {}
+
+    /// access a specific blade
+    ///
+    /// @{
+    template <class B>
+    friend constexpr auto
+    get(multivector& x) -> std::enable_if_t<contains<B>, B&>
+    {
+      return static_cast<B&>(x);
+    }
+    template <class B>
+    friend constexpr auto
+    get(const multivector& x) -> std::enable_if_t<contains<B>, const B&>
+    {
+      return static_cast<const B&>(x);
+    }
+    template <class B>
+    friend constexpr auto
+    // NOLINTNEXTLINE(cppcoreguidelines-rvalue-reference-param-not-moved)
+    get(multivector&& x) -> std::enable_if_t<contains<B>, B&&>
+    {
+      return static_cast<B&&>(x);
+    }
+    template <class B>
+    friend constexpr auto
+    get(const multivector&& x) -> std::enable_if_t<contains<B>, const B&&>
+    {
+      return static_cast<const B&&>(x);
+    }
+    /// @}
+
+    /// access a specific blade or return a default value
+    ///
+    /// @{
+    template <class B>
+    friend constexpr auto
+    get_or(const multivector& x, const B& value) -> const B&
+    {
+      if constexpr (contains<B>) {
+        return get<B>(x);
+      } else {
+        return value;
+      }
+    }
+    /// @}
+
+    /// equality comparison
+    ///
+    template <class... B2s>
+    [[nodiscard]]
+    friend constexpr auto
+    operator==(const multivector& x, const multivector<B2s...>& y) -> bool
+    {
+      return ((get_or(x, Bs{}) == get_or(y, Bs{})) and ...) and
+             ((get_or(x, B2s{}) == get_or(y, B2s{})) and ...);
+    }
+
+    /// unary negation
+    ///
+    [[nodiscard]]
+    friend constexpr auto
+    operator-(const multivector& x) -> multivector
+    {
+      return multivector{-get<Bs>(x)...};
+    }
+
+    /// stream insertion
+    ///
+    friend auto operator<<(std::ostream& os, const multivector& x)
+        -> decltype(((os << std::declval<const Bs&>(), true) and ...), os)
+    {
+      auto first = true;
+
+      std::ignore =
+          ((os << (std::exchange(first, false) ? "" : " + "),
+            os << get<Bs>(x),
+            true) and
+           ...);
+
+      return os;
+    }
+  };
+};
+
+}  // namespace geometry
diff --git a/geometry/src/algebra_type.hpp b/geometry/src/algebra_type.hpp
index d45cd8d..db34e4f 100644
--- a/geometry/src/algebra_type.hpp
+++ b/geometry/src/algebra_type.hpp
@@ -5,14 +5,39 @@
 namespace geometry {
 namespace detail {
 
+auto algebra_for() -> void = delete;
+
+inline constexpr class
+{
+  template <class T>
+  static constexpr auto impl(int) -> typename T::algebra_type
+  {
+    return {};
+  }
+  template <class T>
+  static constexpr auto
+  impl(float) -> decltype(algebra_for(typename T::algebra_core_type{}))
+  {
+    return {};
+  }
+
+public:
+  template <class T>
+  constexpr auto operator()(const T&) const -> decltype(impl<T>(int{}))
+  {
+    return impl<T>(int{});
+  }
+
+} algebra_for_fn{};
+
 template <class T, class = void>
 struct algebra_type
 {};
 
 template <class T>
-struct algebra_type<T, std::void_t<typename T::algebra_type>>
+struct algebra_type<T, std::void_t<decltype(algebra_for_fn(std::declval<T>()))>>
 {
-  using type = typename T::algebra_type;
+  using type = decltype(algebra_for_fn(std::declval<T>()));
 };
 
 }  // namespace detail
diff --git a/geometry/src/detail/blade_list.hpp b/geometry/src/detail/blade_list.hpp
index 445ed6e..0d8ea01 100644
--- a/geometry/src/detail/blade_list.hpp
+++ b/geometry/src/detail/blade_list.hpp
@@ -43,7 +43,7 @@ struct grade
 template <template <std::size_t...> class blade>
 struct blade_list
 {
-  static constexpr auto N = blade<>::algebra_type::dimension;
+  static constexpr auto N = blade<>::algebra_core_type::dimension;
 
   using type = tmp::flatten_t<tmp::transform_t<
       tmp::convert_to_list_t<std::make_index_sequence<N + 1>, tmp::list>,
diff --git a/geometry/src/exterior_product.hpp b/geometry/src/exterior_product.hpp
index fbba33e..f0def45 100644
--- a/geometry/src/exterior_product.hpp
+++ b/geometry/src/exterior_product.hpp
@@ -20,12 +20,12 @@ inline constexpr class
       class D1 = std::decay_t<typename T1::eval_type>,
       class D2 = std::decay_t<typename T2::eval_type>,
       class Algebra = common_algebra_type_t<D1, D2>,
-      class L = detail::pending_dimensions_list_t<D1, D2>>
+      class L = detail::pending_dimensions_list_t<D1, D2>,
+      class R =
+          typename tmp::convert_to_sequence_t<L, Algebra::template reified_blade>::type>
   constexpr auto
   operator()(expression_template::op<expression_template::multiplies, T1, T2>)
-      const -> std::bool_constant<(
-          tmp::convert_to_sequence_t<L, Algebra::template reified_blade>::type::
-              grade < D1::grade + D2::grade)>
+      const -> std::bool_constant<R::grade != D1::grade + D2::grade>
   {
     return {};
   }
diff --git a/geometry/src/interior_product.hpp b/geometry/src/interior_product.hpp
new file mode 100644
index 0000000..b63fcab
--- /dev/null
+++ b/geometry/src/interior_product.hpp
@@ -0,0 +1,63 @@
+#pragma once
+
+#include "geometry/expression_template.hpp"
+#include "geometry/src/common_algebra_type.hpp"
+#include "geometry/src/detail/multivector_product.hpp"
+#include "geometry/src/detail/pending_dimensions_list.hpp"
+#include "geometry/type_metaprogramming.hpp"
+
+#include <type_traits>
+
+namespace geometry {
+namespace detail {
+
+inline constexpr class
+{
+public:
+  template <
+      class T1,
+      class T2,
+      class D1 = std::decay_t<typename T1::eval_type>,
+      class D2 = std::decay_t<typename T2::eval_type>,
+      class Algebra = common_algebra_type_t<D1, D2>,
+      class L = detail::pending_dimensions_list_t<D1, D2>,
+      class R =
+          typename tmp::convert_to_sequence_t<L, Algebra::template reified_blade>::type,
+      int k = tmp::convert_to_sequence_t<
+          L,
+          Algebra::template reified_blade_coefficient>::value>
+  constexpr auto
+  operator()(expression_template::op<expression_template::multiplies, T1, T2>)
+      const -> std::bool_constant<
+          (R::grade !=
+           std::max(D1::grade, D2::grade) - std::min(D1::grade, D2::grade)) or
+          (k == 0)>
+  {
+    return {};
+  }
+} non_grade_contracting{};
+
+}  // namespace detail
+
+/// interior product
+///
+/// @{
+
+inline constexpr auto interior_product =  //
+    detail::multivector_product_with(     //
+        detail::non_grade_contracting,
+        [](auto algebra) {
+          return typename decltype(algebra)::template blade<>{};
+        });
+
+template <class T1, class T2, class A = common_algebra_type_t<T1, T2>>
+[[nodiscard]]
+constexpr auto
+operator|(const T1& x, const T2& y)
+{
+  return interior_product(x, y);
+}
+
+/// @}
+
+}  // namespace geometry
diff --git a/test/BUILD.bazel b/test/BUILD.bazel
index 1da6b41..7cbee02 100644
--- a/test/BUILD.bazel
+++ b/test/BUILD.bazel
@@ -79,3 +79,14 @@ cc_test(
         "@skytest",
     ],
 )
+
+cc_test(
+    name = "algebra_interior_product_test",
+    size = "small",
+    srcs = ["algebra_interior_product_test.cpp"],
+    deps = [
+        ":test_same",
+        "//:geometry",
+        "@skytest",
+    ],
+)
diff --git a/test/algebra_interior_product_test.cpp b/test/algebra_interior_product_test.cpp
new file mode 100644
index 0000000..d037097
--- /dev/null
+++ b/test/algebra_interior_product_test.cpp
@@ -0,0 +1,44 @@
+#include "geometry/geometry.hpp"
+#include "skytest/skytest.hpp"
+
+#include "test/test_same.hpp"
+
+#include <cstddef>
+#include <tuple>
+#include <utility>
+
+auto main() -> int
+{
+  using namespace ::skytest::literals;
+  using ::geometry::get;
+  using ::geometry::test::same;
+  using ::skytest::eq;
+  using ::skytest::expect;
+  using ::skytest::ne;
+
+  "line orthogonal to point (2D)"_ctest = [] {
+    using algebra = ::geometry::algebra<double, 2>;
+
+    // https://bivector.net/PROJECTIVE_GEOMETRIC_ALGEBRA.pdf
+    // section 7.1
+
+    constexpr auto x_hat = algebra::e<2, 0>;
+    constexpr auto y_hat = algebra::e<0, 1>;
+    constexpr auto E0 = algebra::e<1, 2>;
+
+    const auto P = 0 * x_hat + 0 * y_hat + E0;
+
+    // a line of the form
+    // a * x + b * y + c = 0
+    // is represented as
+    // a * e1 + b * e2 + c * e0 = 0
+
+    //using e0 = algebra::blade<0>;
+    constexpr auto x = algebra::e<1>;
+    constexpr auto y = algebra::e<2>;
+
+    return expect(eq(x, y | P) and
+                  eq(y, x | P));
+  };
+
+}