[bzoj4555] [Tjoi2016&Heoi2016]求和
Description
在2016年，佳媛姐姐剛剛學習了第二類斯特林數，非常開心。
現在他想計算這樣一個函數的值:\[f(n)=\sum_{i=0}^{n}\sum_{j=0}^{i}S(i,j)*2^j*(j!)\]
S(i, j)表示第二類斯特林數，遞推公式為:
S(i, j) = j ? S(i ? 1, j) + S(i ? 1, j ? 1), 1 <= j <= i ? 1。
邊界條件為：S(i, i) = 1(0 <= i), S(i, 0) = 0(1 <= i)
你能幫幫他嗎?
Input
輸入只有一個正整數
Output
輸出f(n)。由於結果會很大，輸出f(n)對998244353(7 × 17 × 223 + 1)取模的結果即可。1 ≤ n ≤ 100000
Sample Input
3
Sample Output
87
solution
看到模數，上ntt（最近在學）
然後先推式子\[f(n)=\sum_{i=0}^n\sum_{j=0}^iS(i,j)*2^j*j!\] \[=\sum_{i=0}^n\sum_{j=0}^nS(i,j)*2^j*j!\]
枚舉上限變為n是因為i>j時S為0，方便推式子
然後把斯特林數展開，這步不懂的可以上百度，有證明\[=\sum_{i=0}^n\sum_{j=0}^n2^j\sum_{k=0}^j(-1)^k\binom{j}{k}(j-k)^i\]
再把組合數展開
\[=\sum_{i=0}^n\sum_{j=0}^n2^j\sum_{k=0}^j(-1)^k\frac{j!}{k!(j-k)!}(j-k)^i\]
改變一下枚舉順序
\[=\sum_{j=0}^n2^j*j!\sum_{k=0}^{j}\frac{(-1)^k}{k!}*\frac{\sum_{i=0}^{n}(j-k)^i}{(j-k)!}\]
定義
\[g(x)=\frac{(-1)^x}{x!},h(x)=\frac{\sum_{i=0}^nx^i}{x!}\]
然後原式變為
\[f(x)=\sum_{j=0}^n2^j*j!\sum_{k=0}^{j}g(k)*h(j-k)\]
可以發現後面一部分是一個卷積形式
然後直接上ntt
最後
\[c(j)=\sum_{k=0}^{j}g(k)*h(j-k)\]
\[f(x)=\sum_{j=0}^n2^j*j!*c(j)\]
就可以算答案了
note：\[S(n,m)=\frac{1}{m!}\sum_{k=0}^m(-1)^k\binom{m}{k}(m-k)^n\]
這裏貼一下代碼吧。。
#include<bits/stdc++.h>
using namespace std;
#pragma GCC optimize(3)
#define mod 998244353
#define int long long
#define maxn 2000050
#define read(x) scanf("%lld",&x)
#define write(x) printf("%lld\n",x)
int n,c[maxn],f[maxn],g[maxn],h[maxn],m,nn,pos[maxn],k,Inv,invn;
int qpow(int a,int x){
int res=1;a%=mod;
for(;x;x>>=1,a=1ll*a*a%mod)if(x&1)res=1ll*res*a%mod;
return res;
}
int inv(int x){return qpow(x,mod-2);}
void ntt(int *a,int op){
for(int i=1;i<n;i++) if(pos[i]>i) swap(a[pos[i]],a[i]);
for(int i=2;i<=n;i<<=1){
int wn=qpow(op==1?3:Inv,(mod-1)/i);
for(int j=0;j<n;j+=i){
int w=1;
for(int k=0;k<(i>>1);k++){
int x=a[k+j]%mod,y=a[(i>>1)+k+j]*w%mod;
a[k+j]=(x+y)%mod,a[(i>>1)+k+j]=(x-y+mod)%mod;
w=w*wn%mod;
}
}
}
if(op==1) return ;
for(int i=0;i<n;i++) a[i]=a[i]*invn%mod;
}
void solve(){
Inv=qpow(3,mod-2);
read(n);g[0]=1;int x=1;h[0]=1,h[1]=n+1;
for(int i=1;i<=n;i++) g[i]=(-g[i-1]*inv(i)%mod+mod)%mod;
for(int i=2;i<=n;i++) x=x*i%mod,h[i]=(inv(x)*(1-qpow(i,n+1))%mod*inv(1-i+mod)%mod+mod)%mod;
nn=n;m=n*2;for(n=1;n<=m;n<<=1)k++;
for(int i=1;i<=n;i++) pos[i]=(pos[i>>1]>>1)|((i&1)<<(k-1));invn=inv(n);
ntt(g,1),ntt(h,1);
for(int i=0;i<n;i++) g[i]=g[i]*h[i]%mod;
ntt(g,-1);
x=1;int q=1;int ans=1;
for(int i=1;i<=nn;i++) x=x*i%mod,q=q*2%mod,ans=(ans+x*q%mod*g[i]%mod)%mod;
write(ans);
}
#undef int
int main(){solve();}[bzoj4555] [Tjoi2016&Heoi2016]求和