1. 程式人生 > >[bzoj3527][Zjoi2014]力_FFT

[bzoj3527][Zjoi2014]力_FFT

力 bzoj-3527 Zjoi-2014

題目大意:給定長度為$n$的$q$序列,定義$F_i=\sum\limits_{i<j}\frac{q_iq_j}{(i-j)^2}-\sum\limits_{i>j}\frac{q_iq_j}{(i-j)^2}$。求所有的$E_i=\frac{F_i}{q_i}$。

註釋:$1\le n\le 10^5$,$0\le q\le 10^9$。


想法:我們可以把$F_i$中每一項上的$q_i$刪掉因為我們求得$E_i$除掉了。

進而我們考慮如何求解$F$。

先看$j<i$的部分

$F_i=\sum\limits_{j=0}^{i-1} \frac{q_j}{(i-j)^2}$。

設$p(x)=\frac{1}{x^2}$。

所以$F_i=\sum\limits_{j=0}^{i-1} q_j\cdot p_{i-j}$。

緊接著我們強制令$p_0=0$,$F_i=\sum\limits_{j=0}^i q_j\cdot p_{i-j}$,可以用$FFT$加速。

接下來看$i<j$的部分。

此時$F_i=\sum\limits_{j=i+1}^{n-1} q_j\cdot p_{j-i}$。

bzoj2194一樣,這時我們將$p$序列翻轉,仍然可以用$FFT$加速。

之後把這兩部分加一起即可。

Code

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#define N 100010 
using namespace std; typedef double db;
const db pi=acos(-1);
db E[N<<2],q[N<<2],p[N<<2];
struct cp
{
	db x,y;
	cp() {x=y=0;}
	cp(db x_,db y_) {x=x_,y=y_;}
	cp operator + (const cp &a) const {return cp(x+a.x,y+a.y);}
	cp operator - (const cp &a) const {return cp(x-a.x,y-a.y);}
	cp operator * (const cp &a) const {return cp(x*a.x-y*a.y,x*a.y+y*a.x);}
}a[N<<2],b[N<<2],c[N<<2],d[N<<2];
void fft(cp *a,int len,int flg)
{
	int i,j,k,t;
	cp tmp,w,wn;
	for(i=k=0;i<len;i++)
	{
		if(i>k) swap(a[i],a[k]);
		for(j=len>>1;(k^=j)<j;j>>=1);
	}
	for(k=2;k<=len;k<<=1)
	{
		wn=cp(cos(2*pi*flg/k),sin(2*pi*flg/k));
		t=k>>1;
		for(i=0;i<len;i+=k)
		{
			w=cp(1,0);
			for(j=i;j<i+t;j++)
			{
				tmp=a[j+t]*w;
				a[j+t]=a[j]-tmp;
				a[j]=a[j]+tmp;
				w=w*wn;
			}
		}
	}
	if(flg==-1) for(i=0;i<len;i++) a[i].x/=len;
}
int main()
{
	int n; cin >> n ; for(int i=0;i<n;i++) scanf("%lf",&q[i]);
	for(int i=1;i<=n;i++) p[i]=(double)(1)/(1ll*i*i); p[0]=0;
	for(int i=0;i<n;i++) a[i].x=c[i].x=q[i];
	for(int i=0;i<n;i++) b[i].x=d[n-i-1].x=p[i];
	int len=1; while(len<=(n<<1)) len<<=1;
	fft(a,len,1); fft(b,len,1);
	for(int i=0;i<len;i++) a[i]=a[i]*b[i];
	fft(a,len,-1);
	for(int i=0;i<n;i++) E[i]=a[i].x;
	fft(c,len,1); fft(d,len,1);
	for(int i=0;i<len;i++) c[i]=c[i]*d[i];
	fft(c,len,-1);
	for(int i=0;i<n;i++) E[i]-=c[n+i-1].x;
	for(int i=0;i<n;i++) printf("%.3lf\n",E[i]);
	return 0;
}

小結:對於這兩種形式可以用$FFT$加速應該熟練掌握。