rational_reconstruction.hpp¶
从模意义下的剩余值，在给定分子分母界内唯一恢复最简分数。
\[
\displaystyle \frac{x}{y}\equiv r\pmod m
\]
Complexity: Time: O(log m). Space: O(1).
Implementation¶
当前头文件，省略 include guard；依赖见 #include。
/// @complexity Time: O(log m).
/// Space: O(1).
#include "noya/fraction.hpp"
#include <cassert>
#include <cstdlib>
#include <cstdint>
#include <numeric>
#include <optional>
namespace noya {
/// @brief Reconstruct a reduced fraction numerator/denominator congruent to a
/// residue res modulo mod within supplied bounds. Uniqueness is guaranteed when
/// 2*nb*db < mod.
inline std::optional<fraction> rational_reconstruction(std::int64_t res,
std::int64_t mod,
std::int64_t nb,
std::int64_t db) {
assert(mod >= 2);
assert(nb >= 0 && db >= 1);
res %= mod;
if (res < 0) {
res += mod;
}
std::int64_t r0 = mod;
std::int64_t rem = res;
std::int64_t d0 = 0;
std::int64_t den = 1;
while (std::llabs(rem) > nb) {
if (rem == 0) {
return std::nullopt;
}
std::int64_t quo = r0 / rem;
__int128 r1 = __int128(r0) - __int128(quo) * rem;
__int128 d1 = __int128(d0) - __int128(quo) * den;
r0 = rem;
rem = std::int64_t(r1);
d0 = den;
den = std::int64_t(d1);
}
if (den < 0) {
den = -den;
rem = -rem;
}
if (den == 0 || den > db || std::gcd(rem, den) != 1 ||
std::gcd(den, mod) != 1) {
return std::nullopt;
}
__int128 con = __int128(res) * den - rem;
if (con % mod != 0) {
return std::nullopt;
}
return fraction(rem, den);
}
} // namespace noya
#ifndef NOYA_RATIONAL_RECONSTRUCTION_HPP
#define NOYA_RATIONAL_RECONSTRUCTION_HPP 1
/// @complexity Time: O(log m).
/// Space: O(1).
#include "noya/fraction.hpp"
#include <cassert>
#include <cstdlib>
#include <cstdint>
#include <numeric>
#include <optional>
namespace noya {
/// @brief Reconstruct a reduced fraction numerator/denominator congruent to a
/// residue res modulo mod within supplied bounds. Uniqueness is guaranteed when
/// 2*nb*db < mod.
inline std::optional<fraction> rational_reconstruction(std::int64_t res,
std::int64_t mod,
std::int64_t nb,
std::int64_t db) {
assert(mod >= 2);
assert(nb >= 0 && db >= 1);
res %= mod;
if (res < 0) {
res += mod;
}
std::int64_t r0 = mod;
std::int64_t rem = res;
std::int64_t d0 = 0;
std::int64_t den = 1;
while (std::llabs(rem) > nb) {
if (rem == 0) {
return std::nullopt;
}
std::int64_t quo = r0 / rem;
__int128 r1 = __int128(r0) - __int128(quo) * rem;
__int128 d1 = __int128(d0) - __int128(quo) * den;
r0 = rem;
rem = std::int64_t(r1);
d0 = den;
den = std::int64_t(d1);
}
if (den < 0) {
den = -den;
rem = -rem;
}
if (den == 0 || den > db || std::gcd(rem, den) != 1 ||
std::gcd(den, mod) != 1) {
return std::nullopt;
}
__int128 con = __int128(res) * den - rem;
if (con % mod != 0) {
return std::nullopt;
}
return fraction(rem, den);
}
} // namespace noya
#endif // NOYA_RATIONAL_RECONSTRUCTION_HPP
#include <cassert>
#include <compare>
#include <cstdint>
#include <cstdlib>
#include <numeric>
#include <optional>
/// @complexity Time: O(log m).
/// Space: O(1).
/// @complexity Time: O(log max(|num|,den)) normalization; arithmetic itself is O(1).
/// Space: O(1).
namespace noya {
/// @brief Normalized signed 64-bit rational number; arithmetic results must fit
/// int64, while comparison and intermediate products use signed 128-bit.
struct fraction {
std::int64_t num = 0;
std::int64_t den = 1;
fraction() = default;
fraction(std::int64_t nx, std::int64_t dx = 1) : num(nx), den(dx) {
normalize();
}
void normalize() {
assert(den != 0);
if (den < 0) {
num = -num;
den = -den;
}
std::int64_t div = std::gcd(num, den);
num /= div;
den /= div;
}
friend bool operator==(const fraction &, const fraction &) = default;
friend std::strong_ordering operator<=>(const fraction &a,
const fraction &b) {
__int128 l = __int128(a.num) * b.den;
__int128 r = __int128(b.num) * a.den;
return l < r ? std::strong_ordering::less
: l > r ? std::strong_ordering::greater
: std::strong_ordering::equal;
}
friend fraction operator+(const fraction &a, const fraction &b) {
std::int64_t div = std::gcd(a.den, b.den);
__int128 num =
__int128(a.num) * (b.den / div) + __int128(b.num) * (a.den / div);
__int128 den = __int128(a.den / div) * b.den;
return {std::int64_t(num), std::int64_t(den)};
}
friend fraction operator-(const fraction &a, const fraction &b) {
return a + fraction(-b.num, b.den);
}
friend fraction operator*(fraction a, fraction b) {
std::int64_t g1 = std::gcd(a.num < 0 ? -a.num : a.num, b.den);
std::int64_t g2 = std::gcd(b.num < 0 ? -b.num : b.num, a.den);
a.num /= g1;
b.den /= g1;
b.num /= g2;
a.den /= g2;
return {std::int64_t(__int128(a.num) * b.num),
std::int64_t(__int128(a.den) * b.den)};
}
friend fraction operator/(const fraction &a, const fraction &b) {
assert(b.num != 0);
return a * fraction(b.den, b.num);
}
friend fraction operator-(const fraction &val) { return {-val.num, val.den}; }
};
} // namespace noya
namespace noya {
/// @brief Reconstruct a reduced fraction numerator/denominator congruent to a
/// residue res modulo mod within supplied bounds. Uniqueness is guaranteed when
/// 2*nb*db < mod.
inline std::optional<fraction> rational_reconstruction(std::int64_t res,
std::int64_t mod,
std::int64_t nb,
std::int64_t db) {
assert(mod >= 2);
assert(nb >= 0 && db >= 1);
res %= mod;
if (res < 0) {
res += mod;
}
std::int64_t r0 = mod;
std::int64_t rem = res;
std::int64_t d0 = 0;
std::int64_t den = 1;
while (std::llabs(rem) > nb) {
if (rem == 0) {
return std::nullopt;
}
std::int64_t quo = r0 / rem;
__int128 r1 = __int128(r0) - __int128(quo) * rem;
__int128 d1 = __int128(d0) - __int128(quo) * den;
r0 = rem;
rem = std::int64_t(r1);
d0 = den;
den = std::int64_t(d1);
}
if (den < 0) {
den = -den;
rem = -rem;
}
if (den == 0 || den > db || std::gcd(rem, den) != 1 ||
std::gcd(den, mod) != 1) {
return std::nullopt;
}
__int128 con = __int128(res) * den - rem;
if (con % mod != 0) {
return std::nullopt;
}
return fraction(rem, den);
}
} // namespace noya