数论基础
数论函数
数论函数是指这样一类函数：其定义域是正整数，值域是一个数集。
定义两个数论函数的加法，为逐项相加，即 \((f + g)(n) = f(n) + g(n)\) 。
定义数乘这个数和每一项都相乘，即 \((xf)(n) = x \times f(n)\) 。
常见数论函数
\[1(x) = 1 \\ \epsilon(x) = [x = 1] \\ id(x) = x \\ id^k(x) = x^k \\ \sigma_0(x) = \sum_{d | x} 1 \\ \sigma(x) = \sum_{d | x} d \\ \sigma_k(x) = \sum_{d | x} d^k \]常见数论函数的性质
\[\sum_{i = 1}^n \sigma(i) = \sum_{i = 1}^n \lfloor \dfrac{n}{i} \rfloor \]\[\sigma(ij) = \sum_{x | i} \sum_{y | j} [\gcd(x, y) = 1] \]证明：考虑每个因数的贡献即可。
证明：考虑把每个因子一一映射，如果 \(ij\) 的因子 \(k\) 中有一个因子 \(p^c\)，\(i\) 中有因子 \(p^a\)，\(j\) 中有因子 \(p^b\)。我们规定：
- 如果 \(c\leqslant a\)，那么在 \(i\) 中选择。
- 如果 \(c>a\)，那么我们把 \(c\) 减去 \(a\)，在 \(j\) 中选择 \(p^{c-a}\)（在 \(j\) 中选择 \(p^e\) 表示的是 \(p^{a+e}\)）
对于 \(ij\) 的因子 \(k\) 的其他因子同理。于是对于任何一个 \(k\) 有一个唯一的映射，且每一个选择对应着唯一的 \(k\)。
于是对于 \(ij\) 的因子 \(k=\prod {p_i}^{c_i}\)，我们不可能同时在 \(i\) 和 \(j\) 中选择 \(p_i\)（优先在 \(i\) 中选择，如果不够就只在 \(j\) 中选择不够的指数），故 \(x\) 和 \(y\) 必须互质。
推广：
\[\sigma(ijk) = \sum_{a | i} \sum_{b | j} \sum_{c | k} [\gcd(a, b) = 1] [\gcd(b, c) = 1] [\gcd(a, c) = 1] \]积性函数
定义
积性函数
对于函数 \(f(x)\) ，若 \(\forall \gcd(x, y) = 1, f(x y) = f(x) f(y)\) ，则称 \(f(x)\) 为积性函数。
常见的积性函数：
- 约数函数：\(\sigma_k(n) = \sum_{d | n} d^k\) 。
- 欧拉函数：\(\varphi(n)\) 。
- 莫比乌斯函数：\(\mu(n)\)
完全积性函数
对于函数 \(f(x)\) ，若 \(\forall x, y, f(x y) = f(x) f(y)\) ，则称 \(f(x)\) 为完全积性函数。不难发现完全积性函数属于积性函数。
常见的完全积性函数：
- 常数函数： \(1(n) = 1\) 。
- 元函数： \(\epsilon(n) = [n = 1]\) 。
- 单位函数：\(id(n) = n\) 。
- 幂函数： \(id^x(n) = n^x\) 。
性质
一个显然的事实：\(f(1) = 1\) 。
设 \(n\) 的质因数分解为 \(\prod p_i^{c_i}\) ，则：
- 对于积性函数 \(f(x)\) ，有 \(f(n) = \prod f(p_i^{c_i})\) 。
- 对于完全积性函数 \(f(x)\) ，有 \(f(n) = \prod f(p_i)^{k_i}\) 。
若 \(f(n), g(n)\) 均为积性函数，则以下函数也为积性函数：
\[h(n) = f(n^p) \\ h(n) = f^p(n) \\ h(n) = f(n) g(n) \\ h(n) = \sum_{d | n} f(d) g(\dfrac{n}{d}) \]线性筛
线性筛的时间复杂度是由每个数都会被最小的质因子标记的。
一般的，对于积性函数 \(f(x)\) ，若 \(f(p^c)\) 的值便于分析，则 \(f(1 \sim n)\) 的值可以被线性筛筛出来。
记 \(low_i\) 表示：若 \(i\) 的最小质因子为 \(p_1\) ，其次数为 \(c_1\) ，则 \(low_i = p_1^{c_1}\) 。
在线性筛时，若当前筛出了一个质数 \(p\) ，则将 \(f(p), f(p^2), \cdots, f(p^c)\) 都算出来，否则可以得到 \(f(i) = f(\dfrac{i}{low_i}) \times f(low_i)\) 。
inline void sieve(int n) {
f[1] = 1;
for (int i = 2; i <= n; ++i) {
if (!low[i]) {
pri[++pcnt] = i;
ll x = i;
int c = 1;
while (x <= n) {
low[x] = x;
// calculate f(x)
x *= i, c ++;
}
}
for (int j = 1; j <= pcnt && i * pri[j] <= n; ++j) {
int x = i * pri[j];
low[x] = i % pri[j] ? pri[j] : low[i] * pri[j];
f[x] = f[x / low[x]] * f[low[x]];
if (!(i % pri[j]))
break;
}
}
}
Dirichlet卷积
定义：
-
Dirichlet卷积：\(f(n), g(n)\) 的狄利克雷卷积为 \((f * g) (n) = \sum_{d | n} f(d) g(\dfrac{n}{d})\) 。
-
逆：满足 \(\epsilon = g * f\) 时，\(g, f\) 互逆，对于每个 \(f(1) \neq 0\) 的函数都存在逆元，可以构造：
\[g(n) = \dfrac{1}{f(1)} ([n = 1] - \sum_{i | n, i \neq 1} f(i) g(\dfrac{n}{i})) \]
性质：
- 交换律：\(f * g = g * f\) 。
- 结合律：\((f * g) * h = f * (g * h)\) 。
- 分配律：\((f + g) * h = f * h + g * h\) 。
- 数乘：\((xf) * g = x(f * g)\) 。
- 单位元：\(\epsilon * f = f\) 。
- 若 \(h(1) \neq 0\) ，则 \(f = g \iff f * h = g * h\) 。
- 积性函数的卷积、逆仍为积性函数。
常见Dirichlet卷积式：
\[\sigma_0 = 1 * 1 \\ id = \varphi * 1 \\ \sigma = 1 * id = \sigma_0 * \varphi \\ 1 = \sigma_0 * \mu \\ \epsilon = \mu * 1 \\ \varphi = \mu * id \]线性筛
线性筛的时间复杂度是由每个数都会被最小的质因子标记的。
一般的，对于积性函数 \(f(x)\) ，若 \(f(p^c)\) 的值便于分析，则 \(f(1 \sim n)\) 的值可以被线性筛筛出来。
记 \(low_i\) 表示：若 \(i\) 的最小质因子为 \(p_1\) ，其次数为 \(c_1\) ，则 \(low_i = p_1^{c_1}\) 。
在线性筛时，若当前筛出了一个质数 \(p\) ，则将 \(f(p), f(p^2), \cdots, f(p^c)\) 都算出来，否则可以得到 \(f(i) = f(\dfrac{i}{low_i}) \times f(low_i)\) 。
筛素数
inline void sieve(int n) {
memset(isp, true, sizeof(isp));
isp[1] = false;
for (int i = 2; i <= n; ++i) {
if (isp[i])
pri[++pcnt] = i;
for (int j = 1; j <= pcnt && i * pri[j] <= n; ++j) {
isp[i * pri[j]] = false;
if (!(i % pri[j]))
break;
}
}
}
筛约数个数
\(d_i\) 表示 \(i\) 的约数个数， \(c_i\) 表示 \(i\) 的最小质因子出现次数
inline void prework(int n) {
memset(isp, true, sizeof(isp));
pcnt = 0;
isp[1] = false, d[1] = 1;
for (int i = 2; i <= n; ++i) {
if (isp[i])
pri[++pcnt] = i, d[i] = 2, c[i] = 1;
for (int j = 1; j <= pcnt && i * pri[j] <= n; ++j) {
isp[i * pri[j]] = false;
if (i % pri[j]) {
c[i * pri[j]] = 1;
d[i * pri[j]] = d[i] * 2;
} else {
c[i * pri[j]] = c[i] + 1;
d[i * pri[j]] = d[i] / c[i * pri[j]] * (c[i * pri[j]] + 1);
break;
}
}
}
}
筛约数和
\(f_i\) 表示 \(i\) 的约数和， \(g_i\) 表示 \(i\) 的最小质因子的 \(p^0 + p^1 + \cdots + p^k\)
inline void init(int n) {
memset(IsPrime, true, sizeof(IsPrime));
tot = 0;
g[1] = f[1] = 1;
for (int i = 2; i <= n; ++i) {
if (IsPrime[i])
prime[++tot] = i, g[i] = i + 1, f[i] = i + 1;
for (int j = 1; j <= tot && i * prime[j] <= n; ++j) {
IsPrime[i * prime[j]] = false;
if (i % prime[j]) {
f[i * prime[j]] = f[i] * f[prime[j]];
g[i * prime[j]] = prime[j] + 1;
}
else {
g[i * prime[j]] = g[i] * prime[j] + 1;
f[i * prime[j]] = f[i] / g[i] * g[i * prime[j]];
break;
}
}
}
}
扩展欧几里得
扩展欧几里得算法常用于求解方程 \(ax + by = \gcd(a, b)\) 的一组可行解。由裴蜀定理可知，该方程一定有解。
实现
设：
\[ax_1 + by_1 = \gcd(a,b) \\ bx_2 + (a \bmod b)y_2 = \gcd(b, a \bmod b) \]由欧几里得定理 \(\gcd(a, b) = \gcd(b, a \bmod b)\) 可知：
\[ax_1 + by_1 = bx_2 + (a \bmod b)y_2 \]将 \(a \bmod b\) 用 \(a - \lfloor \dfrac{a}{b} \rfloor \times b\) 带入，得：
\[ax_1 + by_1 = bx_2 + (a - \lfloor \dfrac{a}{b} \rfloor \times b)y_2 \]展开得到：
\[ax_1 + by_1 = ay_2 + bx_2 - \lfloor \dfrac{a}{b} \rfloor \times by_2 = ay_2 + b(x_2 - \lfloor \dfrac{a}{b} \rfloor y_2) \]所以：
\[\begin{cases} x_1 = y_2 \\ y_1 = x_2 - \lfloor \dfrac{a}{b} \rfloor y_2 \end{cases} \]实现：
int exgcd(int a, int b, int &x, int &y) {
b ? (exgcd(b, a % b, y, x), y -= a / b * x) : (x = 1, y = 0);
}
应用
对于不定方程：
\[ax + by = c \]若 \(c\) 不是 \(\gcd(a, b)\) 的倍数，则无解。已知exgcd可以求出
\[ax + by = \gcd(a,b) \]的一组整数特解，记为 \(\begin{cases} x = x_0 \\ y = y_0 \end{cases}\) ，可以得到该方程的一特解：
\[\begin{cases} x_1 = \dfrac{x_0 c}{\gcd(a, b)} \\ y_1 = \dfrac{y_0 c}{\gcd(a, b)} \end{cases} \]考虑构造原方程整数通解形式，设 \(d \in \mathbb{Q}\) ，则：
\[a(x_1 + db) + b(y_1 - da) = c \]同时，通解需保证 \(x_1 + db, y_1 - da\) 均为整数。因为 \(x_1, y_1\) 均为整数，故只需保证 \(db, da\) 为整数即可。
令当 \(d\) 取到最小可能的正值时的 \(d_x = db, d_y = da\) ，那么有通解形式：
\[\begin{cases} x = x_1 + s d_x \\ y = y_1 - s d_y \end{cases} \]其中 \(s\) 为任意整数。
中国剩余定理
中国剩余定理（CRT）是用于求解形如：
\[\begin{aligned} \begin{cases} x &\equiv a_1 \pmod{m_1} \\ x &\equiv a_2 \pmod{m_2} \\ &\vdots \\ x &\equiv a_k \pmod{m_k} \end{cases} \end{aligned} \]的一元线性同余方程组，其中 \(m\) 两两互质。
实现
- 计算 \(M = \prod_{i = 1}^k m_i\) 。
- 对于第 \(i\) 个方程
- 计算 \(p_i = \dfrac{M}{m_i}\) 。
- 计算 \(p_i\) 在模 \(m_i\) 意义下的逆元 \(p_i^{-1}\) 。
- 计算 \(c_i = p_i p_i^{-1}\) （不要对 \(m_i\) 取模）。
- 方程组在模 \(M\) 意义下的唯一解即为 \(x = \sum_{i = 1}^k a_i c_i \pmod{M}\) 。
inline ll CRT() {
ll M = 1, res = 0;
for (int i = 1; i <= n; ++i)
M *= m[i];
for (int i = 1; i <= n; ++i)
res = (res + a[i] * (M / m[i]) % M * inv(M / m[i], m[i]) % M) % M;
return res;
}
证明：只需证明上述算法解出的 \(x\) 均满足每一个方程即可。
当 \(i \not = j\) 时，显然有 \(p_j \equiv 0 \pmod{m_i}\) ，故 \(c_j \equiv p_j \equiv 0 \pmod{m_i}\) ，又因为 \(c_i \equiv p_i (p_i^{-1} \bmod m_i) \equiv 1 \pmod{m_i}\) ，所以有：
\[\begin{aligned} x &\equiv \sum_{j = 1}^k a_j c_j &\pmod{m_i} \\ &\equiv a_i c_i &\pmod{m_i} \\ &\equiv a_i p_i (p_i^{-1} \bmod m_i) &\pmod{m_i} \\ &\equiv a_i &\pmod{m_i} \end{aligned} \]所以上述算法解出的 \(x\) 均满足每一个方程。
另外，若 \(y \not = x\) ，则总存在 \(i\) 使得 \(x, y\) 在模 \(m_i\) 意义下不同余。
Garner算法
CRT的另一个用途是用一组比较小的质数表示一个大整数。
例如，若 \(a\) 满足：
\[\begin{aligned} \begin{cases} A &\equiv a_1 \pmod{p_1} \\ A &\equiv a_2 \pmod{p_2} \\ &\vdots \\ A &\equiv a_k \pmod{p_k} \end{cases} \end{aligned} \]且 \(A < \prod_{i = 1}^k p_i\) ， \(p_i\) 均为质数，则有
\[A = x_1 + x_2 p_1 + x_3p_1p_2 + \cdots + x_kp_1 \cdots p_{k - 1} \]可以用 Garner 算法求出系数 \(x_1, \cdots, x_k\)
令 \(r_{i, j}\) 表示 \(p_i\) 在模 \(p_j\) 意义下的逆元，第一个方程带入 \(A\) 可以得到
\[a_1 \equiv x_1 \pmod{p_1} \]带入第二个方程可以得到
\[\begin{aligned} a_2 &\equiv x_1 + x_2 p_1 &\pmod{p_2} \\ a_2 - x_1 &\equiv x_2 p_1 &\pmod{p_2} \\ (a_2 - x_1) r_{1, 2} &\equiv x_2 &\pmod{p_2} \\ x_2 &\equiv (a_2 - x_1) r_{1, 2} &\pmod{p_2} \end{aligned} \]类似地，我们可以得到
\[x_k = (\cdots((a_k - x_1)r_{1, k} - x_2)r_{2, k} - \cdots) r_{k - 1, k} \pmod{p_k} \]代码实现：
for (int i = 1; i <= k; ++i) {
x[i] = a[i];
for (int j = 1; j < i; ++j)
x[i] = (r[j][i] * (x[i] - x[j]) % p[i] + p[i]) % p[i];
}
某些问题出于某些原因，给出的模数不是质数。但对其质因数分解会发现它没有平方因子，也就是该模数是由一些不重复的质数相乘得到。那么我们可以分别对这些模数进行计算，最后用 CRT 合并答案。
优化
观察：
\[\begin{aligned} x_1 &= a_1 &\pmod{p_1} \\ x_2 &= (a_2 - x_1) r_{1, 2} \\ &= a_2r_{1, 2} - x_1 r_{1, 2} &\pmod{p_2} \\ x_3 &= ((a_3 - x_1) r_{1, 3} - x_2) r_{2, 3} \\ &= a_3 r_{1, 3} r_{2, 3} - (x_1 r_{1, 3} r_{2, 3} + x_2 r_{2, 3}) &\pmod{p_3} \\ x_4 &= (((a_3 - x_1) r_{1, 4} - x_2) r_{2, 4} - x_3) r_{3, 4} \\ &= a_3 r_{1, 4} r_{2, 4} r_{3, 4} - (x_1 r_{1, 4} r_{2, 4} r_{3, 4} + x_2 r_{2, 4} r_{3, 4} + x_3 r_{3, 4}) &\pmod{p_4} \end{aligned} \]归纳：
\[\begin{aligned} x_k &= a_k \prod_{i = 1}^{k - 1} r_{i, k} - \sum_{i = 1}^{k - 1} (x_i \prod_{j = i}^{k - 1} r_{j, k}) &\pmod{p_k} \\ &= \prod_{i = 1}^{k - 1} r_{i, k} (a_k - \sum_{i = 1}^{k - 1} \dfrac{x_i}{\prod_{j = 1}^{i - 1} r_{j, k}}) &\pmod{p_k} \\ &= (\prod_{i = 1}^{k - 1} p_i)^{-1} (a_k - \sum_{i = 1}^{k - 1} (x_i \prod_{j = 1}^{i - 1} p_j)) &\pmod{p_k} \end{aligned} \]那么有：
\[\begin{aligned} A &= \sum_{i = 1}^{k} (x_i \prod_{j = 1}^{i - 1} p_j) &\pmod{M} \\ A &= \sum_{i = 1}^{k} ((\prod_{j = 1}^{i - 1} p_j)^{-1} (a_i - \sum_{j = 1}^{i - 1} (x_j \prod_{u = 1}^{j - 1} p_u)) \prod_{j = 1}^{i - 1} p_j) &\pmod{M} \\ A &= \sum_{i = 1}^{k} (a_i - \sum_{j = 1}^{i - 1} (x_j \prod_{u = 1}^{j - 1} p_u)) &\pmod{M} \end{aligned} \]for (int i = 1; i <= k; ++i)
A = (A + 1ll * ((a[i] - A) % p[i] + p[i]) % p[i] * inv(mod % p[i], p) % p[i] * Mod) % Mod;
此时时间复杂降为 \(O(n)\) 。
exCRT
模数不互质的CRT。
考虑两个方程的情况，设两个方程分别为 \(x \equiv a_1 \pmod{m_1}\) 与 \(x \equiv a_2 \pmod{m_2}\) 。转化为不定方程 \(x = m_1 p + a_1 = m_2 q + a_2\) ，其中 \(p, q \in \mathbb{Z}\) ，则有 $ m_1 p - m_2 q = a_2 - a_1$
由裴蜀定理，当 \(\gcd(m_1, m_2) \nmid (a_2 - a_1)\) 时方程无解。否则我们用exgcd求出一组可行解 \(p, q\) ，则 \(x \equiv m_1 p + a_1 \pmod{\operatorname{lcm}(m_1, m_2)}\) 。
拓展到多个方程，类似地，我们按照上述过程两两合并即可。
inline pair<s128, s128> merge(pair<s128, s128> dt1, pair<s128, s128> dt2) {
s128 m1 = dt1.first, a1 = dt1.second;
s128 m2 = dt2.first, a2 = dt2.second;
s128 c = ((a2 - a1) % m2 + m2) % m2;
s128 x, y, g = exgcd(m1, m2, x, y), bg = m2 / g;
if (c % g)
return make_pair((s128) -1, (s128) -1);
x = x * c / g % bg, a1 += x * m1;
m1 *= bg, a1 = (a1 % m1 + m1) % m1;
return make_pair(m1, a1);
}
inline s128 exCRT() {
pair<s128, s128> now = make_pair(m[1], a[1]);
for (int i = 2; i <= n; ++i)
now = merge(now, make_pair(m[i], a[i]));
return now.second;
}
数论分块
注意到 \(\lfloor \dfrac{n}{x} \rfloor\) 在 \(x \in [1, n]\) 时的取值不超过 \(2 \sqrt{n}\) 种，对于某一些问题，考虑对一段 \(\lfloor \dfrac{n}{x} \rfloor\) 均相等的区间，我们可以直接得出该区间里的所有 \(x\) 对答案的贡献。
如果要实现整块一起统计，我们需要求出每一块的块头 \(l\) 和块尾 \(r\)，则：
\[Ans = \sum_{i=1}^{k} \lfloor \dfrac{n}{i} \rfloor = \sum_{[l,r]} (r-l+1)(\lfloor \dfrac{n}{l} \rfloor) \]注意到每一块的 \(l\) 都可以由上一块的 \(r\) 推出，故我们继续讨论如何在已知 \(l\) 的情况下推出 \(r\) 。令 \(t = \lfloor \dfrac{n}{l} \rfloor\) ，容易得到
\[\begin{cases} r = \min(\lfloor \dfrac{n}{t} \rfloor,n) \ \ \ &(t \neq 0) \\ r = n \ \ \ &(t=0) \end{cases} \]直接计算即可，时间复杂度 \(O(\sqrt{n})\) 。
for (int l = 1, r; l <= n; l = r + 1) {
r = n / (n / l);
// Calculate the contribution of [l, r]
}
扩展到二维甚至多维也是类似的，时间复杂度 \(O(k \sqrt{n})\) ，其中 \(k\) 为维度。
二维数论分块：
for (int l = 1, r; l <= min(n, m); l = r + 1) {
r = min(n / (n / l), m / (m / l));
// Calculate the contribution of [l, r]
}