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

題面

傳送門

思路

首先，咱們發現這個式子中大部分的項都和$j$有關（尤爲是後面的$2^j\ast j!$），因此咱們更換一下枚舉方式，把這道題的枚舉方式變成先$j$再$i$

$f(n)=\sum_{j=0}^n2^j\ast j!\sum_{i=0}^nS_i^j$

第二類斯特林數有一個基於組合意義的公式：

$S_i^j=\frac1{j!}\sum_{k=0}^j(-1)^kC_j^k(j-k)^i=\sum_{k=0}^j\frac{(-1)^k(j-k)^i}{k!(j-k)!}$

把這個公式代回原式中，獲得：

$f(n)=\sum_{j=0}^n2^j\ast j!\sum_{i=0}^n\sum_{k=0}^j\frac{(-1)^k(j-k)^i}{k!(j-k)!}$

再次更換一下枚舉方式，變成：

$f(n)=\sum_{j=0}^n2^j\ast j!\sum_{k=0}^j\frac{(-1)^k}{k!}\sum_{i=0}^n\frac{(j-k)^i}{(j-k)!}$

$f(n)=\sum_{j=0}^n2^j\ast j!\sum_{k=0}^j\frac{(-1)^k}{k!}\ast\frac{\sum_{i=0}^n(j-k)^i}{(j-k)!}$

此時，設兩個函數$a$和$b$，令：

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

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

那麼，

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

其中(a\ast 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\n",(ans+1)%MOD);//這裏不知道爲何，必定要加個1，我也沒有搞明白
}