题目
题目背景
猪八戒调戏嫦娥,玉皇大帝震怒,说:“此罪当诛!”
——然后他就投胎当了猪。
题目描述
有一个
n
n
n 个珠子的手环,其中
m
m
m 个是金色的,其余是白色的。
为了产生交错的美感,要求金色的连续段长度不超过 k k k 。即,任意 k + 1 k+1 k+1 个连续的珠子,都不全是金色。
现在请问,不同的手环有多少种?补充一句,手环是可以旋转的,但是不可以翻转;珠子只有颜色的区别。形式化地说,记金色为 1 1 1,白色为 0 0 0,两个手环形成的颜色序列(从任意一个珠子开始,顺时针写下每个珠子的颜色)分别是 ⟨ a 0 , a 1 , … , a n − 1 ⟩ , ⟨ b 0 , b 1 , … , b n − 1 ⟩ \langle a_0,a_1,\dots,a_{n-1}\rangle,\;\langle b_0,b_1,\dots,b_{n-1}\rangle ⟨a0,a1,…,an−1⟩,⟨b0,b1,…,bn−1⟩,如果存在一个正整数 t t t 使得 a ( i + t ) m o d n = b i a_{(i+t)\bmod n}=b_i a(i+t)modn=bi,那么这两个手环是相同的。
输出方案数对 998244353 998244353 998244353 取模的结构。
数据范围与提示
0
≤
k
≤
m
≤
n
≤
1
0
6
0\le k\le m\le n\le 10^6
0≤k≤m≤n≤106,最多
5
5
5 组数据。
思路
显然是 b u r n s i d e \tt burnside burnside 嘛。我其实也忘了许多,但是记得证明的思路,现场重推了一下。
顺时针旋转 x x x 的操作下,不动点是 a ( i + x ) m o d n = a i a_{(i+x)\bmod n}=a_i a(i+x)modn=ai,显然就是周期 T = gcd ( x , n ) T=\gcd(x,n) T=gcd(x,n) 嘛。只需要确定前 T T T 个珠子的颜色,剩下的就都确定了。
不妨设 m < n m<n m<n,那么至少有一个白色珠子,所以 T T T 中也含有一个白色珠子。真正的手环是 T T T 复制若干次得到的,不难发现,其金色连续段要么完全在 T T T 中,要么是 T T T 首尾相接。
问题转化为,对于一个 T T T 个点的环,求金色连续段长度不超过 k k k 的方案数。不过这个环是不能旋转的(即,旋转相同视为不同方案)。
如果用
d
p
\tt dp
dp 来做,大概就是
f
(
i
,
j
)
=
∑
t
=
0
k
f
(
i
−
1
,
j
−
t
−
1
)
f(i,j)=\sum_{t=0}^{k}f(i-1,j-t-1)
f(i,j)=∑t=0kf(i−1,j−t−1) 的样子。然后你发现它很像母函数。事实上的确如此——将一个金色连续段与它右边的一个白色珠子看成一个物品,每个物品的选择数量不同,则方案不同。用
x
x
x 的指数代表 金色珠子的数量,则
(
∑
i
=
0
k
x
i
)
n
−
m
−
1
\left(\sum_{i=0}^{k}x^i\right)^{n-m-1}
(i=0∑kxi)n−m−1
为什么是 n − m − 1 n-m-1 n−m−1 次方?因为一共有 n − m n-m n−m 个白色珠子,而最后一个是首尾相接,特殊处理,所以剩下的是 n − m − 1 n-m-1 n−m−1 次方。(毕竟每个物品恰好提供一个白色珠子。)
首尾相接其实很好处理。对于长度为
i
+
1
i+1
i+1(即,有
i
i
i 的情况,可以在末尾出现
[
1
,
i
+
1
]
[1,i+1]
[1,i+1] 个珠子(剩下的就放在开头),所以方案数就是
i
+
1
i+1
i+1 了。即
G
(
x
)
=
∑
i
=
0
k
(
i
+
1
)
x
i
G(x)=\sum_{i=0}^{k}(i+1)x^i
G(x)=i=0∑k(i+1)xi
二者乘起来的第 m m m 项系数就是答案。
那么这里有一个推式子的做法,把 G ( x ) G(x) G(x) 解出来,然后 ∑ i = 0 k x i = 1 − x k + 1 1 − x \sum_{i=0}^{k}x^i=\frac{1-x^{k+1}}{1-x} ∑i=0kxi=1−x1−xk+1 的 n − m − 1 n-m-1 n−m−1 次幂就分别计算分子分母。化简式子之后,大概是 ( ( ( 三项的多项式 ) ) ) × \times × ( ( ( 二项式展开式 ) ) ) × \times × ( ( ( 二项式展开式 ) ) ) 。枚举三项多项式中的一项,再 m k = ⌊ m k ⌋ m_k=\lfloor\frac{m}{k}\rfloor mk=⌊km⌋ 枚举 1 − x k + 1 1-x^{k+1} 1−xk+1 的二项式展开式中的一项,直接就找到了 1 − x 1-x 1−x 的二项式展开式中的一项。
最神奇的是,暴力计算非常快。因为每一项要么取 1 1 − x \frac{1}{1-x} 1−x1,要么取 − x k + 1 1 − x \frac{-x^{k+1}}{1-x} 1−x−xk+1,直接枚举第二类的数量(其不超过 m k m_k mk 个),除以了 x k + 1 x^{k+1} xk+1 后全部变为 1 1 − x {1\over 1-x} 1−x1,可以组合数直接计算。——从现实意义来看,就是 容斥原理 计算不超过 k k k 个的方案数。由于 G ( x ) G(x) G(x) 只有 k k k 次,我们只需要求出 ( m − k ) (m-k) (m−k) 次到 m m m 次项,这是 O ( k ⋅ m k ) = O ( m ) \mathcal O(k\cdot m_k)=\mathcal O(m) O(k⋅mk)=O(m) 的。
然后枚举 G ( x ) G(x) G(x) 中的一项,乘上上面算出的一项。我们就已经得到了答案。
但是为何我不敢这么打?因为我不会算复杂度。每次计算是
O
(
m
)
\mathcal O(m)
O(m) 的,但是这里的
m
m
m 不是题目中直接输入的
m
m
m,而是
T
T
T 里面有多少个金色的珠子。如果枚举一个
d
=
n
/
T
d=n/T
d=n/T(即周期的个数),可知复杂度为
∑
d
∣
n
,
d
∣
m
m
d
\sum_{d\mid n,\;d\mid m}\frac{m}{d}
d∣n,d∣m∑dm
糟糕透顶的情况是
m
∣
n
m\mid n
m∣n,此时复杂度是
m
m
m 的因数和
σ
1
(
m
)
=
∏
p
i
t
i
+
1
−
1
p
i
−
1
\sigma_1(m)=\prod{p_i^{t_i+1}-1\over p_i-1}
σ1(m)=∏pi−1piti+1−1,其中
m
=
∏
p
i
t
i
m=\prod p_i^{t_i}
m=∏piti 表示质因数分解。这东西是近乎线性的(打表证明)。
代码
#include <cstdio>
#include <iostream>
#include <cstring>
#include <vector>
using namespace std;
typedef long long int_;
inline int readint(){
int a = 0; char c = getchar(), f = 1;
for(; c<'0'||c>'9'; c=getchar())
if(c == '-') f = -f;
for(; '0'<=c&&c<='9'; c=getchar())
a = (a<<3)+(a<<1)+(c^48);
return a*f;
}
const int Mod = 998244353;
inline int qkpow(int_ b,int_ q){
int_ a = 1;
for(; q; q>>=1,b=b*b%Mod)
if(q&1) a = a*b%Mod;
return a;
}
int n, m, k;
namespace MyMath{
const int MaxN = 3000005;
int jc[MaxN], inv[MaxN];
void importC(){
jc[1] = inv[1] = 1;
for(int i=2; i<MaxN; ++i){
jc[i] = 1ll*jc[i-1]*i%Mod;
inv[i] = inv[Mod%i]*1ll
*(Mod-Mod/i)%Mod;
}
for(int i=2; i<MaxN; ++i)
inv[i] = 1ll*inv[i-1]*inv[i]%Mod;
jc[0] = inv[0] = 1;
}
int getC(int n,int m){
if(m < 0 || n < m) return 0;
return 1ll*jc[n]*inv[m]%Mod*inv[n-m]%Mod;
}
int phi[MaxN];
vector< int > primes;
void importPhi(){
phi[1] = 1; // important
for(int i=2; i<MaxN; ++i)
phi[i] = i-1;
for(int i=2,len=0; i<MaxN; ++i){
if(phi[i] == i-1) ++ len,
primes.push_back(i);
for(int j=0; j<len&&primes[j]<=(MaxN-1)/i; ++j){
phi[i*primes[j]] = phi[i]*(primes[j]-1);
if(i%primes[j] == 0){
phi[i*primes[j]] = phi[i]*primes[j];
break; // till next time
}
}
}
}
int getPhi(int n){ return phi[n]; }
}
namespace HardCalc{
int getPhi(int n){
return MyMath::getPhi(n);
}
// n category, m in total
inline int getDP(int n,int m){
if(!n && !m) return 1;
int res = 0;
for(int i=0; i<=m/(k+1); ++i){
int t = MyMath::getC(n,i)*1ll
*MyMath::getC(m-i*(k+1)+n-1,n-1)%Mod;
if(i&1) res = (res+Mod-t)%Mod;
else res = (res+t)%Mod;
}
return res;
}
int query(int d){
if((1ll*m*d)%n != 0) return 0;
int res = 0, t = 1ll*m*d/n;
for(int i=0; i<=k&&i<=t; ++i)
res = (res+(i+1ll)*
getDP(d-t-1,t-i))%Mod;
return res;
}
void solve(){
int ans = 0, i;
for(i=1; i<n/i; ++i)
if(n%i == 0){
ans = (ans+1ll*getPhi(n/i)*query(i))%Mod;
ans = (ans+1ll*getPhi(i)*query(n/i))%Mod;
}
if(1ll*i*i == n) // sqrt = int
ans = (ans+1ll*getPhi(i)*query(i))%Mod;
ans = 1ll*ans*qkpow(n,Mod-2)%Mod;
printf("%d\n",ans);
}
}
int main(){
freopen("gift.in","r",stdin);
freopen("gift.out","w",stdout);
MyMath::importC();
MyMath::importPhi();
for(int T=readint(); T; --T){
n = readint(), m = readint();
k = readint();
if(!m){ puts("1"); continue; }
if(m == n){
putchar('0' + (k == m));
putchar('\n'); continue;
}
HardCalc::solve();
}
return 0;
}