题意:有$n$堆石子,每堆石子数量相同,以质因数分解给出,不停地从$1$到$n$依次拿石子,使得取完后石子个数为原来的因数(不能不取),当一堆只剩$1$个时结束,问在每堆石子结束的方案数
记石子个数为$\prod\limits_{i=1}^mp_i^{e_i}$,拿石头就是选一些$e_i$把它们减小
令一堆石子取$x$次取完的方案数为$f_x$,取$x-1$次还未取完的方案数为$g_x$,那么$g_x=f_x$,因为对于每个取$x-1$次未取完的方案,把它取完就对应一个$x$次取完的方案,对于每个$x$次取完的方案,它的前$x-1$次操作一定不完全相同,所以两种方案数相同
记$s=\sum e_i$,那么取$j$次取完并在$i$结束的方案数为$g_{j+1}^{i-1}f_jg_j^{n-i}=f_{j+1}^{i-1}f_j^{n-i+1}$,所以答案为$ans_i=\sum\limits_{j=1}^sf_{j+1}^{i-1}f_j^{n-i+1}$
现在我们来算$f$,直接算将$e_i$分成$x$份的方案数就是$h_x=\prod\limits_{i=1}^m\binom{e_i+x-1}{x-1}$,但$h_x\ne f_x$,因为$h_x$中可能包含那些某一次一个石子都没取的方案,所以要容斥
取$x$次中,至少有$y$次一个石子都没取的方案数为$\binom xyh_{x-y}$,那么$f_x=\sum\limits_{y=0}^x(-1)^y\binom xyh_{x-y}$,也就是说求出$h$后我们就可以$O(s\log s)$求$f$了
原题中$h$可以暴力求,但毒瘤代爷稍微加强了一下数据,注意到因为$\sum e_i=s$,所以不同的$e_i$最多会有$O(\sqrt s)$个,所以我们对相同的$e_i$一起处理,就可以在$O(s\sqrt s\log s)$的时间内递推求得$h$
总时间复杂度$O(s\sqrt s\log s)$
#include<stdio.h>
#include<algorithm>
using namespace std;
typedef long long ll;
const int mod=998244353;
int mul(int a,int b){return(ll)a*b%mod;}
int ad(int a,int b){return(a+b)%mod;}
int de(int a,int b){return(a-b)%mod;}
void inc(int&a,int b){(a+=b)%=mod;}
int pow(int a,int b){
int s=1;
while(b){
if(b&1)s=mul(s,a);
a=mul(a,a);
b>>=1;
}
return s;
}
int rev[1048576],N,iN;
void pre(int n){
int i,k=0;
for(N=1,k=0;N<=n;N<<=1)k++;
for(i=0;i<N;i++)rev[i]=(rev[i>>1]>>1)|((i&1)<<(k-1));
iN=pow(N,mod-2);
}
void ntt(int*a,int on){
int i,j,k,t,w,wn;
for(i=0;i<N;i++){
if(i<rev[i])swap(a[i],a[rev[i]]);
}
for(i=2;i<=N;i<<=1){
wn=pow(3,on==1?(mod-1)/i:mod-1-(mod-1)/i);
for(j=0;j<N;j+=i){
w=1;
for(k=0;k<i>>1;k++){
t=mul(a[i/2+j+k],w);
a[i/2+j+k]=de(a[j+k],t);
a[j+k]=ad(a[j+k],t);
w=mul(w,wn);
}
}
}
if(on==-1){
for(i=0;i<N;i++)a[i]=mul(a[i],iN);
}
}
int e[10010],inv[300010],fac[300010],rfac[300010],h[300010],f[300010],c[300010],pe[300010],pc[300010];
int C(int n,int k){return mul(fac[n],mul(rfac[k],rfac[n-k]));}
int a[1048576],b[1048576];
int main(){
int m,n,s,i,j,t,M;
scanf("%d%d",&m,&n);
s=0;
for(i=1;i<=m;i++){
scanf("%d%d",&t,e+i);
s+=e[i];
c[e[i]]++;
}
M=0;
for(i=1;i<=s;i++){
if(c[i]){
M++;
pe[M]=i;
pc[M]=c[i];
}
}
fac[0]=1;
for(i=1;i<=s;i++)fac[i]=mul(fac[i-1],i);
rfac[s]=pow(fac[s],mod-2);
for(i=s;i>0;i--)rfac[i-1]=mul(rfac[i],i);
inv[1]=1;
for(i=2;i<=s;i++)inv[i]=-mul(mod/i,inv[mod%i]);
h[1]=1;
for(i=1;i<s;i++){
h[i+1]=h[i];
for(j=1;j<=M;j++)h[i+1]=mul(h[i+1],pow(mul(pe[j]+i,inv[i]),pc[j]));
}
pre(s<<1);
for(i=0;i<=s;i++){
a[i]=(i&1?-1:1)*rfac[i];
b[i]=mul(h[i],rfac[i]);
}
ntt(a,1);
ntt(b,1);
for(i=0;i<N;i++)a[i]=mul(a[i],b[i]);
ntt(a,-1);
for(i=0;i<=s;i++)f[i]=mul(a[i],fac[i]);
for(i=1;i<=n;i++){
t=0;
for(j=1;j<=s;j++)inc(t,mul(pow(f[j+1],i-1),pow(f[j],n-i+1)));
inc(t,mod);
printf("%d\n",t);
}
}