●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;

}

巴特西

●BZOJ 3527 [Zjoi2014]力

最新文章

热门文章