商の列挙
(number-theory/enumerate-quotients.hpp)
- View this file on GitHub
- Last update: 2025-10-17 21:43:09+09:00
- Include:
#include "number-theory/enumerate-quotients.hpp"
商の列挙
正整数 $N$ に対し,集合 $Q_N$ を
\[Q_N=\left\{\left\lfloor\frac{N}{x}\right\rfloor : x\in\mathbb{Z},1\leq x\leq N\right\}\]| で定める.$ | Q_N | =O(\sqrt{N})$ である. |
table(N)
$Q_N$ の元を昇順に並べた列を求める.
get_range(N, y)
$y\in Q_N$ に対し,$\left\lfloor\frac{N}{x}\right\rfloor=y$ となる $x$ の区間 $(l,r]$ を求める.
左側が開であることに注意.
iterate(N, f)
$q\in Q_N$ について昇順に,$(l,r]=\operatorname{get_range}(q)$ として関数 f(q,l,r) を呼び出す.
Depends on
Verified with
Code
#pragma once
#include "math/util.hpp"
namespace EnumerateQuotients {
using i64 = int64_t;
i64 div(i64 a, i64 b) { return double(a) / b; };
vector<i64> table(i64 N) {
i64 sq = Math::isqrt(N);
vector<i64> xs(sq);
iota(xs.begin(), xs.end(), 1);
if (N <= 1e12) {
for (i64 i = div(N, sq + 1); i > 0; i--) xs.push_back(div(N, i));
} else {
for (i64 i = N / (sq + 1); i > 0; i--) xs.push_back(N / i);
}
return xs;
}
pair<i64, i64> get_range(i64 N, i64 q) {
return N <= 1e12 ? pair<i64, i64>{div(N, q + 1), div(N, q)} : pair<i64, i64>{N / (q + 1), N / q};
}
template <class F>
void iterate(i64 N, F f) {
i64 sq = Math::isqrt(N);
vector<i64> xs;
if (N <= 1e12) {
i64 x = N;
for (i64 q = 1; x <= sq; q++) {
i64 y = div(N, q + 1);
f(q, y, x);
x = y;
}
for (; x > 0; x--) f(div(N, x), x - 1, x);
} else {
i64 x = N;
for (i64 q = 1; x <= sq; q++) {
i64 y = N / (q + 1);
f(q, y, x);
x = y;
}
for (; x > 0; x--) f(N / x, x - 1, x);
}
}
}; // namespace EnumerateQuotients
/**
* @brief 商の列挙
* @docs docs/number-theory/enumerate-quotients.md
*/#line 2 "number-theory/enumerate-quotients.hpp"
#line 2 "math/util.hpp"
namespace Math {
template <class T>
T safe_mod(T a, T b) {
assert(b != 0);
if (b < 0) a = -a, b = -b;
a %= b;
return a >= 0 ? a : a + b;
}
template <class T>
T floor(T a, T b) {
assert(b != 0);
if (b < 0) a = -a, b = -b;
return a >= 0 ? a / b : (a + 1) / b - 1;
}
template <class T>
T ceil(T a, T b) {
assert(b != 0);
if (b < 0) a = -a, b = -b;
return a > 0 ? (a - 1) / b + 1 : a / b;
}
long long isqrt(long long n) {
if (n <= 0) return 0;
long long x = sqrt(n);
while ((x + 1) * (x + 1) <= n) x++;
while (x * x > n) x--;
return x;
}
// return g=gcd(a,b)
// a*x+b*y=g
// - b!=0 -> 0<=x<|b|/g
// - b=0 -> ax=g
template <class T>
T ext_gcd(T a, T b, T& x, T& y) {
T a0 = a, b0 = b;
bool sgn_a = a < 0, sgn_b = b < 0;
if (sgn_a) a = -a;
if (sgn_b) b = -b;
if (b == 0) {
x = sgn_a ? -1 : 1;
y = 0;
return a;
}
T x00 = 1, x01 = 0, x10 = 0, x11 = 1;
while (b != 0) {
T q = a / b, r = a - b * q;
x00 -= q * x01;
x10 -= q * x11;
swap(x00, x01);
swap(x10, x11);
a = b, b = r;
}
x = x00, y = x10;
if (sgn_a) x = -x;
if (sgn_b) y = -y;
if (b0 != 0) {
a0 /= a, b0 /= a;
if (b0 < 0) a0 = -a0, b0 = -b0;
T q = x >= 0 ? x / b0 : (x + 1) / b0 - 1;
x -= b0 * q;
y += a0 * q;
}
return a;
}
template <class T>
T inv_mod(T x, T m) {
x %= m;
if (x < 0) x += m;
T a = m, b = x;
T y0 = 0, y1 = 1;
while (b > 0) {
T q = a / b;
swap(a -= q * b, b);
swap(y0 -= q * y1, y1);
}
if (y0 < 0) y0 += m / a;
return y0;
}
template <class T>
T pow_mod(T x, T n, T m) {
x = (x % m + m) % m;
T y = 1;
while (n) {
if (n & 1) y = y * x % m;
x = x * x % m;
n >>= 1;
}
return y;
}
constexpr long long pow_mod_constexpr(long long x, long long n, int m) {
if (m == 1) return 0;
unsigned int _m = (unsigned int)(m);
unsigned long long r = 1;
unsigned long long y = x % m;
if (y >= m) y += m;
while (n) {
if (n & 1) r = (r * y) % _m;
y = (y * y) % _m;
n >>= 1;
}
return r;
}
}; // namespace Math
#line 4 "number-theory/enumerate-quotients.hpp"
namespace EnumerateQuotients {
using i64 = int64_t;
i64 div(i64 a, i64 b) { return double(a) / b; };
vector<i64> table(i64 N) {
i64 sq = Math::isqrt(N);
vector<i64> xs(sq);
iota(xs.begin(), xs.end(), 1);
if (N <= 1e12) {
for (i64 i = div(N, sq + 1); i > 0; i--) xs.push_back(div(N, i));
} else {
for (i64 i = N / (sq + 1); i > 0; i--) xs.push_back(N / i);
}
return xs;
}
pair<i64, i64> get_range(i64 N, i64 q) {
return N <= 1e12 ? pair<i64, i64>{div(N, q + 1), div(N, q)} : pair<i64, i64>{N / (q + 1), N / q};
}
template <class F>
void iterate(i64 N, F f) {
i64 sq = Math::isqrt(N);
vector<i64> xs;
if (N <= 1e12) {
i64 x = N;
for (i64 q = 1; x <= sq; q++) {
i64 y = div(N, q + 1);
f(q, y, x);
x = y;
}
for (; x > 0; x--) f(div(N, x), x - 1, x);
} else {
i64 x = N;
for (i64 q = 1; x <= sq; q++) {
i64 y = N / (q + 1);
f(q, y, x);
x = y;
}
for (; x > 0; x--) f(N / x, x - 1, x);
}
}
}; // namespace EnumerateQuotients
/**
* @brief 商の列挙
* @docs docs/number-theory/enumerate-quotients.md
*/