Description
求
∑
i
=
1
n
i
k
⋅
r
i
\sum_{i=1}^n{i^k\cdot r^i}
i=1∑nik⋅ri
其中
n
,
r
≤
1
0
18
,
k
≤
2000
n,r\le10^{18},k\le2000
n,r≤1018,k≤2000
Solution
不会做(艹皿艹 ),怒orz题解
记
S
(
k
)
=
∑
i
=
1
n
i
k
⋅
r
i
S(k)=\sum\limits_{i=1}^{n}i^k\cdot r^i
S(k)=i=1∑nik⋅ri,那么
(
r
−
1
)
S
(
k
)
=
r
n
+
1
⋅
n
k
−
r
+
∑
i
=
0
k
−
1
r
i
[
(
i
−
1
)
k
−
i
k
]
(r-1)S(k)=r^{n+1}\cdot n^k-r+\sum\limits_{i=0}^{k-1}{r^i\left[{(i-1)}^k-i^k\right]}
(r−1)S(k)=rn+1⋅nk−r+i=0∑k−1ri[(i−1)k−ik]
我们把
(
i
−
1
)
k
{(i-1)}^k
(i−1)k二项式定理展开,那么就是
(
r
−
1
)
S
(
k
)
=
r
n
+
1
⋅
n
k
−
r
+
∑
j
=
0
k
−
1
(
k
j
)
(
−
1
)
k
−
j
∑
i
=
2
n
r
i
⋅
i
j
(r-1)S(k)=r^{n+1}\cdot n^k-r+\sum\limits_{j=0}^{k-1}{\binom{k}{j}{\left(-1\right)}^{k-j}\sum\limits_{i=2}^{n}r^i\cdot i^j}
(r−1)S(k)=rn+1⋅nk−r+j=0∑k−1(jk)(−1)k−ji=2∑nri⋅ij
然后发现最右边就是
S
(
j
)
−
r
S(j)-r
S(j)−r,于是我们就可以愉快地k2递推这个东西了。。
然后当r=1的时候上面就不成立了,这个时候求自然数幂和就行。
拉格朗日插值:
F
(
x
)
=
∑
i
=
0
k
+
1
y
i
∏
j
=
0
,
j
=
̸
i
k
+
1
x
−
x
j
x
i
−
x
j
F(x)=\sum_{i=0}^{k+1}y_i\prod_{j=0,j =\not i}^{k+1}\frac{x-x_j}{x_i-x_j}
F(x)=i=0∑k+1yij=0,j≠i∏k+1xi−xjx−xj
对于自然数幂和可以发现x都是连续的,那么就变成了下面这个样子
F
(
x
)
=
∑
i
=
0
k
+
1
y
i
∏
j
=
0
i
−
1
(
n
−
j
)
∏
j
=
i
+
1
k
+
1
(
n
−
j
)
i
!
(
k
+
1
−
i
)
!
(
−
1
)
k
+
1
−
i
F(x)=\sum_{i=0}^{k+1}y_i\frac{\prod_{j=0}^{i-1}{(n-j)}\prod_{j=i+1}^{k+1}{(n-j)}}{i!(k+1-i)!{(-1)}^{k+1-i}}
F(x)=i=0∑k+1yii!(k+1−i)!(−1)k+1−i∏j=0i−1(n−j)∏j=i+1k+1(n−j)
分子刚好就是前缀后缀积,分母刚好就是阶乘了。自然数幂和的前k项可以线性筛,阶乘和逆元可以线性求,前缀后缀积同理,那么整个做下来都是O(n)的
Code
#include <stdio.h>
#include <string.h>
#include <algorithm>
#define rep(i,st,ed) for (int i=st;i<=ed;++i)
#define drp(i,st,ed) for (int i=st;i>=ed;--i)
typedef long long LL;
const int MOD=1e9+7;
const int N=200005;
LL s[N],fac[N],inv[N],pre[N],nex[N];
LL ksm(LL x,LL dep) {
LL res=1;
x=(x%MOD+MOD)%MOD;
for (;dep;dep>>=1) {
(dep&1)?(res=res*x%MOD):0;
x=x*x%MOD;
}
return res;
}
void solve1(LL n,int k) {
LL ans=0; n%=MOD;
rep(i,1,k+1) s[i]=(s[i-1]+ksm(i,k))%MOD;
pre[0]=n; nex[k+2]=1;
rep(i,1,k+1) pre[i]=pre[i-1]*(n+MOD-i)%MOD;
drp(i,k+1,0) nex[i]=nex[i+1]*(n+MOD-i)%MOD;
rep(i,0,k+1) {
LL fz=((i>0)?(pre[i-1]):1)*(nex[i+1])%MOD;
LL fm=fac[i]*fac[k+1-i]%MOD;
LL tmp=fz*ksm(fm,MOD-2)%MOD*s[i]%MOD;
if ((k+1-i)&1) tmp=(MOD-tmp);
ans=(ans+tmp)%MOD;
}
printf("%lld\n", ans);
}
void solve2(LL n,LL k,LL r) {
s[0]=(ksm(r,n+1)+MOD-r)*ksm(r-1,MOD-2)%MOD;
rep(i,1,k) {
s[i]=(ksm(r,n+1)*ksm(n,i)%MOD+MOD-r)%MOD;
rep(j,0,i-1) {
LL tmp=fac[i]*inv[j]%MOD*inv[i-j]%MOD*(s[j]+MOD-r)%MOD;
if ((i-j)&1) tmp=(MOD-tmp)%MOD;
s[i]=((s[i]+tmp)%MOD+MOD)%MOD;
}
s[i]=s[i]*ksm(r-1,MOD-2)%MOD;
}
s[k]=(s[k]%MOD+MOD)%MOD;
printf("%lld\n", s[k]);
}
int main(void) {
freopen("data.in","r",stdin);
fac[0]=fac[1]=inv[0]=inv[1]=1;
rep(i,2,N-1) {
fac[i]=fac[i-1]*i%MOD;
inv[i]=inv[MOD%i]*(MOD-MOD/i)%MOD;
}
rep(i,2,N-1) inv[i]=inv[i-1]*inv[i]%MOD;
int T; scanf("%d",&T);
for (;T--;) {
LL n,r,k; scanf("%lld%lld%lld",&n,&k,&r);
r%=MOD;
if (r==1) solve1(n,k);
else solve2(n,k,r);
}
return 0;
}