「学习笔记」单位根反演

单位根反演

公式推导

[n ∣ k] = \frac{1}{n} \sum_{i = 0}^{n - 1} ω_{n}^{i k}

$[ n \mid k ] = \frac{1}{n}\sum\limits_{ i = 0 }^{ n - 1 } \omega_{ n }^{ ik }$

证明：

当 $n\mid k$ 时， $\omega_{n}^{ik}=1$ ，显然成立；

当 $n\nmid k$ 时，对右边项进行等比数列求和得到 $\frac{1}{n}\frac{\omega_{n}^{nk}-1}{\omega_{n}^{k}-1}=0$ 。

单位根反演常用于求某个多项式特定倍数的系数和。

\begin{aligned} \sum_{i = 0}^{⌊ \frac{n}{k} ⌋} [x^{i k}] f (x) & = \sum_{i = 0}^{n} [k ∣ i] [x^{i}] f (x) \\ = \sum_{i = 0}^{n} [x^{i}] f (x) \frac{1}{k} \sum_{j = 0}^{k - 1} ω_{k}^{j i} \\ = \frac{1}{k} \sum_{i = 0}^{n} a_{i} \sum_{j = 0}^{k - 1} ω_{k}^{i j} \\ = \frac{1}{k} \sum_{i = 0}^{n} \sum_{j = 0}^{k - 1} a_{i} (ω_{k}^{j})^{i} \\ = \frac{1}{k} \sum_{j = 0}^{k - 1} f (ω_{k}^{j}) \end{aligned}

$\begin{aligned} \sum_{i=0}^{\lfloor\frac{n}{k}\rfloor}[x^{ik}]f(x)&=\sum_{i=0}^{n}[k\mid i][x^i]f(x)\\ &=\sum_{i=0}^{n}[x^i]f(x)\frac{1}{k}\sum_{j=0}^{k-1}\omega_{k}^{ji}\\ &=\frac{1}{k}\sum_{i=0}^{n}a_i\sum_{j=0}^{k-1}\omega_{k}^{ij}\\ &=\frac{1}{k}\sum_{i=0}^{n}\sum_{j=0}^{k-1}a_i(\omega_{k}^{j})^i\\ &=\frac{1}{k}\sum_{j=0}^{k-1}f(\omega_{k}^{j}) \end{aligned}$

同理可推出

\sum_{i = 0}^{⌊ \frac{n}{k} ⌋} [x^{i k + r}] f (x) = \frac{1}{k} \sum_{i = 0}^{k - 1} f (ω_{k}^{i}) ω_{k}^{- i r}

$\sum_{i=0}^{\lfloor\frac{n}{k}\rfloor}[x^{ik+r}]f(x)=\frac{1}{k}\sum_{i=0}^{k-1}f(\omega_{k}^{i})\omega_{k}^{-ir}\\$

例题

$\color{blue}{\text{[LOJ6485] LJJ } 学二项式定理}$

求 $\big{[}\sum\limits_{i=0}^{n}\big{(}\binom{n}{i}\cdot s^{i}\cdot a_{i\operatorname{mod} 4}\big{)} \big{]}\operatorname{mod} 998244353$ 。
$1\le T\le10^5,1\le n\le 10^{18},1\le s,a_0,a_1,a_2,a_3\le 10^8$ 。

\begin{aligned} \sum_{i = 0}^{n} ((\binom{n}{i}) \cdot s^{i} \cdot a_{i \mod 4}) & = \sum_{j = 0}^{3} \sum_{i = 0}^{n} (\binom{n}{i}) \cdot s^{i} \cdot a_{j} [4 | i - j] \\ = \sum_{j = 0}^{3} \sum_{i = 0}^{n} (\binom{n}{i}) \cdot s^{i} \cdot a_{j} \cdot \frac{1}{4} \sum_{k = 0}^{3} ω_{4}^{k (i - j)} \\ = \frac{1}{4} \sum_{k = 0}^{3} \sum_{j = 0}^{3} a_{j} \cdot ω_{4}^{- j k} \sum_{i = 0}^{n} (\binom{n}{i}) \cdot s^{i} \cdot ω_{4}^{i k} \\ = \frac{1}{4} \sum_{k = 0}^{3} \sum_{j = 0}^{3} a_{j} \cdot ω_{4}^{- j k} (s ω_{4}^{k} + 1)^{n} \end{aligned}

