bzoj4555(数学推导+画柿子+NTT）

题面
今天学了多项式逆元，找到了这题。但是本蒟蒻觉得这题和逆元并没有什么关系，和逆元有关的题解也看不懂。。。
题目要你求这样一个函数，结果对一个费马素数取模。

S(i,j)为第二类斯特林数，题目还恶意地给出了递推公式。。。

经过一顿百度，知道了第二类斯特林数是这样的

含义是将n个不同的元素拆分成m个集合的方案数，显然当m>n时，S(n,m)=0。

所以原式首先可以写成

ans=∑i=0n∑j=0n∗2j∗(j!)∗S(i,j)

然后代第二类斯特林数的公式

ans=∑i=0n∑j=0n∗2j∗∑k=0j(−1)k∗C(j,k)∗(j−k)i

接着把组合数也代进去

ans=∑i=0n∑j=0n∗2j∗∑k=0j(−1)k∗j!k!∗(j−k)!∗(j−k)i

ans=∑i=0n∑j=0n∗2j∗(j!)∑k=0j(−1)kk!∗(j−k)i(j−k)!

这个最内层的求和已经是一个卷积了,但是最外层还有个碍眼的i，难道要做n次卷积吗？

整条式子只有(j−k)i里面有i，根据乘法分配律或者画柿子的套路，把i的枚举扔到最内层就可以了。

ans=∑j=0n∗2j∗(j!)∑k=0j(−1)kk!∗∑ni=0(j−k)i(j−k)!

这样看起来就很舒服，预处理阶乘和阶乘逆元，等比数列求和，然后上NTT就可以了。（998244353的原根是3）。

#include <iostream>#include <fstream>#include <algorithm>#include <cmath>#include <ctime>#include <cstdio>#include <cstdlib>#include <cstring>using namespace std;#define mmst(a, b) memset(a, b, sizeof(a))#define mmcp(a, b) memcpy(a, b, sizeof(b))typedef long long LL;const int N=(1<<18)+5;const LL p=998244353,g=3;LL cheng(LL a,LL b){    LL res=1ll;    for(;b;b>>=1,a=a*a%p)    if(b&1)    res=res*a%p;    return res;}int n,rev[N];void init(int lim){    int k=-1;    n=1;    while(n<lim)    n<<=1,k++;    for(int i=0;i<n;i++)    rev[i]=(rev[i>>1]>>1) | ((i&1)<<k);}void ntt(LL *a,int ops){    for(int i=0;i<n;i++)    if(i<rev[i])    swap(a[i],a[rev[i]]);    for(int l=2;l<=n;l<<=1)    {        int m=l>>1;        LL wn;        if(ops)        wn=cheng(g,(p-1)/l);        else        wn=cheng(g,p-1-(p-1)/l);        for(int i=0;i<n;i+=l)        {            LL w=1;            for(int k=0;k<m;k++)            {                LL t=a[i+k+m]*w%p;                a[i+k+m]=(a[i+k]-t+p)%p;                a[i+k]=(a[i+k]+t)%p;                w=w*wn%p;            }        }    }    if(!ops)    {        LL Inv=cheng(n,p-2);        for(int i=0;i<n;i++)        a[i]=a[i]*Inv%p;    }}LL nn,ans,er[N],jc[N],a[N],b[N];int main(){    cin>>nn;    jc[0]=1ll;    for(LL i=1;i<=nn;i++)    jc[i]=jc[i-1]*i%p;    er[0]=1ll;    for(int i=1;i<=nn;i++)    er[i]=er[i-1]*2%p;    for(LL i=0;i<=nn;i++)    a[i]=cheng(jc[i],p-2);    b[0]=1;    b[1]=nn+1;    for(LL i=2;i<=nn;i++)    b[i]=(cheng(i,nn+1)-1+p)%p*cheng(i-1,p-2)%p*a[i]%p;    for(int i=0;i<=nn;i++)    if(i&1)    a[i]=p-a[i];    init(2*nn+5);    ntt(a,1);    ntt(b,1);    for(int i=0;i<n;i++)    a[i]=a[i]*b[i]%p;    ntt(a,0);    for(int i=0;i<=nn;i++)    ans=(ans+jc[i]*er[i]%p*a[i]%p)%p;    cout<<ans<<endl;    return 0;}