自动搬运

查看原文

来自洛谷，原作者为

	Kubic Go straight ahead 'til we've lost it all.

搬运于2025-08-24 22:55:51，当前版本为作者最后更新于2024-03-03 13:02:51，作者可能在搬运后再次修改，您可在原文处查看最新版

自动搬运只会搬运当前题目点赞数最高的题解，您可前往洛谷题解查看更多

以下是正文

---

一种不需要高级线性代数知识的低复杂度做法。

我们只考虑 $m=0$ 的情况。实际上，通过之后的分析可以知道 $m\neq 0$ 的情况可以归约为若干个 $\sum n$ 不超过原来的 $n$ 且 $m=0$ 的子问题。

依次加入每种颜色的边。加入第一种颜色时，我们得到了若干个环。

令点 $u$ 所在的环大小为 $size(u)$。加入第二种颜色时，对于一条边 $u\rightarrow v$，此时要求 $size(u)=size(v)$，否则与条件四矛盾。

进一步利用条件四，我们可以刻画第二种颜色的边：将环划分成若干个集合，要求每个集合中所有环大小相等。对于一个集合，令其中的环为 $C_1\dots C_p$，大小均为 $L$，令 $C_{i,j}$ 表示环 $C_i$ 中（从任意点开始）顺时针编号的第 $j$ 个点。则 $\forall i\in [1,p),j\in [0,L)$，存在第二种颜色的边 $C_{i,j}\rightarrow C_{i+1,j}$。特殊地，$\exist d\in [0,L),\forall j\in [0,L)$，存在边 $C_{p,j}\rightarrow C_{1,(j+d)\bmod L}$。

相当于将一个集合中的环合并成了一个“柱”。

再加入第三种颜色。通过一些分析可以得到：将“柱”划分为若干个集合，要求每个集合中所有“柱”同构。连边方式也与之前类似，即为每个“柱”向下一个“柱”的某一种同构方式的对应点连边。

以此类推，我们可以猜测出 $k$ 种颜色的情况：之前的“环”与“柱”都可以被推广为“连通块”，每一轮都是对等价的“连通块”进行合并，形成更大的“连通块”。

实际上，我们可以通过归纳得到一个较为关键的性质：对于一个“连通块”，若看作无标号，那么其中任意两点都是等价的，即存在一种自同构使得这两个点对应，并且这种自同构是唯一的。

因此我们令 $f_{i,j}$ 表示进行 $j$ 步操作能够合并出多少种不同的大小为 $i$ 的无标号“连通块”。有转移：

