| // Copyright 2019 The Fuchsia Authors. All rights reserved. |
| // Use of this source code is governed by a BSD-style license that can be |
| // found in the LICENSE file. |
| |
| #include <lib/affine/ratio.h> |
| #include <type_traits> |
| #include <zircon/assert.h> |
| |
| namespace affine { |
| namespace { |
| |
| // Calculates the greatest common denominator (factor) of two values. |
| template <typename T> |
| T BinaryGcd(T a, T b) { |
| internal::DebugAssert(a && b); |
| |
| // Remove and count the common factors of 2. |
| uint8_t twos; |
| for (twos = 0; ((a | b) & 1) == 0; ++twos) { |
| a >>= 1; |
| b >>= 1; |
| } |
| |
| // Get rid of the non-common factors of 2 in a. a is non-zero, so this |
| // terminates. |
| while ((a & 1) == 0) { |
| a >>= 1; |
| } |
| |
| do { |
| // Get rid of the non-common factors of 2 in b. b is non-zero, so this |
| // terminates. |
| while ((b & 1) == 0) { |
| b >>= 1; |
| } |
| |
| // Apply the Euclid subtraction method. |
| if (a > b) { |
| std::swap(a, b); |
| } |
| |
| b = b - a; |
| } while (b != 0); |
| |
| // Multiply in the common factors of two. |
| return a << twos; |
| } |
| |
| enum class RoundDirection { Down, Up }; |
| |
| // Scales a uint64_t value by the ratio of two uint32_t values. If round_up is |
| // true, the result is rounded up rather than down. overflow is set to indicate |
| // overflow. |
| template <RoundDirection ROUND_DIR, uint64_t OVERFLOW_LIMIT_64> |
| uint64_t ScaleUInt64(uint64_t value, uint32_t numerator, uint32_t denominator) { |
| constexpr uint64_t kLow32Bits = 0xffffffffu; |
| |
| // high and low are the product of the numerator and the high and low halves |
| // (respectively) of value. |
| uint64_t high = numerator * (value >> 32u); |
| uint64_t low = numerator * (value & kLow32Bits); |
| |
| // Ignoring overflow and remainder, the result we want is: |
| // ((high << 32) + low) / denominator. |
| |
| // Move the high end of low into the low end of high. |
| high += low >> 32u; |
| low = low & kLow32Bits; |
| |
| // Ignoring overflow and remainder, the result we want is still: |
| // ((high << 32) + low) / denominator. |
| |
| // Compute the divmod of high/D |
| uint64_t high_q = high / denominator; |
| uint64_t high_r = high % denominator; |
| |
| // If high_q is larger than the overflow limit, then we can just get out now. |
| // The overflow limit will be different depending on whether we are scaling |
| // a non-negative number (0x7FFFFFFF) or a negative number (0x80000000) |
| constexpr uint64_t OVERFLOW_LIMIT_32 = OVERFLOW_LIMIT_64 >> 32; |
| if (high_q > OVERFLOW_LIMIT_32) { |
| return OVERFLOW_LIMIT_64; |
| } |
| |
| // The remainder of high/D are the high bits of low. Or them in, and do the |
| // divmod for the low portion |
| low |= high_r << 32u; |
| |
| uint64_t low_q = low / denominator; |
| __UNUSED uint64_t low_r = low % denominator; |
| |
| uint64_t result = (high_q << 32u) | low_q; |
| |
| if constexpr (ROUND_DIR == RoundDirection::Up) { |
| if (result >= OVERFLOW_LIMIT_64) { |
| return OVERFLOW_LIMIT_64; |
| } |
| if (low_r) { |
| ++result; |
| } |
| } |
| |
| return result; |
| } |
| |
| constexpr bool FitsIn32Bits(uint64_t numerator, uint64_t denominator) { |
| return ((numerator <= std::numeric_limits<uint32_t>::max()) && |
| (denominator <= std::numeric_limits<uint32_t>::max())); |
| } |
| |
| } // namespace |
| |
| template <typename T> |
| void Ratio::Reduce(T* numerator, T* denominator) { |
| internal::DebugAssert(numerator != nullptr); |
| internal::DebugAssert(denominator != nullptr); |
| internal::Assert(*denominator != 0); |
| |
| if (*numerator == 0) { |
| *denominator = 1; |
| return; |
| } |
| |
| T gcd = BinaryGcd(*numerator, *denominator); |
| internal::DebugAssert(gcd != 0); |
| |
| if (gcd == 1) { |
| return; |
| } |
| |
| *numerator = *numerator / gcd; |
| *denominator = *denominator / gcd; |
| } |
| |
| void Ratio::Product(uint32_t a_numerator, uint32_t a_denominator, uint32_t b_numerator, |
| uint32_t b_denominator, uint32_t* product_numerator, |
| uint32_t* product_denominator, Exact exact) { |
| internal::DebugAssert(product_numerator != nullptr); |
| internal::DebugAssert(product_denominator != nullptr); |
| |
| uint64_t numerator = static_cast<uint64_t>(a_numerator) * b_numerator; |
| uint64_t denominator = static_cast<uint64_t>(a_denominator) * b_denominator; |
| |
| Ratio::Reduce(&numerator, &denominator); |
| |
| if (!FitsIn32Bits(numerator, denominator)) { |
| internal::Assert(exact == Exact::No); |
| |
| // Try to find the best approximation of the ratio that we can. Our |
| // approach is as follows. Figure out the number of bits to the right |
| // we need to shift the numerator and denominator, rounding up or down |
| // in the process, such that the result can be reduced to fit into 32 |
| // bits. |
| // |
| // This approach tends to beat out a just-shift-until-it-fits approach, |
| // as well as an always-shift-then-reduce approach, but _none_ of these |
| // approaches always finds the best solution. |
| // |
| // TODO(johngro) : figure out if it is reasonable to actually compute |
| // the best solution. Alternatively, consider implementing a "just |
| // shift until it fits" solution if the approximate results are good |
| // enough. |
| // |
| for (uint32_t i = 1; i <= 32; ++i) { |
| // Produce a version of the numerator and denominator which have |
| // each been divided by 2^i, rounding up/down as appropriate |
| // (instead of truncating). |
| uint64_t rounded_numerator = (numerator + (static_cast<uint64_t>(1) << (i - 1))) >> i; |
| uint64_t rounded_denominator = (denominator + (static_cast<uint64_t>(1) << (i - 1))) >> i; |
| |
| if (rounded_denominator == 0) { |
| // Product is larger than we can represent. Return the largest value we |
| // can represent. |
| *product_numerator = std::numeric_limits<uint32_t>::max(); |
| *product_denominator = 1; |
| return; |
| } |
| |
| if (rounded_numerator == 0) { |
| // Product is smaller than we can represent. Return 0. |
| *product_numerator = 0; |
| *product_denominator = 1; |
| return; |
| } |
| |
| Ratio::Reduce(&rounded_numerator, &rounded_denominator); |
| if (FitsIn32Bits(rounded_numerator, rounded_denominator)) { |
| *product_numerator = static_cast<uint32_t>(rounded_numerator); |
| *product_denominator = static_cast<uint32_t>(rounded_denominator); |
| return; |
| } |
| } |
| } |
| |
| *product_numerator = static_cast<uint32_t>(numerator); |
| *product_denominator = static_cast<uint32_t>(denominator); |
| } |
| |
| int64_t Ratio::Scale(int64_t value, uint32_t numerator, uint32_t denominator) { |
| internal::Assert(denominator != 0u); |
| |
| if (value >= 0) { |
| // LIMIT == 0x7FFFFFFFFFFFFFFF |
| constexpr uint64_t LIMIT = static_cast<uint64_t>(std::numeric_limits<int64_t>::max()); |
| return static_cast<int64_t>(ScaleUInt64<RoundDirection::Down, LIMIT>( |
| static_cast<uint64_t>(value), numerator, denominator)); |
| } else { |
| // LIMIT == 0x8000000000000000 |
| // |
| // Note: We are attempting to pass the unsigned distance from zero into |
| // our ScaleUInt64 function. In the case of negative numbers, we pass |
| // the twos compliment into the scale function, and then flip the sign |
| // again on the way out. |
| // |
| // We are taking the advantage of the fact that the twos compliment of |
| // MIN is itself for any signed integer type, and that casting this |
| // value to an unsigned integer of the same size properly produces the |
| // original value's distance from zero. Clamping the limit to the |
| // distance of MIN from zero means that saturated results will likewise |
| // get properly flipped back to MIN during the return. |
| // |
| constexpr uint64_t LIMIT = static_cast<uint64_t>(std::numeric_limits<int64_t>::min()); |
| return static_cast<int64_t>(-ScaleUInt64<RoundDirection::Up, LIMIT>( |
| -static_cast<uint64_t>(value), numerator, denominator)); |
| } |
| } |
| |
| // Explicit instantiation of the two supported types of Ratio |
| template void Ratio::Reduce<uint32_t>(uint32_t*, uint32_t*); |
| template void Ratio::Reduce<uint64_t>(uint64_t*, uint64_t*); |
| |
| } // namespace affine |