diff --git a/src/core/include/units/magnitude.h b/src/core/include/units/magnitude.h index a2d33c08..f757f020 100644 --- a/src/core/include/units/magnitude.h +++ b/src/core/include/units/magnitude.h @@ -26,37 +26,49 @@ #include #include -namespace units::mag -{ +namespace units::mag { -namespace detail -{ -// Helpers to perform prime factorization at compile time. -template - requires requires { N > 0; } -struct prime_factorization; -template -static constexpr auto prime_factorization_v = prime_factorization::value; - -// A way to check whether a number is prime at compile time. -constexpr bool is_prime(std::intmax_t n); -} // namespace detail - -// Integer rep is for prime numbers; long double is for any irrational base we permit. +/** + * @brief Any type which can be used as a basis vector in a BasePower. + * + * We have two categories. + * + * The first is just an `int`. This is for prime number bases. These can always be used directly as NTTPs. + * + * The second category is a _custom type_, which has a static member variable of type `long double` that holds its + * value. We choose `long double` to get the greatest degree of precision; users who need a different type can convert + * from this at compile time. This category is for any irrational base we admit into our representation (on which, more + * details below). + * + * The reason we can't hold the value directly for floating point bases is so that we can support some compilers (e.g., + * GCC 10) which don't yet permit floating point NTTPs. + */ template concept BaseRep = std::is_same_v || std::is_same_v, long double>; /** * @brief A basis vector in our magnitude representation, raised to some rational power. * - * The set of basis vectors must be linearly independent: that is, no product of basis powers can ever equal 1, unless - * all exponents are zero. To achieve this, we use the following kinds of basis vectors. + * The public API is that there is a `power` member variable (of type `ratio`), and a `get_base()` member function (of + * type either `int` or `long double`, as appropriate), for any specialization. + * + * These types exist to be used as NTTPs for the variadic `magnitude<...>` template. We represent a magnitude (which is + * a positive real number) as the product of rational powers of "basis vectors", where each "basis vector" is a positive + * real number. "Addition" in this vector space corresponds to _multiplying_ two real numbers. "Scalar multiplication" + * corresponds to _raising_ a real number to a _rational power_. Thus, this representation of positive real numbers + * maps them onto a vector space over the rationals, supporting the operations of products and rational powers. + * + * (Note that this is the same representation we already use for dimensions.) + * + * As in any vector space, the set of basis vectors must be linearly independent: that is, no product of basis powers + * can ever give the null vector (which in this case represents the number 1), unless all scalars (again, in this case, + * _powers_) are zero. To achieve this, we use the following kinds of basis vectors. * - Prime numbers. (These are the only allowable values for `int` base.) * - Certain selected irrational numbers, such as pi. * * Before allowing a new irrational base, make sure that it _cannot_ be represented as the product of rational powers of - * _existing_ bases, including both prime numbers and any other irrational bases. For example, even though sqrt(2) is - * irrational, we must not ever use it as a base + * _existing_ bases, including both prime numbers and any other irrational bases. For example, even though `sqrt(2)` is + * irrational, we must not ever use it as a base; instead, we would use `base_power{2, ratio{1, 2}}`. */ template struct base_power { @@ -66,9 +78,12 @@ struct base_power { constexpr long double get_base() const { return T::value; } }; +/** + * @brief Specialization for prime number bases. + */ template<> struct base_power { - // The value of the basis "vector". + // The value of the basis "vector". Must be prime to be used with `magnitude` (below). int base; // The rational power to which the base is raised. @@ -77,103 +92,88 @@ struct base_power { constexpr int get_base() const { return base; } }; -template U> -base_power(T, U) -> base_power; - -template -base_power(T) -> base_power; - -template -constexpr bool operator==(base_power t, base_power u) { - return std::is_same_v && (t.get_base() == u.get_base()) && (t.power == u.power); -} - -template -constexpr auto inverse(base_power bp) { - bp.power = -bp.power; - return bp; -} - -namespace detail -{ -template -constexpr bool is_valid_base_power(const base_power &bp) { - if (bp.power == 0) { return false; } - - if constexpr (std::is_same_v) { return is_prime(bp.get_base()); } - else if constexpr (std::is_same_v) { return bp.get_base() > 0; } - else { return false; } // Unreachable. -} +/** + * @brief Deduction guides for base_power: only permit deducing integral bases. + */ +template U> +base_power(T, U) -> base_power; +template +base_power(T) -> base_power; +// Implementation for BasePower concept (below). +namespace detail { template struct is_base_power : std::false_type {}; template struct is_base_power> : std::true_type {}; } // namespace detail +/** + * @brief Concept to detect whether a _type_ is a valid base power. + * + * Note that this is somewhat incomplete. We must also detect whether a _value_ of that type is valid for use with + * `magnitude<...>`. We will defer that second check to the constraints on the `magnitude` template. + */ template concept BasePower = detail::is_base_power::value; +/** + * @brief Equality detection for two base powers. + */ +template +constexpr bool operator==(T t, U u) { + return std::is_same_v && (t.get_base() == u.get_base()) && (t.power == u.power); +} + +/** + * @brief The (multiplicative) inverse of a BasePower. + */ +template +constexpr auto inverse(BP bp) { + bp.power = -bp.power; + return bp; +} + +// Implementation helpers for `magnitude<...>` constraints (below). +namespace detail { +constexpr bool is_valid_base_power(const BasePower auto &bp); + template constexpr bool strictly_increasing(Ts&&... ts); +} // namespace detail /** * @brief A representation for positive real numbers which optimizes taking products and rational powers. * * Magnitudes can be treated as values. Each type encodes exactly one value. Users can multiply, divide, and compare - * for equality using this value API. + * for equality. */ template requires requires { - (detail::is_valid_base_power(BPs) && ... && strictly_increasing(BPs.get_base()...)); + (detail::is_valid_base_power(BPs) && ... && detail::strictly_increasing(BPs.get_base()...)); } struct magnitude {}; -template -constexpr bool operator==(magnitude, magnitude) { - if constexpr (sizeof...(LeftBPs) == sizeof...(RightBPs)) { return ((LeftBPs == RightBPs) && ...); } - else { return false; } -} - +// Implementation for Magnitude concept (below). +namespace detail { +template +struct is_magnitude : std::false_type {}; template -constexpr auto inverse(magnitude) { return magnitude{}; } +struct is_magnitude> : std::true_type {}; +} // namespace detail -constexpr auto operator*(magnitude<>, magnitude<>) { return magnitude<>{}; } +/** + * @brief Concept to detect whether T is a valid Magnitude. + */ +template +concept Magnitude = detail::is_magnitude::value; -template -constexpr auto operator*(magnitude<>, magnitude m) { return m; } - -template -constexpr auto operator*(magnitude m, magnitude<>) { return m; } - -template -constexpr auto operator*(magnitude, magnitude) { - // Shortcut for prepending, which makes it easier to implement some of the other cases. - if constexpr ((sizeof...(T1) == 0) && H1.get_base() < H2.get_base()) { return magnitude{}; } - - if constexpr (H1.get_base() == H2.get_base()) { - constexpr auto partial_product = magnitude{} * magnitude{}; - - // Make a new base_power with the common base of H1 and H2, whose power is their powers' sum. - constexpr auto new_head = [&](auto head) { - head.power = H1.power + H2.power; - return head; - }(H1); - - if constexpr (new_head.power == 0) { - return partial_product; - } else { - return magnitude{} * partial_product; - } - } else if constexpr(H1.get_base() < H2.get_base()){ - return magnitude

