此题bzoj未提供题意--
把定义式直接除以qi,右边的东西再弄成逆序,就是两个卷积,直接套FFT
第一次写FFT把不理解的地方就补在这里了
感觉Rader算法求倒位序很麻烦。。。网上很多人都是用预处理rev数组,个人感觉简单暴力
然后FFT中使用的多项式要补齐到2^k项,这样才能保证蝴蝶操作的正确性
此处2^k取一个刚好>=2*n的东西
以及代码中的h数组不能多赋值,,该是n项就是n项,不然点值表达式是错的
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
using namespace std;
const int maxn = 4E5 + 10;
typedef double DB;
const DB PI = acos(-1.0);
struct Virt{
DB r,i;
Virt(){}
Virt(DB r,DB i): r(r),i(i){}
Virt operator + (const Virt &b) {return Virt(r + b.r,i + b.i);}
Virt operator - (const Virt &b) {return Virt(r - b.r,i - b.i);}
Virt operator * (const Virt &b) {return Virt(r*b.r - i*b.i,r*b.i + i*b.r);}
Virt operator * (const DB &t) {return Virt(r*t,i*t);}
}q[maxn],p[maxn],h[maxn],e1[maxn],e2[maxn],a[maxn];
int n,m,N,L,rev[maxn],dig[maxn];
DB Now[maxn];
void FFT(Virt *A,int on)
{
for (int i = 0; i < N; i++) a[i] = A[rev[i]];
for (int i = 0; i < N; i++) A[i] = a[i];
DB t1 = on*2.00;
for (int k = 2; k <= N; k <<= 1) {
DB t2 = k;
Virt wn = Virt(cos(PI*t1/t2),sin(PI*t1/t2));
for (int i = 0; i < N; i += k) {
Virt w = Virt(1.00,0.00);
for (int j = i; j < i + (k>>1); j++) {
Virt u = A[j];
Virt t = w*A[j+(k>>1)];
A[j] = u + t;
A[j+(k>>1)] = u - t;
w = w*wn;
}
}
}
if (on == -1)
for (int i = 0; i < N; i++)
A[i].r /= (DB)(N);
}
int main()
{
#ifdef DMC
freopen("DMC.txt","r",stdin);
#endif
cin >> n;
for (N = 1; N < n; N <<= 1,++L);
N <<= 1; ++L;
for (int i = 0; i < N; i++) {
for (int x = i,len = 0; x; x >>= 1)
dig[len++] = x&1;
for (int j = 0; j < L; j++)
rev[i] = rev[i]*2 + dig[j];
}
for (int i = 0; i < n; i++) scanf("%lf",&Now[i]);
for (int i = 0; i < N; i++)
p[N-i-1] = q[i] = Virt(Now[i],0);
for (int i = 1; i < n; i++)
h[i] = Virt(1.00/(DB)(i)/(DB)(i),0);
FFT(h,1); FFT(q,1); FFT(p,1);
for (int i = 0; i < N; i++)
e1[i] = q[i]*h[i];
for (int i = 0; i < N; i++)
e2[i] = p[i]*h[i];
FFT(e1,-1); FFT(e2,-1);
for (int i = 0; i < n; i++)
printf("%.8lf\n",e1[i].r - e2[N-i-1].r);
return 0;
}