[ZJOI2014] 力解题报告

Description

F_j~=~\sum_{i = 1}^{j - 1} \frac{q_i \times q_j}{(i - j)^2}~-~\sum_{i = j + 1}^{n} \frac{q_i \times q_j}{(i - j)^2}

E_i~=~\frac{F_i}{q_i}

对 $1 \leq i \leq n$ ，求 $E_i$ 的值。

当两边同时加上j时，仍不变。

F_j~=~\sum_{i = 1}^{j } \frac{q_i \times q_j}{(i - j)^2}~-~\sum_{i = j}^{n} \frac{q_i \times q_j}{(i - j)^2}

E_j=\sum_{i=1}^j\frac {q_i} {(i-j)^2}-\sum_{i=j}^n\frac {q_i} {(i-j)^2}

我们规定

f(n)= q_n,g(n)=\frac 1 {n^2}

因为 $(i-j)^2=(j-i)^2$

E_j=\sum_{i=1}^jf(i)\times g(j-i)-\sum_{i=j}^nf(i)\times g(i-j)

不妨令 $f(0)=0,g(0)=0$ ，那么左边变成了卷积形式。现在来看右边。

\sum_{i=j}^nf(i)\times g(i-j)=\sum_{i=0}^{n-j}f(i+j)\times g(i)

设 $k=n-j$ ，我们引入一个 $f$ 的反转后的多项式 $f′$ ，即让 $f'(i)=f(n-i)$ 。

\sum_{i=0}^{k}f'(n-i-j)\times g(i)=\sum_{i=0}^kf(k-i)\times g(i)

至此，两边都是卷积形式，可以直接FFT。

Code：

#include<cmath>
#include<algorithm>
#include<cstdio>
#include<cstring>
using namespace std;
#define R register
const int N=1e5+9;
const double pi=3.1415926535898;
struct complex
{
	complex(double xx=0,double yy=0){x=xx,y=yy;}
	double x,y;
};
complex operator + (complex a,complex b){return complex(a.x+b.x,a.y+b.y);}
complex operator - (complex a,complex b){return complex(a.x-b.x,a.y-b.y);}
complex operator * (complex a,complex b){return complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}
complex a[N],b[N],c[N];
int r[N],l,m=1;
inline int read()
{
	int res=0,f=1;
	char ch=getchar();
	while(ch<'0'||ch>'9'){
		if(ch=='-') f=-1;
		ch=getchar();
	}
	while(ch>='0'&&ch<='9'){
		res=(res<<1)+(res<<3)+ch-'0';
		ch=getchar();
	}
	return res*f;
}

	int n;
inline void FFT(complex *f,int op)
{
	for(R int i=0;i<m;++i) if(i<r[i]) swap(f[i],f[r[i]]);
	for(R int mid=1;mid<m;mid<<=1)
	{
		complex Wn=complex(cos(pi/mid),op*sin(pi/mid));
		for(R int rr=mid<<1,i=0;i<m;i+=rr)
		{
			complex w=complex(1,0);
			for(R int k=0;k<mid;k++,w=w*Wn)
			{
				complex X=f[k+i],Y=f[k+mid+i]*w;
				f[k+i]=X+Y;
				f[k+i+mid]=X-Y;
			}
		}
	}
	if(op==-1) for(int i=0;i<m;i++) f[i].x/=m;
}

int main()
{
	n=read();
	for(R int i=1;i<=n;++i)
	{
		scanf("%lf",&a[i].x);
		c[n-i].x=a[i].x;
		b[i].x=(double)(1.0/i/i);
	}
	while(m<=(n<<1))  ++l,m<<=1;
	for(R int i=0;i<m;++i) r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
	FFT(a,1);
	FFT(b,1);
	FFT(c,1);
	for(R int i=0;i<m;i++) a[i]=a[i]*b[i],c[i]=c[i]*b[i];
	FFT(a,-1);FFT(c,-1);
	for(R int i=1;i<=n;++i)
		printf("%.3lf\n",a[i].x-c[n-i].x);
	return 0;
}