{} * (magnitude{} * magnitude{}); - } else { // We know H2.get_base() < H1.get_base() - return magnitude

{} * (magnitude{} * magnitude{}); - } -} - -template -constexpr auto operator/(magnitude l, magnitude r) { return l * inverse(r); } +/** + * @brief Convert any positive integer to a Magnitude. + */ +template + requires requires { N > 0; } +constexpr Magnitude auto as_magnitude(); /** * @brief Make a Magnitude that is a rational number. @@ -182,7 +182,7 @@ constexpr auto operator/(magnitude l, magnitude r) { re * manually adding base powers. */ template -constexpr auto make_ratio() { return detail::prime_factorization_v / detail::prime_factorization_v; } +constexpr auto make_ratio() { return as_magnitude() / as_magnitude(); } /** * @brief A base to represent pi. @@ -191,6 +191,9 @@ struct pi_base { static constexpr long double value = std::numbers::pi_v; }; +/** + * @brief A simple way to create a Magnitude representing a rational power of pi. + */ template constexpr auto pi_to_the() { return magnitude{Power}>{}; } @@ -198,42 +201,9 @@ constexpr auto pi_to_the() { return magnitude{Power}>{}; } // Implementation details below. //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// -//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// -// Magnitude concept implementation. - -template -struct pairwise_all { - Predicate predicate; - - template - constexpr bool operator()(Ts&&... ts) const { - // Carefully handle different sizes, avoiding unsigned integer underflow. - constexpr auto num_comparisons = [](auto num_elements) { - return (num_elements > 1) ? (num_elements - 1) : 0; - }(sizeof...(Ts)); - - // Compare zero or more pairs of neighbours as needed. - return [this](std::tuple &&t, std::index_sequence) { - return (predicate(std::get(t), std::get(t)) && ...); - }(std::make_tuple(std::forward(ts)...), std::make_index_sequence()); - } -}; - -template -pairwise_all(T) -> pairwise_all; - -// Check whether a tuple of (possibly heterogeneously typed) values are strictly increasing. -template -constexpr bool strictly_increasing(Ts&&... ts) { - return pairwise_all{std::less{}}(std::forward(ts)...); -} - namespace detail { -//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// -// Prime factorization implementation. - // Find the smallest prime factor of `n`. constexpr std::intmax_t find_first_factor(std::intmax_t n) { @@ -265,11 +235,7 @@ constexpr std::intmax_t remove_power(std::intmax_t base, std::intmax_t pow, std: return n; } -// Specialization for the prime factorization of 1 (base case). -template<> -struct prime_factorization<1> { static constexpr magnitude<> value{}; }; - -// Specialization for the prime factorization of larger numbers (recursive case). +// Helpers to perform prime factorization at compile time. template requires requires { N > 0; } struct prime_factorization { @@ -278,10 +244,122 @@ struct prime_factorization { static constexpr std::intmax_t remainder = remove_power(first_base, first_power, N); static constexpr auto value = magnitude{} - * prime_factorization_v; + * prime_factorization::value; }; +template +static constexpr auto prime_factorization_v = prime_factorization::value; + +// Specialization for the prime factorization of 1 (base case). +template<> +struct prime_factorization<1> { static constexpr magnitude<> value{}; }; + +// A way to check whether a number is prime at compile time. constexpr bool is_prime(std::intmax_t n) { return (n > 1) && (find_first_factor(n) == n); } +constexpr bool is_valid_base_power(const BasePower auto &bp) { + if (bp.power == 0) { return false; } + + if constexpr (std::is_same_v) { return is_prime(bp.get_base()); } + else { return bp.get_base() > 0; } +} + +// A function object to apply a predicate to all consecutive pairs of values in a sequence. +template +struct pairwise_all { + Predicate predicate; + + template + constexpr bool operator()(Ts&&... ts) const { + // Carefully handle different sizes, avoiding unsigned integer underflow. + constexpr auto num_comparisons = [](auto num_elements) { + return (num_elements > 1) ? (num_elements - 1) : 0; + }(sizeof...(Ts)); + + // Compare zero or more pairs of neighbours as needed. + return [this](std::tuple &&t, std::index_sequence) { + return (predicate(std::get(t), std::get(t)) && ...); + }(std::make_tuple(std::forward(ts)...), std::make_index_sequence()); + } +}; + +// Deduction guide: permit constructions such as `pairwise_all{std::less{}}`. +template +pairwise_all(T) -> pairwise_all; + +// Check whether a sequence of (possibly heterogeneously typed) values are strictly increasing. +template +constexpr bool strictly_increasing(Ts&&... ts) { + return pairwise_all{std::less{}}(std::forward(ts)...); +} + } // namespace detail + +//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// +// Magnitude equality implementation. + +template +constexpr bool operator==(magnitude, magnitude) { + if constexpr (sizeof...(LeftBPs) == sizeof...(RightBPs)) { return ((LeftBPs == RightBPs) && ...); } + else { return false; } +} + +//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// +// Magnitude inverse implementation. + +template +constexpr auto inverse(magnitude) { return magnitude{}; } + +//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// +// Magnitude product implementation. + +// Base cases, for when either (or both) inputs are the identity. +constexpr auto operator*(magnitude<>, magnitude<>) { return magnitude<>{}; } +constexpr auto operator*(magnitude<>, Magnitude auto m) { return m; } +constexpr auto operator*(Magnitude auto m, magnitude<>) { return m; } + +// Recursive case for the product of any two non-identity Magnitudes. +template +constexpr auto operator*(magnitude, magnitude) { + // Shortcut for the "pure prepend" case, which makes it easier to implement some of the other cases. + if constexpr ((sizeof...(T1) == 0) && H1.get_base() < H2.get_base()) { return magnitude{}; } + + // "Same leading base" case. + if constexpr (H1.get_base() == H2.get_base()) { + constexpr auto partial_product = magnitude{} * magnitude{}; + + // Make a new base_power with the common base of H1 and H2, whose power is their powers' sum. + constexpr auto new_head = [&](auto head) { + head.power = H1.power + H2.power; + return head; + }(H1); + + if constexpr (new_head.power == 0) { + return partial_product; + } else { + return magnitude{} * partial_product; + } + } + // Case for when H1 has the smaller base. + else if constexpr(H1.get_base() < H2.get_base()){ + return magnitude

