寒鸽儿のBlog

view

1. 学会了fft
2. 过了fft板子

detail

fft板子

talk

终于学会了fft , 在某位L姓带佬的指导下。
话说能给我这种理解能力和水平巨差的人全程讲下来真的不容易qaq
稍微记一点关键的东西备忘 , 实际上多数还是源于该带佬的blog

note

没写完 , 咕咕咕~

code

void dft(cp *a, int n) {
	if(n == 1) return;
	int hf = n >> 1;
	cp a0[hf], a1[hf];
	for(int i = 0; i < hf; ++i)	{
		a0[i] = a[i<<1];
		a1[i] = a[i<<1|1];
	}
	dft(a0, hf);
	dft(a1, hf);
	for(int i = 0; i < hf; ++i) {
		cp omg = cp(cos(2*PI*i/n), sin(2*PI*i/n));
		a[i] = a0[i] + omg * a1[i];
		a[i+hf] = a0[i] - omg * a1[i];	
	}
}

void idft(cp *a, int n) {
	reverse(a + 1, a + n);
	dft(a, n);
	for(int i = 0; i < n; ++i) a[i] = a[i] / n;
}

P3803 fft板子

事项 :

1. 需要使用二进制翻转进行卡常 , 通过观察容易发现在fft递归到底层之后的下标值为原下标的二进制翻转(长度为序列长度) , 该过程可以使用 $\mathcal{O(n)}$ 递推求解 ，即 t[i] = t[i>>1]>>1 + (i&n?n>>1:0)
2. 复数单位根应当采用递推的形式去做 , 即 $\omega_{n}^i = \omega_{n}^{i-1}*\omega_{n}$
   好的操作
   
   不优的操作
   

code

#include <bits/stdc++.h>

using namespace std;

const double PI = acos(-1);
const int maxn = 2345678;

struct cp {
	double x, y;
	cp() { x = y = 0; }
	cp(double X, double Y) : x(X), y(Y) {}
	cp operator + (const cp& o) const {
		return cp(x + o.x, y + o.y);
	}
	cp operator - (const cp& o) const {
		return cp(x - o.x, y - o.y);
	}
	cp operator * (const cp& o) const {
		return cp(x*o.x-y*o.y, x*o.y+y*o.x);	
	}
	void operator *= (const cp& o) {
		double tx = x * o.x - y * o.y;
		y = x * o.y + y * o.x;
		x = tx;
	}
	void operator /= (const int& n) {
		x /= n; y /= n;
	}
}a[maxn], b[maxn];

inline void flip(int &x) {
	int t = 1;
	while(t < x) t <<= 1;
	x = t;
}

int target[maxn];

void dft(cp a[], int n) {
	for(int i = 0; i < n; ++i) if(i <= target[i]) swap(a[i], a[target[i]]);
	for(int m = 2; m <= n; m <<= 1) {
		int hf = m >> 1;
		for(int i = 0; i < n; i += m) {
			cp w0 = cp(cos(2*PI/m), sin(2*PI/m));
			cp omg = cp(1, 0);
			for(int k = 0; k < hf; ++k, omg *= w0) {
				cp t = omg * a[i+k+hf];
				a[i+k+hf] = a[i+k] - t;
				a[i+k] = a[i+k] + t;
			}
		}
	}
}

void idft(cp a[], int n) {
	reverse(a + 1, a + n);
	dft(a, n);
	for(int i = 0; i < n; ++i) a[i] /= n;
}

int main() {
	ios_base::sync_with_stdio(false);
	cin.tie(NULL);
	int n, m;
	cin >> n >> m;
	for(int i = 0; i <= n; ++i) cin >> a[i].x;
	for(int i = 0; i <= m; ++i) cin >> b[i].x;
	int len = n + m + 1;
	flip(len);
	for(int i = 0; i < len; ++i)
		target[i] = target[i>>1]>>1 | (i&1?len>>1:0);
	dft(a, len);
	dft(b, len);
	for(int i = 0; i < len; ++i) a[i] *= b[i];
	idft(a, len);
	for(int i = 0; i <= n + m; ++i) cout << (int)(a[i].x + 0.5) << ' ';
	cout << endl;
	return 0;
}

thought

今天被一些事情缠住了 , 没怎么做题qaq
明天下午好像还要去兵役体检啊QaQ , 明明我近视500°的说 , 一下午又要没了(心疼)
明天上午有ACM集训的白嫖课 , 开心~