LOJ6519. 魔力环(莫比乌斯反演+生成函数)
题目链接
题解
这里给出的解法基于莫比乌斯反演。可以用群论计数的相关方法代替莫比乌斯反演,但两种方法的核心部分是一样的。
环计数的常见套路就是将环视为序列。我们统计所有不同的序列,那么对于最小循环节长度为 \(d\) 的序列对应的环就会被统计 \(d\) 次。因此假设最小循环节长度为 \(x\) 的合法序列数为 \(f(x)\),那么答案即为 \(\sum_\limits{d | {\rm gcd}(n, m)} \frac{1}{d}f(d)\)。
直接求 \(f(x)\) 不好求,使用莫比乌斯反演转化为求 \(g(x)\), \(g(x)\) 表示最小循环节长度为 \(x\) 的因子的合法序列个数。那么有 \(g(x) = \sum_\limits{d | x} f(d)\),即 \(f(x) = \sum_\limits{d | x} \mu(d) g(\frac{x}{d})\)。
我们假设整个长度为 \(n\) 的序列被分成了 \(l\) 段,每一段都是原序列的一个循环节。令 \(a = \frac{n}{l}, b = \frac{m}{l}\)(\(a\) 即为每一段的珠子数,\(b\) 即为每一段的黑色珠子数量),那么我们求 \(g(a)\),等价于求下面方程的整数解的数量:
即被 \(a - b\) 颗白色珠子划分开的 \(a - b + 1\) 段黑色珠子的和为 \(b\),且满足每连续一段长度不超过 \(k\) 的限制条件。运用生成函数的知识,求上面方程的解的数量等同于求如下多项式 \(h(x)\) 中 \(x^b\) 的系数:
进一步地,有:
由于 \(\sum_\limits{i = 0}^k x^i = \frac{1 - x^{k + 1}}{1 - x}\),故:
再展开右侧的式子 \(\sum_\limits{i = 0}^k(i + 1)x^i\):
因此,我们得到了:
其中,\((1 - x^{k + 1})^{a - b - 1}\) 可化为 \(\sum_\limits{i = 0}^{\infty}\binom{a - b - 1}{i}(-1)^ix^{(k + 1)i}\),而 \(\frac{1}{(1 - x)^{a - b + 1}}\) 即 \((1 - x)^{-(a - b + 1)}\),可通过负整数次幂的二项式定理化为 \(\sum_\limits{i = 0}^{\infty}\binom{a - b + i}{i}x^i\),那么可得:
当把 \(h(x)\) 化成该形式后,要求 \(h(x)\) 中 \(x^b\) 的系数就变得非常简单了。记 \(s_1 = \sum_\limits{(k + 1)i + j = b}(-1)^i\binom{a - b - 1}{i}\binom{a - b+ j}{j}\),\(s_2 = (k + 2)\sum_\limits{(k + 1)i + j = b - k - 1}(-1)^i\binom{a - b - 1}{i}\binom{a - b+ j}{j}\),\(s_3 = (k + 1)\sum_\limits{(k + 1)i + j = b - k - 2}(-1)^i\binom{a - b - 1}{i}\binom{a - b+ j}{j}\),\(x^b\) 的系数为 \(w\),那么有 \(w = s_1 - s_2 + s_3\)。
求 \(s_1, s_2, s_3\) 只需按照 \(s_1, s_2, s_3\) 的式子枚举 \(i\) 即可,因为 \(i\) 确定 \(j\) 也就确定了。因此,我们可以在 \(\frac{b}{k + 1}\) 的时间内求出 \(h(x)\) 中 \(x^b\) 的系数。
除去反演部分,我们就能够在 \(\frac{\sigma(n)}{k + 1}\) 的时间内解决此题,其中,\(\sigma(n)\) 表示 \(n\) 的约数和。由于 \(\sigma(n)\) 可近似看作 \(n\ {\rm log}\ {\rm log}\ n\),接近线性,因此时间复杂度是非常优秀的。
代码
#include<bits/stdc++.h>
using namespace std;
#define X first
#define Y second
#define mp make_pair
#define pb push_back
#define debug(...) fprintf(stderr, __VA_ARGS__)
typedef long long ll;
typedef long double ld;
typedef unsigned int uint;
typedef pair<int, int> pii;
typedef unsigned long long ull;
template<typename T> inline void read(T& x) {
char c = getchar();
bool f = false;
for (x = 0; !isdigit(c); c = getchar()) {
if (c == '-') {
f = true;
}
}
for (; isdigit(c); c = getchar()) {
x = x * 10 + c - '0';
}
if (f) {
x = -x;
}
}
template<typename T, typename... U> inline void read(T& x, U&... y) {
read(x), read(y...);
}
template<typename T> inline bool checkMax(T& a, const T& b) {
return a < b ? a = b, true : false;
}
template<typename T> inline bool checkMin(T& a, const T& b) {
return a > b ? a = b, true : false;
}
const int N = 1e5 + 10, mod = 998244353, G = 3;
int n, m, k;
inline void add(int& x, int y) {
if (y < 0) {
y += mod;
}
x = (x + y) % mod;
}
inline void mul(int& x, int y) {
x = 1ll * x * y % mod;
}
inline int qpow(int v, int p) {
int res = 1;
for (; p; p >>= 1, mul(v, v)) {
if (p & 1) {
mul(res, v);
}
}
return res;
}
inline int gcd(int x, int y) {
return !y ? x : gcd(y, x % y);
}
int p[N], pri[N], mu[N], fac[N], invfac[N], c;
inline int binom(int n, int m) {
return n < 0 || m < 0 || n < m ? 0 : 1ll * fac[n] * invfac[m] % mod * invfac[n - m] % mod;
}
void sieve(int n) {
mu[1] = 1;
for (register int i = 2; i <= n; ++i) {
p[i] = 1;
}
for (register int i = 2; i <= n; ++i) {
if (p[i]) {
mu[i] = -1, pri[++c] = i;
}
for (register int j = 1, d; j <= c && (d = pri[j] * i) <= n; ++j) {
p[d] = 0;
if (i % pri[j] == 0) {
mu[d] = 0; break;
} else {
mu[d] = -mu[i];
}
}
}
fac[0] = invfac[0] = 1;
for (register int i = 1; i <= n; ++i) {
mul(fac[i] = fac[i - 1], i);
}
invfac[n] = qpow(fac[n], mod - 2);
for (register int i = n - 1; i; --i) {
mul(invfac[i] = invfac[i + 1], i + 1);
}
}
inline int calc(int n, int m) {
int res = 0;
for (register int i = 0; i * (k + 1) <= m; ++i) {
int j = m - i * (k + 1);
add(res, 1ll * binom(n - m - 1, i) * binom(n - m + j, j) % mod * (i & 1 ? mod - 1 : 1) % mod);
j = m - i * (k + 1) - k - 1;
if (j >= 0) {
add(res, 1ll * (k + 2) * binom(n - m - 1, i) % mod * binom(n - m + j, j) % mod * (i & 1 ? 1 : mod - 1) % mod);
}
j = m - i * (k + 1) - k - 2;
if (j >= 0) {
add(res, 1ll * (k + 1) * binom(n - m - 1, i) % mod * binom(n - m + j, j) % mod * (i & 1 ? mod - 1 : 1) % mod);
}
}
return res;
}
int f[N], g[N];
int main() {
read(n, m, k);
sieve(n);
int d = gcd(n, m);
for (register int i = 1; i <= d; ++i) {
if (d % i) {
continue;
}
g[n / i] = calc(n / i, m / i);
}
for (register int i = 1; i <= n; ++i) {
for (register int j = i; j <= n; j += i) {
add(f[j], mu[i] * g[j / i]);
}
}
int ans = 0;
for (register int i = 1; i <= n; ++i) {
if (n % i == 0) {
add(ans, 1ll * f[i] * qpow(i, mod - 2) % mod);
}
}
printf("%d\n", ans);
return 0;
}