$\begin{aligned} \sum\limits_{i=0}^{n}\big{(}\binom{n}{i}\cdot s^{i}\cdot a_{i\operatorname{mod} 4}\big{)}&=\sum_{j=0}^{3}\sum_{i=0}^{n}\binom{n}{i}\cdot s^{i}\cdot a_j[4|i-j]\\ &=\sum_{j=0}^{3}\sum_{i=0}^{n}\binom{n}{i}\cdot s^i \cdot a_j\cdot\frac{1}{4}\sum_{k=0}^{3}\omega_{4}^{k(i-j)}\\ &=\frac{1}{4}\sum_{k=0}^{3}\sum_{j=0}^{3}a_j\cdot\omega_{4}^{-jk}\sum_{i=0}^{n}\binom{n}{i}\cdot s^i\cdot \omega_{4}^{ik}\\ &=\frac{1}{4}\sum_{k=0}^{3}\sum_{j=0}^{3}a_j\cdot\omega_{4}^{-jk}(s\omega_{4}^{k}+1)^n \end{aligned}$

直接计算即可，时间复杂度 $O(\log n)$ 。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int mod = 998244353;
int T, s, ans, a[5]; ll n;
inline ll qpow(ll a, ll b) {
	ll ans = 1;
	for (; b; b >>= 1, a = a * a % mod) if (b & 1) ans = ans * a % mod;
	return ans;
}
const int omega = qpow(3, (mod - 1) / 4);
const int omega_ = qpow(omega, mod - 2);
const int inv4 = qpow(4, mod - 2);
int main() {
	scanf("%d", &T);
	while (T --) {
		int ans = 0;
		scanf("%lld%d", &n, &s);
		for (int i = 0; i < 4; ++ i) scanf("%d", a + i);
		for (int k = 0; k < 4; ++ k)
			for (int j = 0; j < 4; ++ j)
				(ans += a[j] * qpow(omega_, j * k) % mod * (qpow((s * qpow(omega, k) + 1) % mod, n)) % mod) %= mod;
		printf("%d\n", (ll)ans * inv4 % mod);
	}
	return 0;
}

$[\color{blue}{\text{AtCoder Regular Contest 145 F - Modulo Sum of Increasing Sequences]}}$

对于每个 $K=0,1,\ldots,\text{MOD}-1$ ，计数长度为 $N$ ，值域为 $[0,M]$ ，总和对 $\text{MOD}$ 取模后结果为 $K$ 的不降序列 $A$ 的数量，答案对 $998244353$ 取模。

$1\le N,M\le10^6,1\le\text{MOD}\le 500$ 。

令 $b_i=a_i+i-1$ ，则有 $\forall i<n, b_i<b_{i+1}$ 。问题转化成了在 $[0,M+N-1]$ 中选取 $N$ 个数使得总和对 $\text{MOD}$ 取模后结果为 $K$ 。下令 $n=M+N$ 。为了方便，记 $p=\text{MOD}$ 。

先不考虑长度的限制，列出 OGF

F (x) = \prod_{i = 0}^{n - 1} (1 + x^{i})

$F(x)=\prod_{i=0}^{n-1}(1+x^i)$

答案即为

\frac{1}{p} \sum_{i = 0}^{p - 1} F (ω_{p}^{i}) ω_{p}^{- i k} = \frac{1}{p} \sum_{i = 0}^{p - 1} \prod_{j = 0}^{n - 1} (1 + ω_{p}^{i j}) ω_{p}^{- i k}

$\frac{1}{p}\sum_{i=0}^{p-1}F(\omega_{p}^{i})\omega_{p}^{-ik}=\frac{1}{p}\sum_{i=0}^{p-1}\prod_{j=0}^{n-1}(1+\omega_{p}^{ij})\omega_{p}^{-ik}$

$\color{blue}{\textbf{Special Case: } p\mid n}$

假设 $p\mid n$ ，记 $d=\gcd (p,i)$ ，则

\frac{1}{p} \sum_{i = 0}^{p - 1} \prod_{j = 0}^{n - 1} (1 + ω_{p}^{i j}) ω_{p}^{- i k} = \frac{1}{p} \sum_{i = 0}^{p - 1} (\prod_{j = 0}^{\frac{p}{d} - 1} (1 + ω_{\frac{p}{d}}^{i j})^{d})^{\frac{n}{p}} ω_{p}^{- i k}

$\frac{1}{p}\sum_{i=0}^{p-1}\prod_{j=0}^{n-1}(1+\omega_{p}^{ij})\omega_{p}^{-ik}=\frac{1}{p}\sum_{i=0}^{p-1}(\prod_{j=0}^{\frac{p}{d}-1}(1+\omega_{\frac{p}{d}}^{ij})^{d})^{\frac{n}{p}}\omega_{p}^{-ik}$

用分圆多项式作简化

x^{n} - 1 = \prod_{i = 0}^{n - 1} (x - ω_{n}^{i})

$x^n-1=\prod_{i=0}^{n-1}(x-\omega_{n}^{i})$

