可以先令 gi,j g i , j 表示 i i 个人里分个弹果的欢乐度, h(i) h ( i ) 表示给一个小朋友 i i 个糖果的欢乐度
则有
发现是个卷积的形式,则
gi=gi−1∗h
g
i
=
g
i
−
1
∗
h
而本题要求的是所有
g
g
的和,即
这么怎么搞呢?
理性思考,可以发现
gi+j=gi∗gj
g
i
+
j
=
g
i
∗
g
j
令
f(i)=∑ij=1gj
f
(
i
)
=
∑
j
=
1
i
g
j
则有
因为 gi+j=gi∗gj g i + j = g i ∗ g j ,则 ∑ij=i2+1gj ∑ j = i 2 + 1 i g j 和 f(i/2)∗g(i/2) f ( i / 2 ) ∗ g ( i / 2 ) 是等价的
可以理解为 f(i/2) f ( i / 2 ) 是左半边的和,要求出右半边需要将其向右平移 i/2 i / 2 个长度,而 gi+j=gi∗gj g i + j = g i ∗ g j ,那么将其乘上 g(i/2) g ( i / 2 ) 就可以了
那么式子可以变成这个样子:
这样就可以运用分治FFT来解决了,递归时同时维护 f f 和,顺便转移即可
blog参照于https://www.cnblogs.com/Skyminer/p/6561689.html
代码如下:
#include<algorithm>
#include<cstdio>
#include<cmath>
#define N 50050
using namespace std;
const double DFT=2.0,IDFT=-2.0;
const double pi=acos(-1);
struct Complex{
double x,y;
Complex(double _=0.0,double __=0.0):x(_),y(__){}
Complex operator + (Complex b) {return Complex(x+b.x,y+b.y);}
Complex operator - (Complex b) {return Complex(x-b.x,y-b.y);}
Complex operator * (Complex b) {return Complex(x*b.x-y*b.y,x*b.y+y*b.x);}
}a[N],b[N];
int m,n,len,O,S,U,MOD;
int pos[N],f[N],g[N],h[N],p[N];
inline void FFT(Complex x[],double mode){
for(int i=0;i<len;i++)
if(i<pos[i])
swap(x[i],x[pos[i]]);
for(int i=2,mid=1;i<=len;i<<=1,mid<<=1){
Complex wm(cos(2.0*pi/i),sin(mode*pi/i));
for(int j=0;j<len;j+=i){
int r=j+mid;
Complex w(1,0);
for(int k=j;k<r;k++,w=w*wm){
Complex l=x[k],r=w*x[k+mid];
x[k]=l+r;x[k+mid]=l-r;
}
}
}
if(mode==IDFT)
for(int i=0;i<len;i++)
x[i].x/=len;
return;
}
inline void Mul(int x[],int y[],int z[]){
for(int i=0;i<len;i++)
a[i]=Complex(x[i],0),b[i]=Complex(y[i],0);
FFT(a,DFT);FFT(b,DFT);
for(int i=0;i<len;i++)
a[i]=a[i]*b[i];
FFT(a,IDFT);
for(int i=0;i<=m;i++)
z[i]=int(floor(a[i].x+0.5))%MOD;
return;
}
void solve(int k){
if(k==1) return;
solve(k/2);
Mul(f,g,p);Mul(g,g,g);
for(int i=0;i<len;i++) f[i]=(f[i]+p[i])%MOD;
if(k&1){///奇数时候特判搞一搞
Mul(g,h,g);
for(int i=0;i<len;i++) f[i]=(f[i]+g[i])%MOD;
}
return;
}
int main(){
scanf("%d%d%d%d%d%d",&m,&MOD,&n,&O,&S,&U);
for(len=1;len<=m<<1;len<<=1);
for(int i=0;i<len;i++){
pos[i]=pos[i>>1]>>1;
if(i&1) pos[i]|=len>>1;
}
for(int i=1;i<=m;i++)
f[i]=g[i]=h[i]=(O*i*i+S*i+U)%MOD;
solve(n);
printf("%d",f[m]);
return 0;
}