题目
给出n个数qi，给出Fj的定义如下：
令Ei=Fi/qi，求Ei.
思路
年轻人的第一道FFT。
[E_j=sum_{i<j}frac{q_i}{(i-j)^2}-sum_{i>j}frac{q_i}{(i-j)^2}
]
设：
[f[i]=q_i\
g[i]=frac{1}{i^2}
]
那么原式可以化为：
[E_j=sum_{i<j}f[i]g[j-i]-sum_{i>j}f[i]g[i-j]
]
也就是:
[E_j=sum_{i=0}^{j-1}f[i]g[j-i]-sum_{i=j+1}^{n}f[i]g[i-j]
]
前半部分是卷积的形式，显然可以用FFT维护，重点是后面的部分。
我们来变一下形：
[sum_{i=j+1}^{n}f[i]g[i-j]\
令i=i-j\
sum_{i=0}^{n-j}f[i+j]g[i]
]
然后将(i+j+i)凑成(n-j),也就是令(f'[n-j-i]=f[i+j])
再简单变形：(f'[i]=f[n-i])
于是(f'[])也构造出来了，这样两边都化为了卷积的形式，就可以用FFT求解了。
代码
#include<bits/stdc++.h>
#define M 400005
using namespace std;
int n;
struct Complex{
double x,y;
Complex(){}
Complex(double _x,double _y):x(_x),y(_y){}
Complex operator + (const Complex &res) const{
return (Complex){x+res.x,y+res.y};
}
Complex operator - (const Complex &res) const{
return (Complex){x-res.x,y-res.y};
}
Complex operator * (const Complex &res) const{
return (Complex){x*res.x-y*res.y,x*res.y+y*res.x};
}
};
Complex A[M],B[M];
double pi=acos(-1.0);
void FFT(Complex *y,int n,int f){
if(n==1)return;
Complex L[n>>1],R[n>>1];
for(int i=0;i<n;i+=2)L[i>>1]=y[i],R[i>>1]=y[i+1];
FFT(L,n>>1,f);FFT(R,n>>1,f);
Complex wn(cos(2*pi/n),f*sin(2*pi/n)),w(1,0);
for(int i=0;i<(n>>1);i++,w=w*wn){
y[i]=L[i]+w*R[i];
y[i+(n>>1)]=L[i]-w*R[i];
}
}
double q[M],b[M],C[M];
int main(){
cin>>n;n--;
for(int i=0;i<=n;i++)scanf("%lf",&q[i]);
for(int i=0;i<=n;i++)A[i]=(Complex){q[i],0};
for(int i=1;i<=n;i++){
b[i]=1.0/i/i;
B[i]=Complex(b[i],0);
}
int nn=n,mm=n;
mm+=nn;
for(nn=1;nn<=mm;nn<<=1);
FFT(A,nn,1);FFT(B,nn,1);
for(int i=0;i<=nn;i++)A[i]=A[i]*B[i];
FFT(A,nn,-1);
for(int i=0;i<=n;i++)C[i]=A[i].x/nn;
for(int i=0;i<=nn;i++)A[i]=B[i]=Complex(0,0);
for(int i=0;i<=n;i++)A[i]=Complex(q[n-i],0);
for(int i=1;i<=n;i++)B[i]=Complex(b[i],0);
FFT(A,nn,1);FFT(B,nn,1);
for(int i=0;i<=nn;i++)A[i]=A[i]*B[i];
FFT(A,nn,-1);
for(int i=0;i<=n;i++)C[i]-=A[n-i].x/nn;
for(int i=0;i<=n;i++)printf("%.3f
",C[i]);
return 0;
}