模拟赛 斯特林数计数 题解（拉格朗日插值）
题目大意：
求
答案对\(10^9 + 7\)取模。
首先，我们发现
这个东西出现了很多次。
设R=。
R可以矩阵乘法求。
根据斯特林数的
原式可化为
设所求为\(f_k(N)\)。
对于\(R=1\)的情况,自然数幂求和即可。\(O(k+logn)\)。
所以，
由上述暴力展开的情况可以发现，存在一个关于\(N\)的\(m\)次多项式\(F_m(N)\)，
使得\(f_m(N)=F_m(N)R^{N+1}-F_m(0)R\)，利用上述的展开方式也可归纳证明这个结论。
所以，求出\(F_m(N)\)即可，由于它是\(m\)次多项式，所以，求出\(F_m(0)到F_m(m)\)即可插值求出\(F_m(N)\)。
因为\(f_m(N)-f_m(N-1)=F_m(N)R^{N+1}-F_m(N-1)R^N=N^m R^N\)，所以\(F_m(N)R-F_m(N-1)=N^m\)，
即\(F_m(N)=\frac{N^m+F_m(N-1)}{R}\)。
设\(F_m(0)=x\)，那么\(F_m(i)\)就能表示为\(k_ix+b_i\)的形式。
那么，我们只要再找一条\(F_m\)的等量关系即可求出。
代入到上一条式子，即可得出
（其实就是说m次多项式的m+1次差分等于0，这条式子可以保证\(F_m\)形成一个m次多项式）。
也可以这样理解：我们用\(m+1\)个点(0~m)就能确定一个多项式，但要保证确定的多项式代入\(m+1\)的值是正确的。
通过拉格朗日插值也能得出这条式子。
这样，我们就能求出\(F_m(0)到F_m(m)\)，再通过插值，即可求出\(F_m(N)\)，进而求出答案。
通过线性筛,可以在计算质数\(p\)的\(p^m\)后线性地算出所有正整数\(x\)的\(x^m\)。
时间复杂度为\(O(\frac{m}{ln(m)}*log_2(m))=O(m)\)。
总时间复杂度\(O(k+logn)\)。
代码：
#include <stdio.h>
#define md 1000000007
#define ll long long
int ny[200010],jc[200010],jn[200010];
int C(int n,int m)
{
return 1ll*jc[n]*jn[m]%md*jn[n-m]%md;
}
int ksm(int a,int b)
{
int jg=1;
while(b>0)
{
if(b&1)
jg=1ll*jg*a%md;
a=1ll*a*a%md;
b=(b>>1);
}
return jg;
}
int getfib(ll x)
{
int a=1,b=0,c=0,d=1,w[70],s=0;
while(x>0)
{
w[s++]=x%2;
x/=2;
}
for(int i=s-1;i>=0;i--)
{
int bc=1ll*b*c%md,ad=(a+d)%md;
b=1ll*b*ad%md;c=1ll*c*ad%md;
a=(1ll*a*a+bc)%md;
d=(1ll*d*d+bc)%md;
if(w[i])
{
int ta=a,tc=c;
a=(a+b)%md;
c=(c+d)%md;
b=ta;d=tc;
}
}
return c;
}
int qz[200010],hz[200010];
int jisuan(int sz[200010],int n,ll x)
{
for(int i=0;i<=n;i++)
{
if(i==0)qz[i]=1;
else qz[i]=qz[i-1];
qz[i]=1ll*qz[i]*(x%md-i+md)%md;
}
for(int i=n;i>=0;i--)
{
if(i==n)hz[i]=1;
else hz[i]=hz[i+1];
hz[i]=1ll*hz[i]*(x%md-i+md)%md;
}
int ans=0;
for(int i=0;i<=n;i++)
{
int l=1ll*jn[i]*jn[n-i]%md;
if(i>0)l=1ll*l*qz[i-1]%md;
if(i<n)l=1ll*l*hz[i+1]%md;
if((n-i)%2==1)
l=1ll*l*(md-1)%md;
ans=(ans+1ll*l*sz[i])%md;
}
return ans;
}
int mi[200010],k[200010],b[200010],S[200010],F[200010];
int sa[200010],ss[200010],sl=0;
void yucl(int x,int m)
{
mi[1]=1;
for(int i=2;i<=x;i++)
{
if(!sa[i])
{
ss[sl++]=i;
mi[i]=ksm(i,m);
}
for(int j=0;j<sl;j++)
{
if(1ll*i*ss[j]>x)
break;
sa[i*ss[j]]=true;
mi[i*ss[j]]=1ll*mi[i]*mi[ss[j]]%md;
if(i%ss[j]==0)break;
}
}
}
int main()
{
freopen("c.in","r",stdin);
freopen("c.out","w",stdout);
ll n,r;int K;
scanf("%lld%lld%d",&n,&r,&K);
ny[1]=jc[0]=jn[0]=1;
for(int i=2;i<=K+1;i++)
ny[i]=1ll*(md-md/i)*ny[md%i]%md;
for(int i=1;i<=K+1;i++)
{
jc[i]=1ll*jc[i-1]*i%md;
jn[i]=1ll*jn[i-1]*ny[i]%md;
}
yucl(K+1,K);
int R=(getfib(r+2)-1+md)%md;
if(R==1)
{
for(int i=1;i<=K+1;i++)
mi[i]=(mi[i]+mi[i-1])%md;
printf("%d",jisuan(mi,K+1,n));
}
else
{
int rn=ksm(R,md-2);
k[0]=1;b[0]=0;
for(int i=1;i<=K+1;i++)
{
k[i]=1ll*k[i-1]*rn%md;
b[i]=1ll*(b[i-1]+mi[i])*rn%md;
}
int hk=0,hb=0;
for(int i=0;i<=K+1;i++)
{
int t=C(K+1,i);
if(i%2==1)
t=1ll*t*(md-1)%md;
hk=(hk+1ll*t*k[i])%md;
hb=(hb+1ll*t*b[i])%md;
}
F[0]=1ll*(md-hb)*ksm(hk,md-2)%md;
for(int i=1;i<=K;i++)
F[i]=(1ll*k[i]*F[0]+b[i])%md;
int fn=jisuan(F,K,n);
int ans=1ll*fn*ksm(R,(n+1)%(md-1))%md;
ans=(ans-1ll*F[0]*R%md+md)%md;
printf("%d\n",ans);
}
fclose(stdin);
fclose(stdout);
return 0;
}
