Description###
给出n个数qi,给出Fj的定义如下:
[F_j=sumlimits_{i<j} frac{q_iq_j}{(i-j)^2} - sumlimits_{i>j} frac{q_iq_j}{(i-1)^2}
]
令Ei=Fi/qi,求Ei.
Input###
第一行一个整数n。
接下来n行每行输入一个数,第i行表示qi。
n≤100000,0<qi<1000000000
Output###
n行,第i行输出Ei。与标准答案误差不超过1e-2即可。
Sample Input###
5
4006373.885184
15375036.435759
1717456.469144
8514941.004912
1410681.345880
Sample Output###
-16838672.693
3439.793
7509018.566
4595686.886
10903040.872
想法##
对FFT及卷积等还不是很熟,所以这道题还是参考的题解。
首先,原式化简得:
[E_j=sumlimits_{i<j} frac{q_i}{(i-1)^2} - sumlimits_{i>j} frac{q_i}{(i-j)^2}
]
设 (g_i= frac{1}{i^2}) ,则
[E_j=sumlimits_{i<j} q_ig_{j-i} - sumlimits_{i>j} q_i g_{j-i}
]
我们发现前面的求和式子就是一个卷积形式,可用fft
而后面的求和式子,若将下标i反过来,就也是一个标准的卷积形式,可用fft
之后就出来啦。
代码##
#include<cstdio>
#include<iostream>
#include<algorithm>
#include<cmath>
using namespace std;
const int N = 300005;
const double pi = 3.1415926535897932384626433832795;
struct c{
double r,i;
c() { r=i=0.0; }
c(double x,double y) { r=x; i=y; }
c operator + (const c &b) { return c(r+b.r,i+b.i); }
c operator += (const c &b) { return *this=*this+b; }
c operator - (const c &b) { return c(r-b.r,i-b.i); }
c operator -= (const c &b) { return *this=*this-b; }
c operator * (const c &b) { return c(r*b.r-i*b.i,r*b.i+b.r*i); }
c operator *= (const c &b) { return *this=*this*b; }
}a1[N],a2[N],b[N],x[N];
int l;
int r[N];
void fft(c A[],int ty){
for(int i=0;i<l;i++) x[r[i]]=A[i];
for(int i=0;i<l;i++) A[i]=x[i];
for(int i=2;i<=l;i<<=1){ /**/
c wn(cos(pi*2/i),ty*sin(pi*2/i));
for(int j=0;j<l;j+=i){
c w(1,0);
for(int k=j;k<j+i/2;k++){
c t=w*A[k+i/2];
A[k+i/2]=A[k]-t;
A[k]+=t;
w*=wn;
}
}
}
}
int n;
int main()
{
scanf("%d",&n);
for(int i=0;i<n;i++) {
scanf("%lf",&a1[i].r);
a2[n-i-1].r=a1[i].r;
}
for(int i=1;i<n;i++) b[i].r=1.0/((double)i*i);
l=1;
while(l<n*2) l<<=1;
for(int i=0;i<l;i++) r[i]=(r[i>>1]>>1)|((i&1)*(l>>1));
fft(a1,1); fft(a2,1); fft(b,1);
for(int i=0;i<l;i++) {
a1[i]*=b[i];
a2[i]*=b[i];
}
fft(a1,-1); fft(a2,-1);
for(int i=0;i<n;i++)
printf("%.3lf
",a1[i].r/l-a2[n-i-1].r/l);
return 0;
}