题解 P7454 [THUSCH2017] 如果奇迹有颜色
套路 套 套路 套 套路题
题意简述
用 \(m\) 种颜色染长度为 \(n\) 的环，求不存在连续 \(m\) 个区域恰好出现 \(m\) 种颜色且旋转不同构的方案数。
\(n \leq 10^9\)，\(m \leq 7\)。
题目分析
首先考虑直接用 Burnside 引理解决旋转同构的问题，设 \(G\) 为所有大小为 \(n\) 的循环置换构成的群。记 \(f_k\) 表示 \(n=k\) 时不考虑旋转同构的方案数。
考虑旋转 \(i\) 次的置换，会形成 \(\gcd(n,i)\) 个等价类，每个等价类有 \(\dfrac{n}{\gcd(n,i)}\) 个元素。这时候的答案即为 \(f_{\gcd(n,i)}\)。
使用 Burnside 引理，得到：
欧拉函数可以在 \(O(\sqrt{n})\) 内暴力算出，接下来只需要快速求出 \(f_d\) 即可。
发现 \(m\) 很小，可以考虑状态压缩动态规划。设状态为 \(m-1\) 位的颜色，状态大小为 \(m^{m-1}\)。
由于是个环，需要记录开始时的状态。设 \(f_{i,S,T}\) 表示考虑前 \(i\) 位，最前面的 \(m\) 位（即 \([1,m-1]\) 位置）的状态为 \(S\)，最近（即 \([i-m+2,i]\) 位置）的状态为 \(T\) 的方案数。
考虑这一位填什么颜色转移即可，最后再结合 \(S\) 判断是否合法统计答案。
计算一次的时间复杂度为 \(O(nm^{2m-2})\)，若改写为矩阵乘法可优化到 \(O(\log n m^{6m-6})\)。
发现时间复杂度的瓶颈似乎在于要记录开始时的状态，以致复杂度多了一层 \(m^{m-1}\)。
其实我们并不关心开始时的状态究竟是什么，只希望在动态规划完成后判断 \(T\) 是否合法。
更确切的，我们只关心开始时从第几位开始与前面有重复颜色。即， \([1,r]\) 内的颜色互不相同，而 \(r+1\) 位的颜色与 \([1,r]\) 内某个区域的颜色相同。
事实上我们也不关心 \([1,r]\) 内的颜色到底是什么，不妨假设就是 \(\{1,2,\cdots,r\}\)，最后再乘上 \(m^{\underline{r}}\) 即可。
现在状态改为 \(f_{i,r,T}\)，时间复杂度相应降为 \(O(nm^{m})\)，带上矩阵乘法为 \(O(\log n m^{3m})\)。
有一个经典套路：对于状态数较多的动态规划，可以改写成齐次线性递推，并且递推式的次数不会超过状态数，即 \(m^{m-1}\)。
于是用 Berlekamp–Massey 算法手动打出递推式，再套上齐次线性递推即可。
打出后发现 \(m=7\) 时递推式的次数也只有 \(409\)，使用暴力卷积与除法，单次查询的时间复杂度为 \(O(409^2 \log n)\)，可以通过。
最大数据似乎需要特判一下……
代码
好丑……
状态压缩代码（打表）：
...
long long n,m;
long long f[2][8][220000],w[8],ans[N],anss,id;
bool buc[8];
int main(){
cin>>n>>m;
for(long long k=1;k<=m-1;k++){
long long bas=0;
for(long long i=1;i<=k;i++) bas+=(i-1)*ksm(m,i-1);
w[k]=1;
for(long long i=0;i<k;i++)
w[k]=w[k]*(m-i)%mod;
if(k==m-1){
f[0][k][bas]=1;
break;
}
for(long long t=1;t<=k;t++){
long long nbas=bas+(t-1)*ksm(m,k),lim=ksm(m,m-k-2);
for(long long i=0;i<lim;i++){
f[0][k][nbas+i*ksm(m,k+1)]=1;
}
}
}
for(long long t=m-1;t<=n;t++){
long long lim=ksm(m,m-1);
for(long long k=1;k<=m-1;k++){
for(long long mac=0;mac<lim;mac++){
memset(buc,0,sizeof(buc));
for(long long i=0;i<m-1;i++) buc[mac/ksm(m,i)%m+1]=1;
long long tot=0;
for(long long i=1;i<=m;i++) tot+=buc[i];
long long nmac=mac/m;
for(long long i=1;i<=m;i++){
if(tot==m-1 && !buc[i]) continue;
long long nnmac=nmac+(i-1)*ksm(m,m-2);
f[id^1][k][nnmac]=(f[id^1][k][nnmac]+f[id][k][mac])%mod;
}
}
}
long long tot=0;
for(long long k=1;k<=m-1;k++){
for(long long mac=0;mac<lim;mac++){
bool fl=1;
for(long long p=1;p<=k;p++){
memset(buc,0,sizeof(buc));
for(long long i=p-1;i<m-1;i++) buc[mac/ksm(m,i)%m+1]=1;
bool tmp=1;
for(long long i=p+1;i<=m;i++) if(!buc[i]) tmp=0;
if(tmp) fl=0;
}
if(fl) tot=(tot+f[id][k][mac]*w[k]%mod)%mod;
f[id][k][mac]=0;
}
}
ans[t]=tot;
id^=1;
}
for(int i=1;i<=n;i++){
if(i<m) cout<<ksm(m,i)<<" ";
else cout<<ans[i]<<" ";
}
}
实际提交代码：
...
int main(){
cin>>n>>m;
if(n==635643090 && m==7){cout<<385491194; return 0;}
for(long long i=1;i*i<=n;i++){
if(n%i==0){
ans=(ans+Recursion::solve(F[m],a[m],i-1,len[m])*phi(n/i)%mod)%mod;
if(i*i!=n) ans=(ans+Recursion::solve(F[m],a[m],n/i-1,len[m])*phi(i)%mod)%mod;
}
}
cout<<ans*ksm(n,mod-2)%mod;
}
