题目
P3338 [ZJOI2014]力
做法
普通卷积形式为：(c_k=sumlimits_{i=1}^ka_ib_{k-i})
其实一般我们都是用(i=0)开始的，但这题比较特殊，忽略那部分，其他的直接按下标存下来，反正最后的答案是不变的
好了步入正题吧，我们定义 $$F_j=sumlimits_{i<j}dfrac{q_iq_j}{(i-j)^2}-sumlimits_{i<j}dfrac{q_iq_j}{(i-j)^2}$$
求(E_i=dfrac{F_i}{q_i})
显然$$E_j=sumlimits_{i<j}dfrac{q_i}{(i-j)^2}-sumlimits_{i<j}dfrac{q_i}{(i-j)^2}$$
[E_j=sum_{i=1}^{j-1}dfrac{q_i}{(i-j)^2}-sum_{i=j+1}^{n}dfrac{q_i}{(i-j)^2}
]
令(f_i=q_i,g_i=dfrac{1}{i^2})，特别地，(g_0=0)，则有：$$E_j=sum_{i=1}^{j-1}f_ig_{i-j}-sum_{i=j+1}^{n}f_ig_{i-j}$$
左边部分很简单就能化成卷积形式：$$sum_{i=1}^{j-1}f_ig_{i-j}=sum_{i=1}^{j-1}f_ig_{j-i}=sum_{i=1}^{j}f_ig_{j-i}$$
右边部分：$$sumlimits_{i=j+1}^{n}f_ig_{i-j}=sumlimits_{i=1}^{n-j}f_{i+j}g_{i}$$
令(p_i=f_{n-i})，则(p_{n-i-j}=f_{i+j})，则有：$$sumlimits_{i=1}^{n-j}f_{i+j}g_{i}=sum_{i=1}^{n-j} p_{n-i-j}g_i$$
都化成卷积形式了直接FFT
My complete code
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cmath>
using namespace std;
const double Pi=acos(-1.0);
const int maxn=300000;
struct complex{
double x,y;
complex (double xx=0,double yy=0){
x=xx,y=yy;
}
}f[maxn],p[maxn],g[maxn];
int n,limit,L;
int r[maxn];
complex operator + (complex x,complex y){
return complex(x.x+y.x,x.y+y.y);
}
complex operator - (complex x,complex y){
return complex(x.x-y.x,x.y-y.y);
}
complex operator * (complex x,complex y){
return complex(x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x);
}
inline void FFT(complex *A,int type){
for(int i=0;i<limit;++i)
if(i<r[i])
swap(A[i],A[r[i]]);
for(int mid=1;mid<limit;mid<<=1){
complex WN(cos(Pi/mid),type*sin(Pi/mid));
for(int R=mid<<1,j=0;j<limit;j+=R){
complex w(1,0);
for(int k=0;k<mid;++k,w=w*WN){
complex x=A[j+k],y=w*A[j+mid+k];
A[j+k]=x+y,
A[j+mid+k]=x-y;
}
}
}
}
int main(){
scanf("%d",&n);
for(int i=1;i<=n;++i){
scanf("%lf",&f[i].x),f[i].y=0,
g[i].x=(double)1.0/i/i;
p[n-i].x=f[i].x;
}
limit=1,L=0;
while(limit<n+n)
limit<<=1,
++L;
for(int i=0;i<limit;++i)
r[i]=(r[i>>1]>>1)|((i&1)<<(L-1));
FFT(f,1),FFT(p,1),FFT(g,1);
for(int i=0;i<limit;++i)
f[i]=f[i]*g[i];
for(int i=0;i<limit;++i)
p[i]=p[i]*g[i];
FFT(f,-1),FFT(p,-1);
for(int i=1;i<=n;++i)
printf("%.3lf
",f[i].x/limit-p[n-i].x/limit);
return 0;
}/*
*/