带入 $x=-1$ ，得

(- 1)^{n} - 1 = (- 1)^{n} \prod_{i = 0}^{n - 1} (1 + ω_{n}^{i})

$(-1)^{n}-1=(-1)^n\prod_{i=0}^{n-1}(1+\omega_{n}^{i})$

所以

\prod_{i = 0}^{n - 1} (1 + ω_{n}^{i}) = 1 - (- 1)^{n}

$\prod_{i=0}^{n-1}(1+\omega_{n}^{i})=1-(-1)^{n}$

那么

\begin{aligned} \frac{1}{p} \sum_{i = 0}^{p - 1} (\prod_{j = 0}^{\frac{p}{d} - 1} (1 + ω_{\frac{p}{d}}^{i j})^{d})^{\frac{n}{p}} ω_{p}^{- i k} & = \frac{1}{p} \sum_{i = 0}^{p - 1} (\prod_{j = 0}^{\frac{p}{d} - 1} (1 + ω_{\frac{p}{d}}^{i j}))^{\frac{d n}{p}} ω_{p}^{- i k} \\ = \frac{1}{p} \sum_{i = 0}^{p - 1} (1 - (- 1)^{\frac{p}{d}})^{\frac{d n}{p}} ω_{p}^{- i k} \\ = \frac{1}{p} \sum_{d ∣ p} (1 - (- 1)^{\frac{p}{d}})^{\frac{d n}{p}} (\sum_{gcd (i, p) = d} ω_{p}^{- i k}) \\ = \frac{1}{p} \sum_{d ∣ p} (1 - (- 1)^{\frac{p}{d}})^{\frac{d n}{p}} (\sum_{gcd (i, \frac{p}{d}) = 1} ω_{\frac{p}{d}}^{- i k}) \\ = \frac{1}{p} \sum_{d ∣ p} (1 - (- 1)^{\frac{p}{d}})^{\frac{d n}{p}} (\sum_{gcd (i, \frac{p}{d}) = 1} ω_{\frac{p}{d}}^{- i k}) \\ = \frac{1}{p} \sum_{d ∣ p} (1 - (- 1)^{\frac{p}{d}})^{\frac{d n}{p}} (\sum_{j ∣ \frac{p}{d}} μ (j) \sum_{i = 0}^{\frac{p}{d j} - 1} ω_{\frac{p}{d j}}^{- i k}) \\ = \frac{1}{p} \sum_{d ∣ p} (1 - (- 1)^{\frac{p}{d}})^{\frac{d n}{p}} (\sum_{j ∣ \frac{p}{d}} μ (j) [\frac{p}{d j} ∣ k] \frac{p}{d j}) \end{aligned}

$\begin{aligned} \frac{1}{p}\sum_{i=0}^{p-1}(\prod_{j=0}^{\frac{p}{d}-1}(1+\omega_{\frac{p}{d}}^{ij})^{d})^{\frac{n}{p}}\omega_{p}^{-ik}&=\frac{1}{p}\sum_{i=0}^{p-1}\big{(}\prod_{j=0}^{\frac{p}{d}-1}(1+\omega_{\frac{p}{d}}^{ij}) \big{)}^\frac{dn}{p}\omega_{p}^{-ik}\\ &=\frac{1}{p}\sum_{i=0}^{p-1}(1-(-1)^{\frac{p}{d}} )^\frac{dn}{p}\omega_{p}^{-ik}\\ &=\frac{1}{p}\sum_{d\mid p}(1-(-1)^{\frac{p}{d}})^\frac{dn}{p}(\sum_{\gcd(i,p)=d}\omega_{p}^{-ik})\\ &=\frac{1}{p}\sum_{d\mid p}(1-(-1)^{\frac{p}{d}})^\frac{dn}{p}(\sum_{\gcd(i,\frac{p}{d})=1}\omega_{\frac{p}{d}}^{-ik})\\ &=\frac{1}{p}\sum_{d\mid p}(1-(-1)^{\frac{p}{d}})^\frac{dn}{p}(\sum_{\gcd(i,\frac{p}{d})=1}\omega_{\frac{p}{d}}^{-ik})\\ &=\frac{1}{p}\sum_{d\mid p}(1-(-1)^{\frac{p}{d}})^\frac{dn}{p}(\sum_{j\mid\frac{p}{d}}\mu(j)\sum_{i=0}^{\frac{p}{dj}-1}\omega_{\frac{p}{dj}}^{-ik})\\ &=\frac{1}{p}\sum_{d\mid p}(1-(-1)^{\frac{p}{d}})^\frac{dn}{p}(\sum_{j\mid\frac{p}{d}}\mu(j)[\frac{p}{dj}\mid k]\frac{p}{dj}) \end{aligned}$

