XLkxc
Description
神犇LYD虐完HEOI之后给傻×XLk出了一题:
SHY是某国的公主,平时的一大爱好是读诗…(中间略)…结果mod p就可以了
简明题意
给定 k,a,n,d,p
f(i)=1^k+2^k+3^k+……+i^k
g(x)=f(1)+f(2)+f(3)+….+f(x)
求(g(a)+g(a+d)+g(a+2d)+……+g(a+nd))mod p
对于所有数据
1<=k<=123
0<=a,n,d<=123456789
p==1234567891
Input
第一行数据组数,(保证小于6)
以下每行四个整数 k,a,n,d
Output
每行一个结果。
Sample Input
5
1 1 1 1
1 1 1 1
1 1 1 1
1 1 1 1
1 1 1 1
Sample Output
5
5
5
5
5
这个样例啊,亦可赛艇。
思路:
可以发现,
f(x)
f
(
x
)
就是一个自然数幂和。
那么通过差分可知,
f(x)
f
(
x
)
为一个
k+1
k
+
1
次多项式。
同理根据差分可得
g(x)
g
(
x
)
为一个
k+2
k
+
2
次多项式。
令
h(x)=g(a)+g(a+d)+g(a+2∗d)+......+g(a+n∗d)
h
(
x
)
=
g
(
a
)
+
g
(
a
+
d
)
+
g
(
a
+
2
∗
d
)
+
.
.
.
.
.
.
+
g
(
a
+
n
∗
d
)
那么同样根据差分可以得到
h(x)
h
(
x
)
为一个
k+3
k
+
3
次多项式。
于是考虑拉格朗日插值。
首先预处理处
f(x)
f
(
x
)
,从而算出
g(x)
g
(
x
)
。
然后,插出
g(x)
g
(
x
)
以计算
h(x)
h
(
x
)
的前
k+4
k
+
4
项。
最后,插出
h(x)
h
(
x
)
并计算答案。
然后就没有了……
记录一下重心拉格朗日插值的公式:
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll md=1234567891;
const ll N=139;
ll k,a,n,d;
ll g[N],f[N],fac[N],inv[N],tl[N],w[N];
inline ll qpow(ll a,ll b)
{
ll ret=1;
while(b)
{
if(b&1)ret=ret*a%md;
a=a*a%md;b>>=1;
}
return ret;
}
inline void init()
{
fac[0]=1;
for(ll i=1;i<N;i++)
fac[i]=fac[i-1]*i%md;
inv[N-1]=qpow(fac[N-1],md-2);
for(ll i=N-1;i>=1;i--)
inv[i-1]=inv[i]*i%md;
}
inline ll cg(ll x,ll m,ll *y)
{
if(x<=m)return y[x];
ll l=1,sum=0;
for(ll i=1;i<=m;i++)
l=l*(x-i)%md;
for(ll i=1;i<=m;i++)
(sum+=w[i]*qpow(x-i,md-2)%md*y[i]%md)%=md;
return sum*l%md;
}
inline ll initw(ll m)
{
for(ll i=1;i<=m;i++)
w[i]=((((m-i)&1)?-1ll:1ll)*inv[m-i]*inv[i-1]%md+md)%md;
}
inline void mina()
{
scanf("%lld%lld%lld%lld",&k,&a,&n,&d);
for(ll i=0;i<=k+3;i++)
g[i]=qpow(i,k);
for(ll i=1;i<=k+3;i++)
(g[i]+=g[i-1])%=md;
for(ll i=1;i<=k+3;i++)
(g[i]+=g[i-1])%=md;
initw(k+3);
for(ll i=0;i<=k+4;i++)
f[i]=cg((a+i*d)%md,k+3,g);
for(ll i=1;i<=k+4;i++)
f[i]=(f[i]+f[i-1])%md;
initw(k+4);
printf("%lld\n",cg(n,k+4,f));
}
int main()
{
init();
int T;scanf("%d",&T);
while(T--)
mina();
return 0;
}