题意
给定\(n(1\leq n\leq10^9),k(1\leq k\leq5000)\),求\(\sum_{i=1}^n{n\choose i}i^k\)对\(10^9+7\)取模
注意到如果\(i=0\),右边式子为\(0\),所以不妨把答案式子写成:\(\sum_{i=0}^n{n\choose i}i^k\)
一个很\(tricky\)的斯特林数代换:\(i^k=\sum_{i=0}^k\begin{Bmatrix}k\\j\end{Bmatrix}{i\choose j}j!\),左边是\(k\)个元素任意放在\(i\)个集合的方案数,右边是枚举有\(j\)个非空集合,把\(k\)个元素放入其中,而这\(j\)个集合之间可以任意排列,所以乘个排列数。
那么\(\sum_{i=0}^n{n\choose i}i^k=\sum_{i=0}^n{n\choose i}\sum_{j=0}^i\begin{Bmatrix}k\\j\end{Bmatrix}{i\choose j}j!\),交换求和符号得到
我们的斯特林数代换有什么好处呢?可以发现当\(j>k\)时,\(\begin{Bmatrix}k\\j\end{Bmatrix}\)的值就为\(0\)了,所以可以进一步把式子写成:
如果递推求斯特林数\(\begin{Bmatrix}n\\m\end{Bmatrix}=m\begin{Bmatrix}n-1\\m\end{Bmatrix}+\begin{Bmatrix}n-1\\m-1\end{Bmatrix}\),这个题便可以\(O(k^2)\)计算了。
但用卷积求斯特林数可以做到\(O(klog_2k)\),只不过要\(MTT\),非常毒瘤。
不过神鱼说可以推到\(O(k)\),便试了一下发现真的可以。
当\(n>k\)时,把斯特林数的展开式\( \begin{Bmatrix}n\\m\end{Bmatrix}=\frac{1}{m!} \sum_{i=0}^m(-1)^i{m\choose i}(m-i)^n\)带入最后的式子:
设\(f(i)=\sum_{j=0}^{k-i}{n-i\choose j}(-\frac{1}{2})^j\),现在考虑这个是否可以递推。由于随着\(i\)增大,求和式的上界变小,因此考虑按自变量从大到小递推。这里设\(w=(-\frac{1}{2})\)
注意到组合数的递推是\({n\choose m}={n-1\choose m}+{n-1\choose m-1}\),因此考虑用\(f(i)-f(i+1)\)。
移项得\(f(i)=(w+1)f(i+1)+{n-i-1\choose k-i}w^{k-i}\)
但是你可能会说“啊这哪里来的\(O(k)\)啊式子里不还是有个\(i^k\)吗”但你想过没有,\(i^k\)是不是可以线筛?我是没想过。但真的可以。具体就是线筛的思想,用每个\(i\)的最小质因子来求\(i^k\)。由于质数的数量远不足\(\frac{k}{log_2k}\),我们暴力算所有质数的\(k\)次方也无伤大雅。
最后,递推组合数就可以做到严格\(O(k)\)了。
#include<bits/stdc++.h>
#define rg register
#define il inline
#define cn const
#define fp(i,a,b) for(rg int i=(a),ed=(b);i<=ed;++i)
#define fb(i,a,b) for(rg int i=(a),ed=(b);i>=ed;--i)
using namespace std;
typedef cn int cint;
cint maxn=5010,mod=1e9+7;
int n,m,ans,inv[maxn],c[maxn],f[maxn];
int cnt,pri[maxn],pk[maxn];
il int fpow(int a,int b,int ans=1){
for(;b;b>>=1,a=1ll*a*a%mod)if(b&1)
ans=1ll*ans*a%mod;
return ans;
}
int main(){
cin>>n>>m;
inv[1]=1; fp(i,2,m+1)inv[i]=1ll*(mod-mod/i)*inv[mod%i]%mod;
c[0]=1; fp(i,1,m+1)c[i]=1ll*(n-i+1)*c[i-1]%mod*inv[i]%mod;
pk[1]=1; fp(i,2,m+1){
if(!pk[i])pri[++cnt]=i,pk[i]=fpow(i,m);
for(rg int j=1;j<=cnt&&pri[j]*i<=m+1;++j){
pk[i*pri[j]]=1ll*pk[i]*pk[pri[j]]%mod;
if(i%pri[j]==0)break;
}
}
if(n<=m+1)fp(i,1,n)ans=(ans+1ll*c[i]*pk[i]%mod)%mod;
else{
rg int w=mod-inv[2],pw=1,c2=1;
f[m]=1;
fb(i,m-1,1){
c2=1ll*c2*(n-i-1)%mod*inv[m-i]%mod,pw=1ll*pw*w%mod;
f[i]=(1ll*(w+1)*f[i+1]%mod+1ll*c2*pw%mod)%mod;
}
rg int p2=fpow(2,n);
fp(i,1,m){
p2=1ll*p2*inv[2]%mod;
ans=(ans+1ll*p2*c[i]%mod*pk[i]%mod*f[i]%mod)%mod;
}
}
printf("%d\n",ans);
return 0;
}