现在加上长度的限制，在初始的 OGF 中多加一个元 $y$ ，变为

\prod_{i = 0}^{n - 1} (1 + x^{i} y)

$\prod_{i=0}^{n-1}(1+x^iy)$

推导后的结果为

[y^{N}] \frac{1}{p} \sum_{d ∣ p} (1 - (- y)^{\frac{p}{d}})^{\frac{d n}{p}} (\sum_{j ∣ \frac{p}{d}} μ (j) [\frac{p}{d j} ∣ k] \frac{p}{d j})

$[y^{N}]\frac{1}{p}\sum_{d\mid p}(1-(-y)^{\frac{p}{d}})^\frac{dn}{p}(\sum_{j\mid\frac{p}{d}}\mu(j)[\frac{p}{dj}\mid k]\frac{p}{dj})$

可以用二项式定理 $O(p^2)$ 求出。

$\color{blue}{\textbf{General Case: } p\nmid n}$

根据上面的推导，可以 $O(p^2)$ 求出 $p\mid n$ 的情况，我们在此基础上做 $p\nmid n$ 的情况。

首先令 $n'=\lfloor\frac{n}{p}\rfloor\cdot p$ ， $O(p^2)$ 求出 $0\sim n'-1$ 的答案。

考虑添加进 $n'\sim n-1$ 之间的数，也就是求出在 $n'\sim n-1$ 之间选取 $i$ 个数满足这些数的和 $\operatorname{mod}{p}=j$ 的方案数，最后再和之前的答案进行合并。

因为 $n'-n<p$ ，所以这部分直接用 $dp$ 解决，时间复杂度为 $O(p^3)$ 。

至此，本题解决。总时间复杂度 $O(N+M+p^3)$ 。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 505, M = 2e6 + 5, mod = 998244353;
int n, m, p, _n, mu[N], ans[N], dp[N][N];
ll fac[M], inv[M], Inv[M];
vector<int> divisors[N];
vector<pair<int, vector<int>>> gr;
inline void precalc(int n) {
	fac[0] = inv[0] = Inv[0] = fac[1] = inv[1] = Inv[1] = 1;
	for (int i = 2; i < n; ++ i) {
		fac[i] = fac[i - 1] * i % mod;
		Inv[i] = (mod - mod / i) * Inv[mod % i] % mod;
		inv[i] = inv[i - 1] * Inv[i] % mod; 
	}
}
inline ll C(int n, int m) {
	if (n < m || n < 0 || m < 0) return 0;
	return fac[n] * inv[m] % mod * inv[n - m] % mod;
}
int main() {
	precalc(M);
	scanf("%d%d%d", &n, &m, &p), n += m, m = n - m, _n = n / p * p;
	dp[0][0] = 1;
	for (int i = 1; i <= n - _n; ++ i)
	        for (int k = i; k; -- k)
			for (int j = 0; j < p; ++ j)
				(dp[k][j] += dp[k - 1][(j - (i + _n - 1) % p + p) % p]) %= mod;
	for (int i = 1; i <= p; ++ i)
		for (int j = i; j <= p; j += i)
			divisors[j].emplace_back(i);
	mu[1] = 1;
	for (int i = 1; i <= p; ++ i)
		for (int j = i + i; j <= p; j += i)
			mu[j] -= mu[i];
	for (int i = 2; i <= p; ++ i) (mu[i] += mod) %= mod;
	for (auto d : divisors[p]) {
		vector<int> c(p);
		for (int k = 0; k < p; ++ k)
			for (auto j : divisors[p / d])
				if ((k * d * j) % p == 0)
					(c[k] += (ll)p / d / j * mu[j] % mod) %= mod;
		gr.emplace_back(d, c);
	}
	for (int use = 0; use <= min(n - _n, m); ++ use) {
		int leave = m - use;
		vector<int> coef(p);
		for (auto &&[d, c] : gr) {
			int b = p / d; ll co;
			if (leave % b > 0) continue;
			if (b & 1) co = Inv[p];
			else co = (leave / b) & 1 ? mod - Inv[p] : Inv[p];
			(co *= C(d * _n / p, leave / b)) %= mod;
			for (int i = 0; i < p; ++ i)
				(coef[i] += co * c[i] % mod) %= mod;
		}
		for (int i = 0; i < p; ++ i)
			for (int j = 0; j < p; ++ j)
				(ans[(i + j) % p] += (ll)coef[i] * dp[use][j] % mod) %= mod;
	}
	int move = (m - 1LL) * m / 2 % p;
	for (int i = 0; i < p; ++ i) printf("%d ", ans[(i + move) % p]);
	return 0;
}