stern_brocot.hpp¶
Run-length encoded L/R path from 1/1 to a positive reduced fraction.
Verified by stern_brocot_tree.
\[
\displaystyle \frac{a}{b}\oplus\frac{c}{d}=\frac{a+c}{b+d}
\]
Implementation¶
#ifndef NOYA_STERN_BROCOT_HPP
#define NOYA_STERN_BROCOT_HPP 1
/// @complexity Time: O(log max(numerator,denominator)) run steps.
/// Space: O(number of returned runs).
#include <cassert>
#include <cstdint>
#include <limits>
#include <numeric>
#include <optional>
#include <utility>
#include <vector>
namespace noya {
using stern_brocot_step = std::pair<char, std::uint64_t>;
/// @brief Run-length encoded L/R path from 1/1 to a positive reduced fraction.
inline std::vector<stern_brocot_step>
stern_brocot_path(std::uint64_t numerator, std::uint64_t denominator) {
assert(numerator >= 1 && denominator >= 1);
assert(std::gcd(numerator, denominator) == 1);
std::vector<stern_brocot_step> path;
while (numerator != denominator) {
if (numerator < denominator) {
std::uint64_t count = (denominator - 1) / numerator;
path.emplace_back('L', count);
denominator -= count * numerator;
} else {
std::uint64_t count = (numerator - 1) / denominator;
path.emplace_back('R', count);
numerator -= count * denominator;
}
}
return path;
}
/// @brief Decode a run-length encoded Stern-Brocot path to its reduced
/// fraction; asserts when a 64-bit numerator or denominator would overflow.
inline std::pair<std::uint64_t, std::uint64_t>
stern_brocot_fraction(const std::vector<stern_brocot_step> &path) {
using u128 = unsigned __int128;
u128 left_numerator = 0;
u128 left_denominator = 1;
u128 numerator = 1;
u128 denominator = 1;
u128 right_numerator = 1;
u128 right_denominator = 0;
constexpr u128 maximum = std::numeric_limits<std::uint64_t>::max();
for (auto [direction, count] : path) {
assert((direction == 'L' || direction == 'R') && count >= 1);
if (direction == 'L') {
right_numerator = numerator + u128(count - 1) * left_numerator;
right_denominator = denominator + u128(count - 1) * left_denominator;
numerator += u128(count) * left_numerator;
denominator += u128(count) * left_denominator;
} else {
left_numerator = numerator + u128(count - 1) * right_numerator;
left_denominator = denominator + u128(count - 1) * right_denominator;
numerator += u128(count) * right_numerator;
denominator += u128(count) * right_denominator;
}
assert(numerator <= maximum && denominator <= maximum);
}
return {std::uint64_t(numerator), std::uint64_t(denominator)};
}
/// @brief Return the Stern-Brocot node at the given depth on the root-to-node
/// path, or nullopt when the requested depth is greater than the node depth.
inline std::optional<std::pair<std::uint64_t, std::uint64_t>>
stern_brocot_ancestor(std::uint64_t depth, std::uint64_t numerator,
std::uint64_t denominator) {
auto path = stern_brocot_path(numerator, denominator);
std::vector<stern_brocot_step> prefix;
for (auto [direction, count] : path) {
std::uint64_t used = std::min(depth, count);
if (used > 0) {
prefix.emplace_back(direction, used);
}
if (depth <= count) {
return stern_brocot_fraction(prefix);
}
depth -= count;
}
if (depth == 0) {
return stern_brocot_fraction(prefix);
}
return std::nullopt;
}
/// @brief Return the lowest common ancestor of two positive reduced fractions.
inline std::pair<std::uint64_t, std::uint64_t>
stern_brocot_lca(std::uint64_t first_numerator, std::uint64_t first_denominator,
std::uint64_t second_numerator,
std::uint64_t second_denominator) {
auto first = stern_brocot_path(first_numerator, first_denominator);
auto second = stern_brocot_path(second_numerator, second_denominator);
std::vector<stern_brocot_step> prefix;
for (std::size_t index = 0; index < std::min(first.size(), second.size());
index++) {
if (first[index].first != second[index].first) {
break;
}
prefix.emplace_back(first[index].first,
std::min(first[index].second, second[index].second));
if (first[index].second != second[index].second) {
break;
}
}
return stern_brocot_fraction(prefix);
}
/// @brief Return the open rational interval containing precisely the
/// descendants of a positive reduced fraction. The boundaries may be 0/1 or
/// 1/0.
inline std::pair<std::pair<std::uint64_t, std::uint64_t>,
std::pair<std::uint64_t, std::uint64_t>>
stern_brocot_subtree_range(std::uint64_t numerator, std::uint64_t denominator) {
using u128 = unsigned __int128;
auto path = stern_brocot_path(numerator, denominator);
u128 left_numerator = 0;
u128 left_denominator = 1;
u128 right_numerator = 1;
u128 right_denominator = 0;
for (auto [direction, count] : path) {
if (direction == 'L') {
right_numerator += u128(count) * left_numerator;
right_denominator += u128(count) * left_denominator;
} else {
left_numerator += u128(count) * right_numerator;
left_denominator += u128(count) * right_denominator;
}
}
constexpr u128 maximum = std::numeric_limits<std::uint64_t>::max();
assert(left_numerator <= maximum && left_denominator <= maximum);
assert(right_numerator <= maximum && right_denominator <= maximum);
return {{std::uint64_t(left_numerator), std::uint64_t(left_denominator)},
{std::uint64_t(right_numerator), std::uint64_t(right_denominator)}};
}
} // namespace noya
#endif // NOYA_STERN_BROCOT_HPP
#include <cassert>
#include <cstdint>
#include <limits>
#include <numeric>
#include <optional>
#include <utility>
#include <vector>
/// @complexity Time: O(log max(numerator,denominator)) run steps.
/// Space: O(number of returned runs).
namespace noya {
using stern_brocot_step = std::pair<char, std::uint64_t>;
/// @brief Run-length encoded L/R path from 1/1 to a positive reduced fraction.
inline std::vector<stern_brocot_step>
stern_brocot_path(std::uint64_t numerator, std::uint64_t denominator) {
assert(numerator >= 1 && denominator >= 1);
assert(std::gcd(numerator, denominator) == 1);
std::vector<stern_brocot_step> path;
while (numerator != denominator) {
if (numerator < denominator) {
std::uint64_t count = (denominator - 1) / numerator;
path.emplace_back('L', count);
denominator -= count * numerator;
} else {
std::uint64_t count = (numerator - 1) / denominator;
path.emplace_back('R', count);
numerator -= count * denominator;
}
}
return path;
}
/// @brief Decode a run-length encoded Stern-Brocot path to its reduced
/// fraction; asserts when a 64-bit numerator or denominator would overflow.
inline std::pair<std::uint64_t, std::uint64_t>
stern_brocot_fraction(const std::vector<stern_brocot_step> &path) {
using u128 = unsigned __int128;
u128 left_numerator = 0;
u128 left_denominator = 1;
u128 numerator = 1;
u128 denominator = 1;
u128 right_numerator = 1;
u128 right_denominator = 0;
constexpr u128 maximum = std::numeric_limits<std::uint64_t>::max();
for (auto [direction, count] : path) {
assert((direction == 'L' || direction == 'R') && count >= 1);
if (direction == 'L') {
right_numerator = numerator + u128(count - 1) * left_numerator;
right_denominator = denominator + u128(count - 1) * left_denominator;
numerator += u128(count) * left_numerator;
denominator += u128(count) * left_denominator;
} else {
left_numerator = numerator + u128(count - 1) * right_numerator;
left_denominator = denominator + u128(count - 1) * right_denominator;
numerator += u128(count) * right_numerator;
denominator += u128(count) * right_denominator;
}
assert(numerator <= maximum && denominator <= maximum);
}
return {std::uint64_t(numerator), std::uint64_t(denominator)};
}
/// @brief Return the Stern-Brocot node at the given depth on the root-to-node
/// path, or nullopt when the requested depth is greater than the node depth.
inline std::optional<std::pair<std::uint64_t, std::uint64_t>>
stern_brocot_ancestor(std::uint64_t depth, std::uint64_t numerator,
std::uint64_t denominator) {
auto path = stern_brocot_path(numerator, denominator);
std::vector<stern_brocot_step> prefix;
for (auto [direction, count] : path) {
std::uint64_t used = std::min(depth, count);
if (used > 0) {
prefix.emplace_back(direction, used);
}
if (depth <= count) {
return stern_brocot_fraction(prefix);
}
depth -= count;
}
if (depth == 0) {
return stern_brocot_fraction(prefix);
}
return std::nullopt;
}
/// @brief Return the lowest common ancestor of two positive reduced fractions.
inline std::pair<std::uint64_t, std::uint64_t>
stern_brocot_lca(std::uint64_t first_numerator, std::uint64_t first_denominator,
std::uint64_t second_numerator,
std::uint64_t second_denominator) {
auto first = stern_brocot_path(first_numerator, first_denominator);
auto second = stern_brocot_path(second_numerator, second_denominator);
std::vector<stern_brocot_step> prefix;
for (std::size_t index = 0; index < std::min(first.size(), second.size());
index++) {
if (first[index].first != second[index].first) {
break;
}
prefix.emplace_back(first[index].first,
std::min(first[index].second, second[index].second));
if (first[index].second != second[index].second) {
break;
}
}
return stern_brocot_fraction(prefix);
}
/// @brief Return the open rational interval containing precisely the
/// descendants of a positive reduced fraction. The boundaries may be 0/1 or
/// 1/0.
inline std::pair<std::pair<std::uint64_t, std::uint64_t>,
std::pair<std::uint64_t, std::uint64_t>>
stern_brocot_subtree_range(std::uint64_t numerator, std::uint64_t denominator) {
using u128 = unsigned __int128;
auto path = stern_brocot_path(numerator, denominator);
u128 left_numerator = 0;
u128 left_denominator = 1;
u128 right_numerator = 1;
u128 right_denominator = 0;
for (auto [direction, count] : path) {
if (direction == 'L') {
right_numerator += u128(count) * left_numerator;
right_denominator += u128(count) * left_denominator;
} else {
left_numerator += u128(count) * right_numerator;
left_denominator += u128(count) * right_denominator;
}
}
constexpr u128 maximum = std::numeric_limits<std::uint64_t>::max();
assert(left_numerator <= maximum && left_denominator <= maximum);
assert(right_numerator <= maximum && right_denominator <= maximum);
return {{std::uint64_t(left_numerator), std::uint64_t(left_denominator)},
{std::uint64_t(right_numerator), std::uint64_t(right_denominator)}};
}
} // namespace noya