{} * (magnitude{} * magnitude{}); + } + // Case for when H2 has the smaller base. + else { + return magnitude

{} * (magnitude{} * magnitude{}); + } +} + +//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// +// Magnitude quotient implementation. + +constexpr auto operator/(Magnitude auto l, Magnitude auto r) { return l * inverse(r); } + +//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// +// `as_magnitude()` implementation. + +template + requires requires { N > 0; } +constexpr Magnitude auto as_magnitude() { return detail::prime_factorization_v; } + } // namespace units::mag diff --git a/test/unit_test/runtime/magnitude_test.cpp b/test/unit_test/runtime/magnitude_test.cpp index a4a5e74e..cf6f5ecd 100644 --- a/test/unit_test/runtime/magnitude_test.cpp +++ b/test/unit_test/runtime/magnitude_test.cpp @@ -28,28 +28,6 @@ namespace units::mag { -TEST_CASE("strictly_increasing") -{ - SECTION ("Empty input is sorted") - { - CHECK(strictly_increasing()); - } - - SECTION ("Single-element input is sorted") - { - CHECK(strictly_increasing(3)); - CHECK(strictly_increasing(15.42)); - CHECK(strictly_increasing('c')); - } - - SECTION ("Multi-value inputs compare correctly") - { - CHECK(strictly_increasing(3, 3.14)); - CHECK(!strictly_increasing(3, 3.0)); - CHECK(!strictly_increasing(4, 3.0)); - } -} - TEST_CASE("make_ratio performs prime factorization correctly") { SECTION("Performs prime factorization when denominator is 1") @@ -218,6 +196,28 @@ TEST_CASE("pairwise_all evaluates all pairs") } } +TEST_CASE("strictly_increasing") +{ + SECTION ("Empty input is sorted") + { + CHECK(strictly_increasing()); + } + + SECTION ("Single-element input is sorted") + { + CHECK(strictly_increasing(3)); + CHECK(strictly_increasing(15.42)); + CHECK(strictly_increasing('c')); + } + + SECTION ("Multi-value inputs compare correctly") + { + CHECK(strictly_increasing(3, 3.14)); + CHECK(!strictly_increasing(3, 3.0)); + CHECK(!strictly_increasing(4, 3.0)); + } +} + } // namespace detail } // namespace units::mag