给出$n$个数$q_i$,给出$F_i$的定义如下:
$F_i=\sum_{j<i} \frac {q_iq_j}{(i-j)^2}-\sum_{j>i} \frac{q_iq_j}{(i-j)^2}$
求所有$E_i=\frac {F_i} {q_i}$
怎么看这题都是非常的奇怪啊,在上面套了个$q_i$在求$E_i$的时候又除掉了。。
那么$E_i=\sum_{j<i} \frac {q_j}{(i-j)^2}-\sum_{j>i} \frac{q_j}{(i-j)^2}$
考虑构造卷积
揣摩一下$2$年前毛主力的心路历程:
这个卷积可以直接秒掉啊!
新开一个数组$p_j=\begin{cases} -(i-j)^{-2} && \text{j<n-i} \\ (i-j)^{-2} && \text{j>n-i} \\ 0 && \text{j=n-i}\end{cases}$
$E_i=\sum_{j=1}^n p_j q_{n-j}$
哎呀好像每次都要重新算怎么办
哎呀反正和$E_i$没什么关系只要得到那个卷积式子即可
为了让$q_i$不被实时计算,不如无脑加个$i$好了
于是变成了这样的东西:$E_{i}=\sum_{j=1}^{n}p_jq_{n+i-j}$
又当$j>n$时,$p_j=0$
所以式子还可以进一步扩展:
$E_i=\sum_{j=1}^{n+i}p_jq_{n+i-j}$
有点像卷积的形式了,把$E_i$换成$A_{i+n}$
于是就构成了卷积形式$A_{i+n}=\sum_{j=1}^{n+i}p_jq_{n+i-j}$
于是$FFT$一发可以得到$A$咯,再无脑减个$n$得到$E$
浙江OI真是简单
#include<iostream>
#include<string.h>
#include<algorithm>
#include<stdio.h>
#include<math.h>
#include<complex>
using namespace std;
const int M=5e5+5;
typedef double db;
typedef complex<db> comp;
const db twpi=3.14159265358*2;
comp A[M],B[M],C[M];
int rev[M];
inline void DFT(comp *res,int n,int sg){
for(int i=0;i<n;++i)
if(rev[i]>i)swap(res[rev[i]+1],res[i+1]);
for(int s=2;s<=n;s<<=1){
comp wm(cos(twpi/s),sin(twpi/s)*sg);
for(int i=1;i<=n;i+=s){
comp w(1,0);
for(int j=0;j<s>>1;++j,w*=wm){
comp a=res[i+j],b=res[i+j+(s>>1)];
res[i+j]=a+b*w,res[i+j+(s>>1)]=a-b*w;
}
}
}
}
#define bug(x) cout<<#x<<"="<<x<<" "
#define debug(x) cout<<#x<<"="<<x<<endl
int main(){
int n;cin>>n;
db a;
for(int i=1;i<=n;++i)
scanf("%lf",&a),A[i]=a;
int w=1,s=0;
for(;w<=(n<<1);w<<=1,++s);
for(int i=1;i<w;++i)
rev[i]=rev[i>>1]>>1|(i&1)<<s-1;
for(int i=1;i<=w;++i){
if(i<n)B[i]=-1.0/((db)(n-i)*(n-i));
else if(i>n)B[i]=1.0/((db)(n-i)*(n-i));
else B[i]=0;
}
DFT(A,w,1),DFT(B,w,1);
for(int i=1;i<=w;++i)
C[i]=A[i]*B[i];
DFT(C,w,-1);
for(int i=1;i<=w;++i)C[i]/=w;
for(int i=1;i<=n;++i)
printf("%.3lf\n",C[i+n-1].real());
return 0;
}