$$f_{i,j}\times i\rightarrow f_{i',j+1}$$

其中要求 $i\mid i'$。

其意义为这一步合并了 $\dfrac{i'}{i}$ 个两两同构的大小为 $i$ 的“连通块”。因为无标号，所以我们只乘一个 $i$ 表示最后一层到第一层的连边方式。

我们只需要对于每个 $i\in [1,n]$ 计算 $f_{i,k}$，再进行多项式 exp 即可得到答案。

但 $j$ 这一维可能非常巨大，需要使用一些手段对其进行优化。

将 $f$ 写作生成函数，即令 $F_i(x)=\sum f_{i,j}x^j$。转移可以改写为：

$$F_i=\dfrac{1}{1-ix}\sum\limits_{j\mid i}F_j\times jx $$

进一步将 $F_i$ 写作有理分式，容易发现可以令 $F_i$ 分母为 $\prod\limits_{j\mid i}{(1-jx)}$。

只需维护 $F_i$ 的分子（可以证明次数不超过分母），最后用线性递推求出 $k$ 次项系数，即可将时间复杂度降为 $O(\text{poly}(n,\log k))$，使用较好的实现方式可通过本题。

更进一步地，注意到分母形式的特殊性，我们可以将 $F_i$ 写作分式分解形式，即 $F_i=\sum\limits_{j\mid i}\dfrac{w_{i,j}}{1-jx}$。

按照上述转移式直接暴力维护 $w_{i,j}$，即可做到 $O\left(\sum\dfrac{n}{i}d(i)\right)=O(n\log^2 n)$。

这个转移还可以通过类似于“在线高维前缀和”的方式进一步优化，令 $s_{i,j}$ 表示在 $i$ 处只对前 $j$ 个质因子做高维前缀和得到的结果。可以将时间复杂度降为 $O(n\log n\log\log n)$。

参考代码（$O(n\log^2 n)$）：

#include <bits/stdc++.h>
using namespace std;
#define ll long long
#define pb push_back
const int N=2e6+5,MOD=998244353;
int n,pw[N],id[N],z[N],z1[N];ll m;vector<int> d[N],f[N];
int l,lim,lim1,invL,r[N],inv[N],tmp1[N],tmp2[N],tmp3[N],g[2][N];
void W(int &x,int y) {x+=y;if(x>=MOD) x-=MOD;}
void W(int &x,ll y) {x=(x+y)%MOD;}
int add(int x,int y) {x+=y;return x<MOD?x:x-MOD;}
int qPow(int x,int y)
{int res=1;for(;y;y>>=1,x=1ll*x*x%MOD) if(y&1) res=1ll*res*x%MOD;return res;}
void init(int n)
{
	l=0;lim=1;while(lim<n) ++l,lim*=2;invL=qPow(lim,MOD-2);
	for(int i=0;i<lim;++i) r[i]=(r[i>>1]>>1)|((i&1)<<l-1);
	if(lim>lim1)
	{
		for(int i=1;i<lim;++i) inv[i]=i>1?1ll*inv[MOD%i]*(MOD-MOD/i)%MOD:1;
		for(int i=1,t1,t2,t3,t4;i<lim;i*=2)
		{
			t1=qPow(3,(MOD-1)/(i*2));t2=qPow(t1,MOD-2);t3=t4=1;
			for(int j=0;j<i;++j,t3=1ll*t3*t1%MOD,t4=1ll*t4*t2%MOD)
				g[0][i+j]=t3,g[1][i+j]=t4;
		}lim1=lim;
	}
}
void deriv(int n,int a[]) {for(int i=1;i<n;++i) a[i-1]=1ll*a[i]*i%MOD;a[n-1]=0;}
void integ(int n,int a[]) {for(int i=n-1;i;--i) a[i]=1ll*a[i-1]*inv[i]%MOD;a[0]=0;}
void NTT(bool fl,int a[])
{
	for(int i=0;i<lim;++i) if(i<r[i]) swap(a[i],a[r[i]]);
	for(int i=1,t1,t2;i<lim;i*=2) for(int j=0;j<lim;j+=i*2) for(int k=0;k<i;++k)
	{
		t1=a[j+k];t2=1ll*g[fl][i+k]*a[i+j+k]%MOD;
		a[j+k]=add(t1,t2);a[i+j+k]=add(t1,MOD-t2);
	}if(fl) for(int i=0;i<lim;++i) a[i]=1ll*a[i]*invL%MOD; 
}
void polyInv(int n,int a[],int res[])
{
	if(n==1) {res[0]=qPow(a[0],MOD-2);return;}polyInv((n+1)/2,a,res);
	for(int i=0;i<n;++i) tmp1[i]=a[i];for(int i=n;i<lim;++i) tmp1[i]=0;
	init(n*2-1);NTT(0,tmp1);NTT(0,res);
	for(int i=0;i<lim;++i) res[i]=(2-1ll*tmp1[i]*res[i]%MOD+MOD)*res[i]%MOD;
	NTT(1,res);for(int i=n;i<lim;++i) res[i]=0;
}
void polyLn(int n,int a[])
{
	init(n*2-1);for(int i=0;i<lim;++i) tmp2[i]=0;
	polyInv(n,a,tmp2);deriv(n,a);NTT(0,a);NTT(0,tmp2);
	for(int i=0;i<lim;++i) a[i]=1ll*a[i]*tmp2[i]%MOD;
	NTT(1,a);integ(n,a);for(int i=n;i<lim;++i) a[i]=0;
}
void polyExp(int n,int a[],int res[])
{
	if(n==1) {res[0]=1;return;}polyExp((n+1)/2,a,res);
	for(int i=0;i<n;++i) tmp3[i]=res[i];for(int i=n;i<lim;++i) tmp3[i]=0;
	polyLn(n,tmp3);for(int i=0;i<n;++i) tmp3[i]=add(a[i],MOD-tmp3[i]);++tmp3[0];
	NTT(0,tmp3);NTT(0,res);for(int i=0;i<lim;++i) res[i]=1ll*res[i]*tmp3[i]%MOD;
	NTT(1,res);for(int i=n;i<lim;++i) res[i]=0;
}
int main()
{
	scanf("%*d %d %*d %lld",&n,&m);m%=MOD-1;
	for(int i=1;i<=n;++i)
	{
		pw[i]=qPow(i,m);inv[i]=i>1?1ll*inv[MOD%i]*(MOD-MOD/i)%MOD:1;
		for(int j=i;j<=n;j+=i) d[j].pb(i);
	}
	for(int i=1,t;i<=n;++i)
	{
		f[i].resize(d[i].size());if(i==1) {f[i][0]=z[1]=1;continue;}
		for(int j=0;j<d[i].size();++j) id[d[i][j]]=j;
		for(auto j:d[i]) if(j<i) for(int k=0;k<d[j].size();++k)
			W(f[i][id[d[j][k]]],f[j][k]);
		for(int j=0;j<d[i].size();++j)
		{
			if(d[i][j]<i)
			{
				t=1ll*f[i][j]*(MOD-inv[i-d[i][j]])%MOD;
				f[i][j]=t;W(f[i].back(),MOD-t);
			}W(z[i],1ll*f[i][j]*pw[d[i][j]]);
		}z[i]=1ll*z[i]*inv[i]%MOD;for(auto &x:f[i]) x=1ll*x*i%MOD;
	}polyExp(n+1,z,z1);for(int i=1;i<=n;++i) z1[n]=1ll*z1[n]*i%MOD;
	printf("%d\n",z1[n]);return 0;
}