题意

给出一个数组A&＃xff0c;经过一次处理&＃xff0c;生成一个数组S&＃xff0c;数组S中的每个值相当于数组A的累加&＃xff0c;比如&＃xff1a;A &＃61; {1 3 5 6} &＃61;> S &＃61; {1 4 9 15}。如果对生成的数组S再进行一次累加操作&＃xff0c;{1 4 9 15} &＃61;> {1 5 14 29}&＃xff0c;现在给出数组A&＃xff0c;问进行K次操作后的结果。(输出结果 Mod 10^9 &＃43; 7)
2 <&＃61; n <&＃61; 50000, 0 <&＃61; k <&＃61; 10^9, 0 <&＃61; a[i] <&＃61; 10^9

分析

首先有个结论就是ans[n]&＃61;∑ni&＃61;0Cii&＃43;k−1∗a[n−i]ans[n]&＃61;∑i&＃61;0nCi&＃43;k−1i∗a[n−i]。
这条式子是怎么来的呢&＃xff1f;我的理解就是&＃xff0c;把每一次前缀和后的数组放在一起&＃xff0c;初始数组在第0行&＃xff0c;前缀和放在第一行&＃xff0c;如此类推。那么a[i]对a[n]贡献的系数就相当于每次可以从(x,y)走到(x&＃43;1,y&＃43;l)&＃xff08;l>&＃61;0&＃xff09;&＃xff0c;从(0,i)走到(k,n)的不同方案数。根据隔板法不难得到系数就是Cn−in−i&＃43;k−1Cn−i&＃43;k−1n−i。
上面的式子是一个卷积的形式&＃xff0c;可以用FFT来优化。但注意到这里的模数并不是NTT模数&＃xff0c;那么就要用到任意模数FFT。
具体来讲就是取M&＃61;P−−√M&＃61;P&＃xff0c;设f(x)&＃61;k(x)∗M&＃43;r(x)f(x)&＃61;k(x)∗M&＃43;r(x)&＃xff0c;那么有
f1(x)∗f2(x)&＃61;(k1(x)∗M&＃43;r1(x))∗(k2(x)∗M&＃43;r2(x))f1(x)∗f2(x)&＃61;(k1(x)∗M&＃43;r1(x))∗(k2(x)∗M&＃43;r2(x))
&＃61;M2∗k1(x)∗k2(x)&＃43;M∗(k1(x)∗r2(x)&＃43;k2(x)∗r1(x))&＃43;r1(x)∗r2(x)&＃61;M2∗k1(x)∗k2(x)&＃43;M∗(k1(x)∗r2(x)&＃43;k2(x)∗r1(x))&＃43;r1(x)∗r2(x)
总共需要做7次FFT且卷积后元素的数量级为O(nP)&＃xff0c;可以用long double来进行FFT而不是高精度。
一开始打完后被卡了精度&＃xff0c;后面把pi从define变成const就过了。。。
第一次手打复数类型&＃xff0c;发现比自带的要快将近一倍。。。

代码

#include
#include
#include
#include
#include
#includeusing namespace std;typedef long long LL;const int MOD&＃61;1000000007;
const int N&＃61;50005;
const long double pi&＃61;acos((long double)-1);struct com
{long double x,y;com operator &＃43; (const com &d) {return (com){x&＃43;d.x,y&＃43;d.y};}com operator - (const com &d) {return (com){x-d.x,y-d.y};}com operator * (const com &d) {return (com){x*d.x-y*d.y,x*d.y&＃43;d.x*y};}com operator / (const long double &d) {return (com){x/d,y/d};}
};int n,M,ny[N],c[N],a[N],k,L,rev[N*4];
com s1[N*4],s2[N*4],s3[N*4],k1[N*4],r1[N*4],k2[N*4],r2[N*4];int read()
{int x&＃61;0,f&＃61;1;char ch&＃61;getchar();while (ch<&＃39;0&＃39;||ch>&＃39;9&＃39;){if(ch&＃61;&＃61;&＃39;-&＃39;)f&＃61;-1;ch&＃61;getchar();}while (ch>&＃61;&＃39;0&＃39;&&ch<&＃61;&＃39;9&＃39;){x&＃61;x*10&＃43;ch-&＃39;0&＃39;;ch&＃61;getchar();}return x*f;
}void init()
{n&＃61;read();k&＃61;read();for (int i&＃61;0;i0]&＃61;ny[1]&＃61;c[0]&＃61;1;c[1]&＃61;k;for (int i&＃61;2;i1]*(k&＃43;i-1)%MOD;for (int i&＃61;2;i1]%MOD,c[i]&＃61;(LL)c[i]*ny[i]%MOD;
}void fft(com *a,int f)
{for (int i&＃61;0;iif (ifor (int i&＃61;1;i1){com wn&＃61;(com){cos(pi/i),f*sin(pi/i)};for (int j&＃61;0;j1)){com w&＃61;(com){1,0};for (int k&＃61;0;kif (f&＃61;&＃61;-1) for (int i&＃61;0;i}int main()
{init();M&＃61;sqrt(MOD);int lg&＃61;0;for (L&＃61;1;L<&＃61;n*2;L<<&＃61;1,lg&＃43;&＃43;);for (int i&＃61;0;i>1]>>1)|((i&1)<<(lg-1));for (int i&＃61;0;i0.0},r1[i]&＃61;(com){a[i]%M,0.0},k2[i]&＃61;(com){c[i]/M,0.0},r2[i]&＃61;(com){c[i]%M,0.0};fft(k1,1);fft(r1,1);fft(k2,1);fft(r2,1);for (int i&＃61;0;i1);fft(s2,-1);fft(s3,-1);for (int i&＃61;0;iint x1&＃61;(LL)(s1[i].x&＃43;0.5)%MOD*M*M%MOD,x2&＃61;(LL)(s2[i].x&＃43;0.5)%MOD*M%MOD,x3&＃61;(LL)(s3[i].x&＃43;0.5)%MOD;int ans&＃61;x1&＃43;x2;ans-&＃61;ans>&＃61;MOD?MOD:0;ans&＃43;&＃61;x3;ans-&＃61;ans>&＃61;MOD?MOD:0;printf("%d\n",ans);}return 0;
}