思路
颓柿子的题目
要求求这样的一个式子
\[F_j=\sum_{i<j}\frac{q_iq_j}{(i-j)^2}-\sum_{i>j}\frac{q_iq_j}{(i-j)^2}
\]
\]
令\(E_i=\frac{F_i}{q_i}\),求所有的\(E_i\)
对于Ei,显然可以
\[E_i=\sum_{j=0}^{i-1}\frac{q_j}{(i-j)^2}-\sum_{j=i+1}^n\frac{q_j}{(i-j)^2}
\]
\]
前后没什么关联,可以分开考虑,首先考虑前面部分
\[\sum_{j=0}^{i-1}\frac{q_j}{(i-j)^2}
\]
\]
设\(f(i)=q_i\),\(g(i)=\frac{1}{i^2}\)
则前面一部分变为
\[\sum_{j=0}^{i-1}f(j)g(i-j)
\]
\]
变成了卷积的形式
后面一部分是
\[\sum_{j=i+1}^{n}\frac{q_j}{(i-j)^2}
\]
\]
同理变成
\[\sum_{j=i+1}^nf(j)g(j-i)
\]
\]
似乎没什么办法,但是我们可以把\(f\)整个反过来(就是\(f(1)\)与\(f(n)\)交换,与\(f(2)\)与\(f(n-1)\)交换)
设交换之后的\(f\)为\(f'\),\(f'_i=f_{n-i+1}\)
则有
\[\sum_{j=0}^{i-1}f'(j)g(j-i)
\]
\]
则原来的第i项和反转后的第\(n-i+1\)项相同
对应相减即可
代码
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#define int long long
using namespace std;
const double pi = acos(-1.0);
struct complex{//a+bi
double a,b;
complex operator + (complex &bx){
return (complex){a+bx.a,b+bx.b};
}
complex operator - (complex &bx){
return (complex){a-bx.a,b-bx.b};
}
complex operator * (complex &bx){
return (complex){a*bx.a-b*bx.b,b*bx.a+a*bx.b};
}
complex conj(void){
return (complex){a,-b};
}
};
int n;
double q[401000],ans[401000];
complex inv[401000],wnk[401000],a[401000],b[401000];
void init(int len){
for(int i=0;i<len;i++){
wnk[i]=(complex){cos(2*pi*i/len),sin(2*pi*i/len)};
inv[i]=wnk[i].conj();
}
}
void FFT(complex *a,complex *opt,int n){
int lim=0;
while((1<<lim)<n){
lim++;
}
n=(1<<lim);
for(int i=0;i<n;i++){
int t=0;
for(int j=0;j<lim;j++)
if((i>>j)&1)
t|=(1<<(lim-j-1));
if(t<i)
swap(a[i],a[t]);
}
for(int i=2;i<=n;i<<=1){
int len=i/2;
for(complex *j=a;j!=a+n;j+=i){
for(int k=0;k<len;k++){
complex t=j[k+len]*opt[n/i*k];
j[k+len]=j[k]-t;
j[k]=j[k]+t;
}
}
}
}
void Do_FFT(int n){
int lim=0;
while((1<<lim)<n){
lim++;
}
n=(1<<lim);
init(n);
FFT(a,wnk,n);
FFT(b,wnk,n);
for(int i=0;i<n;i++)
a[i]=a[i]*b[i];
FFT(a,inv,n);
for(int i=0;i<n;i++)
a[i].a/=n;
}
signed main(){
// freopen("1.in","r",stdin);
// freopen("test.out","w",stdout);
scanf("%lld",&n);
for(int i=1;i<=n;i++){
scanf("%lf",&q[i]);
}
for(int i=1;i<=n;i++){
a[i].a=q[i];
b[i].a=1.0/(i*i);
a[i].b=0;
b[i].b=0;
}
Do_FFT(2*n+4);
for(int i=1;i<=n;i++)
ans[i]=a[i].a;
memset(a,0,sizeof(a));
memset(b,0,sizeof(b));
for(int i=1;i<=n;i++){
a[i].a=q[i];
b[i].a=1.0/(i*i);
a[i].b=0;
b[i].b=0;
}
for(int i=1,j=n;i<j;i++,j--)
swap(a[i],a[j]);
Do_FFT(2*n+4);
for(int i=1;i<=n;i++)
ans[i]-=a[n-i+1].a;
for(int i=1;i<=n;i++)
printf("%.3lf\n",ans[i]);
return 0;
}