[BZOJ4555][Tjoi2016&Heoi2016]求和(FFT)

来源:互联网 发布:linux 分配权限 编辑:程序博客网 时间:2024/06/05 08:31

=== ===

这里放传送门

=== ===

题解

首先值得吐槽的就是这个题如果不知道第二类斯特林数的通项公式的话是很难做的。。想想万一考试的时候遇到这个题然后不知道通项公式两眼抓瞎不是很GG。。而且题目还偏偏给了你递推公式。。明显就是在说【有本事你自己推嘛】之类的嘛。。。= =

然后ATP看到这个题第一反应就是去查了通项公式。。

S(n,m)=1m!k=1m(1)kCkm(mk)n

没错就决定是它了。
题目要求的式子是
i=1nj=1iS(i,j)×2j×j!

第二重枚举的上限跟i有关让人看起来很不爽,我们可以把枚举上界强行扩大到n。因为斯特林数一个直观意义就是把n个不同的小球放到m个相同的盒子里并且不允许有空的方案数,当n小于m的时候当然就是0。所以多枚举的那些对最后的答案不会有影响。

然后把通项公式代进去:

i=1nj=1n1j!k=1j(1)kCjk(jk)i×2j×j!

显然j!1j!可以约掉了,那我们就把它约掉。然后怎么办呢?似乎看起来有什么(1)k(jk)i让我们能够联想到卷积之类的东西,但那个组合数好像很碍事,那我们把它化成阶乘的形式。
然后式子变成了:
i=1nj=1nk=1j(1)kj!k!(jk)!×(jk)i×2j

那个2jj!好像跟k无关诶那我们把它扔到第三层枚举外面去吧,然后把和k有关的量和跟jk有关的量写到一起去,把式子变成:

i=1nj=1n2j×j!k=1j(1)kk!×(jk)i(jk)!

这是卷积吗?好像不大是。。最大的不和谐因素就是后面那个(jk)i,这个东西跟i有关,但我们又不能做n遍卷积。。。可是当我们固定一个j的时候,i=1..n的所有值都被算了一遍。那么我们可以把i那个东西拿到里面去让它先把所有的都加起来。然后式子就变成了:
j=1n2j×j!k=1j(1)kk!×ni=1(jk)i(jk)!

那么我们设f(i)=(1)ii!g(i)=k=1niki!h(x)=f(x)×g(x),那么式子就变成了

j=1n2j×j!×h(j)

一个卷积就可以搞出来啦!

代码

#include<cstdio>#include<cstring>#include<algorithm>using namespace std;const long long Mod=998244353;const long long G=3;int n,L,R[300010],N,M;long long a[300010],b[300010],mul[300010],ans,inv;long long powww(long long a,int t){    long long ans=1;a%=Mod;    while (t!=0){        if (t&1) ans=(ans*a)%Mod;        a=(a*a)%Mod;t>>=1;    }    return ans;}void NTT(long long *a,int N,int opt){    long long wn,w,x,y;    for (int i=0;i<N;i++)      if (i<R[i]) swap(a[i],a[R[i]]);    for (int k=1;k<N;k<<=1){        wn=powww(G,(Mod-1)/(k<<1));        for (int p=(k<<1),i=0;i<N;i+=p){            w=1;            for (int j=0;j<k;j++){                x=a[i+j];y=a[i+j+k]*w%Mod;                a[i+j]=(x+y)%Mod;                a[i+j+k]=(x-y+Mod)%Mod;                w=(w*wn)%Mod;            }        }    }    if (opt==-1) reverse(a+1,a+N);}int main(){    scanf("%d",&n);mul[0]=1;    for (int i=1;i<=n;i++)      mul[i]=mul[i-1]*i%Mod;    M=2*n;    for (N=1;N<=M;N<<=1) L++;    for (int i=0;i<=N;i++)      R[i]=(R[i>>1]>>1)|((i&1)<<(L-1));    for (int i=0;i<=n;i++){        a[i]=powww(mul[i],Mod-2);        if  (i&1) a[i]=Mod-a[i];        b[i]=(1-powww(i,n+1))%Mod;        b[i]=b[i]*powww(1-i,Mod-2)%Mod;        b[i]=b[i]*powww(mul[i],Mod-2)%Mod;        b[i]=(b[i]+Mod)%Mod;    }    b[1]=n+1;    NTT(a,N,1);NTT(b,N,1);    for (int i=0;i<=N;i++) a[i]=(a[i]*b[i])%Mod;    NTT(a,N,-1);    inv=powww(N,Mod-2)%Mod;    for (int i=0;i<=n;i++){        long long tmp=powww(2,i);        tmp=tmp*mul[i]%Mod;        tmp=tmp*a[i]%Mod*inv%Mod;        ans=(ans+tmp)%Mod;    }    ans=(ans+Mod)%Mod;    printf("%I64d\n",ans);    return 0;}
0 0