●BZOJ 3527 [Zjoi2014]力

题链：

http://www.lydsy.com/JudgeOnline/problem.php?id=3527

题解：

FFT求卷积。

$$\begin{aligned}
E_i&=\frac{F_i}{q_i}\\
&=\sum_{k<i}\frac{q_k}{(i-k)^2}-\sum_{k>i}\frac{q_k}{(i-k)^2}\\
&=\sum_{k=1}^{n}{q_kP(i-k)}
\end{aligned}$$

其中，

$$P(x)=\begin{cases}
x^{-2}\;\;\quad x>0\\
0\;\;\quad\quad x=0\\
-x^{-2}\quad x<0\\
\end{cases}$$

如上，就构造出了一个卷积的形式，可以用FFT求出所有的$E_i$。

另外为了避免下标出现负数，所以E的下标和P的下标都加上n，(同时为了FFT的方便，把k从0开始枚举)即： 

$$
E_i=E'_{i+n}=\sum_{k=0}^{n-1}{q_kP'(i+n-k)}\\
P'(x)=\begin{cases}
(x-n)^{-2}\;\;\quad x>n\\
0\quad \qquad\qquad x=0\\
-(n-x)^{-2}\;\quad x<n\\
\end{cases}
$$

代码：

 

#include<bits/stdc++.h>
#define MAXN 262145 
using namespace std;
const double Pi=acos(-1);
struct Complex{
	double real,image;
	Complex(double _real=0,double _image=0):real(_real),image(_image){}
	Complex operator - () const{return Complex(-real,-image);}
	friend Complex operator + (Complex A,Complex B){return Complex(A.real+B.real,A.image+B.image);}
	friend Complex operator - (Complex A,Complex B){return A+(-B);}
	friend Complex operator * (Complex A,Complex B){return Complex(A.real*B.real-A.image*B.image,A.image*B.real+A.real*B.image);}
};
Complex Q[MAXN],P[MAXN];
int pos[MAXN];
void FFT(Complex *Y,int n,int sign){
	for(int i=0;i<n;i++) if(i<pos[i]) swap(Y[i],Y[pos[i]]);
	for(int d=2;d<=n;d<<=1){
		Complex dw(cos(2*Pi/d),sin(sign*2*Pi/d)),w,tmp;
		for(int i=0;w=Complex(1,0),i<n;i+=d)
			for(int k=i;k<i+d/2;w=w*dw,k++)
				tmp=w*Y[k+d/2],Y[k+d/2]=Y[k]-tmp,Y[k]=Y[k]+tmp;
	}
}
int main(){
	int n,m,len; scanf("%d",&n);
	for(int i=0;i<n;i++) scanf("%lf",&Q[i].real);
	for(int k=1;k<n;k++) P[n-k].real=-1.0/(1.0*k*k),P[n+k].real=1.0/(1.0*k*k);
	m=2*n; for(n=1,len=0;n<m;n<<=1) len++;
	for(int i=0;i<n;i++) pos[i]=(pos[i>>1]>>1)|((i&1)<<(len-1));
	FFT(Q,n,1); FFT(P,n,1);
	for(int i=0;i<n;i++) Q[i]=Q[i]*P[i];
	FFT(Q,n,-1);
	for(int i=m/2;i<m;i++) printf("%.5lf\n",Q[i].real/n);
	return 0;
}