blob: 25c7b7c27a633dc1b5a445f3db14ce7e0d926cff [file] [log] [blame] [edit]
// 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