bzoj 2194: 快速傅立叶之二 （FFT）

2194: 快速傅立叶之二

Time Limit: 10 Sec   Memory Limit: 259 MB
Submit: 1178   Solved: 682
[ Submit][ Status][ Discuss]

Description

请计算C[k]=sigma(a[i]*b[i-k]) 其中 k < = i < n ，并且有 n < = 10 ^ 5。 a,b中的元素均为小于等于100的非负整数。

Input

第一行一个整数N,接下来N行，第i+2..i+N-1行，每行两个数，依次表示a[i],b[i] (0 < = i < N)。

Output

输出N行，每行一个整数，第i行输出C[i-1]。

Sample Input

5
3 1
2 4
1 1
2 4
1 4

Sample Output

24
12
10
6
1

HINT

Source

[ Submit][ Status][ Discuss]

题解：FFT

C(x)=A(x)*B(x)

C(x)=sigma(j=0..2n-2) cj*x^j

cj=sigma(k=0..j)ak*b(j-k) 快速傅里叶变换实际就是对两个高次函数求卷积，两函数的傅里叶变换的乘积等于它们卷积后的傅里叶变换

我们发现两个高次函数的下标和是相等的。

但是这道题中ck=sigma a[i]*b[i-k] k<=i<n

可以这个式子的下标和不相等，但是ck=simga a[i]*b[n-(i-k)] 这样下标和就是n+k是相等的，所以我们可以经b数组翻转，使b[i]中存储原本b[n-i]中的信息。这样就可以用快速傅里叶变换来求解，那么答案就是c[n+k] k<n(原本的n) 

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
#include<cstring>
#define N 300003
#define pi acos(-1)
using namespace std;
struct data{
	double x,y; 
	data (double X=0,double Y=0) {
		x=X,y=Y;
	}
}a[N],b[N],c[N];
data operator +(data a,data b){  return data(a.x+b.x,a.y+b.y); }
data operator -(data a,data b){  return data(a.x-b.x,a.y-b.y); }
data operator *(data a,data b){  return data(a.x*b.x-a.y*b.y,a.y*b.x+a.x*b.y); }
int n,m,L,R[N];
char s[N];
void fft(data a[N],int opt)
{
	for (int i=0;i<n;i++) 
	if (i<R[i]) swap(a[R[i]],a[i]);
	for (int i=1;i<n;i<<=1) {
		data wn=data(cos(pi/i),opt*sin(pi/i));
		for (int p=i<<1,j=0;j<n;j+=p) {
			data w=data(1,0);
			for (int k=0;k<i;k++,w=w*wn) {
				data x=a[j+k],y=w*a[j+k+i];
				a[j+k]=x+y; a[j+k+i]=x-y;
			}
		}
	}
}
int main()
{
	freopen("a.in","r",stdin);
	scanf("%d",&n); n--;
	for (int i=0;i<=n;i++) scanf("%lf%lf",&a[i].x,&b[n-i].x);
	m=2*n; for (n=1;n<=m;n<<=1) L++;
	for (int i=0;i<n;i++)
	 R[i]=(R[i>>1]>>1)|((i&1)<<(L-1));
	fft(a,1); fft(b,1);
	for (int i=0;i<=n;i++) a[i]=a[i]*b[i];
	fft(a,-1);
	for (int i=m/2;i<=m;i++) printf("%d\n",(int)(a[i].x/n+0.5));
}