[HEOI2016/TJOI2016][bzoj4555] 求和 [斯特林数+NTT]

题面

传送门

思路

首先，我们发现这个式子中大部分的项都和$j$有关（尤其是后面的$2^jast j!$），所以我们更换一下枚举方式，把这道题的枚举方式变成先$j$再$i$

$f(n)=sum_{j=0}ⁿ²jast j!sum_{i=0}^nS_ij$

第二类斯特林数有一个基于组合意义的公式：

$S_i^{j=frac1{j!}sum_{k=0}}j(-1)^kC_jk(j-k)^i=sum_{k=0}jfrac{(-1)^k(j-k)i}{k!(j-k)!}$

把这个公式代回原式中，得到：

$f(n)=sum_{j=0}ⁿ²jast j!sum_{i=0}^nsum_{k=0}jfrac{(-1)^k(j-k)i}{k!(j-k)!}$

再次更换一下枚举方式，变成：

$f(n)=sum_{j=0}ⁿ²jast j!sum_{k=0}^jfrac{(-1)k}{k!}sum_{i=0}^nfrac{(j-k)i}{(j-k)!}$

$f(n)=sum_{j=0}ⁿ²jast j!sum_{k=0}^jfrac{(-1)k}{k!}astfrac{sum_{i=0}^n(j-k)i}{(j-k)!}$

此时，设两个函数$a$和$b$，令：

$a(i)=frac{(-1)^i}{i!}$

$b(i)=frac{sum_{j=0}ⁿⁱj}{i!}=frac{i^{n+1}-1}{(i-1)i!}$

那么，

$f(n)=sum_{j=0}^n 2^jast j!ast(aast b)(j)$

其中(aast b)(j)表示$a$和$b$的$0-j$项的卷积

模数为$998244353$，用$NTT$做一遍卷积即可，时间效率为$O(nlog_2n)$

注意事项

$b(0)=1,b(1)=n+1$

这两个要提前保存一下，因为用公式推的话会div 0

还有一个奇怪的问题我没有解决，具体看代码最后面吧

Code

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
inline int read(){
    int re=0,flag=1;char ch=getchar();
    while(ch>'9'||ch<'0'){
        if(ch=='-') flag=-1;
        ch=getchar();
    }
    while(ch>='0'&&ch<='9') re=(re<<1)+(re<<3)+ch-'0',ch=getchar();
    return re*flag;
}
#define ll long long
ll MOD=998244353,g=3,inv[400010],f[400010],finv[400010];
int qpow(ll a,ll b){//快速幂
    ll re=1;
    while(b){
        if(b&1) re=re*a%MOD;
        a=a*a%MOD;b>>=1;
    }
    return re;
}
ll n,A[400010],B[400010],C[400010],r[400010],limit,cnt;
void ntt(ll *a,ll type){
    int i,j,k,mid;ll y,w,wn;
    for(i=0;i<limit;i++) if(i<r[i]) swap(a[i],a[r[i]]);
    for(mid=1;mid<limit;mid<<=1){
        wn=qpow((type==1)?g:inv[g],(MOD-1)/(mid<<1));
        for(j=0;j<limit;j+=(mid<<1)){
            w=1;
            for(k=0;k<mid;k++,w=w*wn%MOD){
                y=a[j+k+mid]*w%MOD;
                a[j+k+mid]=(a[j+k]-y+MOD)%MOD;
                a[j+k]=(a[j+k]+y)%MOD;
            }
        }
    }
    if(type==-1) for(i=0;i<limit;i++) a[i]=a[i]*inv[limit]%MOD;
}
void init(){
    limit=1;cnt=0;int i;
    while(limit<=(n<<1)) limit<<=1,cnt++;
    for(i=0;i<limit;i++) r[i]=((r[i>>1]>>1)|((i&1)<<(cnt-1)));
    inv[1]=A[0]=B[0]=f[1]=finv[1]=1;A[1]=MOD-1;B[1]=n+1;
    for(i=2;i<=limit;i++) inv[i]=(MOD-MOD/i)*inv[MOD%i]%MOD;
    for(i=2;i<=limit;i++){
        f[i]=f[i-1]*i%MOD;
        finv[i]=finv[i-1]*inv[i]%MOD;
    }
}
int main(){
    n=read();
    init();int i;
    for(i=2;i<=n;i++) A[i]=(((i%2)?-1:1)*finv[i]+MOD)%MOD;
    for(i=2;i<=n;i++) B[i]=((qpow(i,n+1)-1)*inv[i-1]%MOD*finv[i])%MOD;
    ntt(A,1);ntt(B,1);
    for(i=0;i<limit;i++) C[i]=A[i]*B[i]%MOD;
    ntt(C,-1);
    ll ans=0;
    for(i=0;i<=n;i++) ans=(ans+qpow(2,i)*f[i]%MOD*C[i]%MOD)%MOD;
    printf("%lld
",(ans+1)%MOD);//这里不知道为什么，一定要加个1，我也没有搞明白
}