BZOJ3527 「ZJOI2014」力 FFT

问题描述

BZOJ3527


题解

(F(i)=q_i,G(i)=frac{1}{i^2}),FFT一下就行了。


(mathrm{Code})

#include<bits/stdc++.h>
using namespace std;

const double Pi=acos(-1);
const int maxn=150000;
int n,m;

struct CP{
	CP(double X=0,double Y=0) {x=X,y=Y;}
	double x,y;
	CP operator + (CP const &a) const
	{return CP(x+a.x,y+a.y);}
	CP operator - (CP const &a) const
	{return CP(x-a.x,y-a.y);}
	CP operator * (CP const &a) const
	{return CP(x*a.x-y*a.y,x*a.y+y*a.x);}
}f[maxn<<1],p[maxn<<1],q[maxn<<1];
int tr[maxn<<1];

void fft(CP *f,int type){
	for(int i=0;i<n;i++) if(i<tr[i]) swap(f[i],f[tr[i]]);
	for(int dlen=2;dlen<=n;dlen<<=1){
		int len=dlen>>1;
		CP w(cos(2*Pi/dlen),type*sin(2*Pi/dlen));
		for(int k=0;k<n;k+=dlen){
			CP buf(1,0);
			for(int l=k;l<k+len;l++){
				CP tt=buf*f[l+len];f[l+len]=f[l]-tt;
				f[l]=f[l]+tt;buf=buf*w;
			}
		}
	}
}

int main(){
	int N;
	scanf("%d",&N);
	for(int i=1;i<=N;i++) scanf("%lf",&f[i].x),q[N-i+1].x=f[i].x;
	for(int i=1;i<=N;i++) p[i].x=1.0/((double)i)/((double)i);
	for(n=1;n<=2*N;n<<=1);
//	for(int i=1;i<=n;i++) printf("**%lf %lf %lf
",f[i].x,q[i].x,p[i].x);
	for(int i=0;i<n;i++) tr[i]=(tr[i>>1]>>1)|((i&1)?n>>1:0);
//	for(int i=0;i<n;i++) tr[i]=(tr[i>>1]>>1)|((i&1)?N>>1:0);//µäÐÍ´íÎó 
	fft(f,1);fft(q,1);fft(p,1);
	for(int i=0;i<=n;i++) f[i]=f[i]*p[i],q[i]=q[i]*p[i];
	fft(f,-1);fft(q,-1);
	for(int i=0;i<=n;i++) f[i].x/=n,q[i].x/=n;
	for(int i=1;i<=N;i++) printf("%.8f
",f[i].x-q[N-i+1].x);
	return 0;
}
原文地址:https://www.cnblogs.com/liubainian/p/12148170.html