# FFT是个啥？

## 总述

$$A(x)=A_0+A_1x+A_2x^2+A_3x^3+...+A_{n-1}x^{n-1}$$

$$B(x)=B_0+B_1x+B_2x^2+B_3x^3+...+B_{m-1}x^{m-1}$$

## 点值表达式

FFT即为取若干个特殊的值进行点值表达

## n次单位复数根

n次单位复数根：即$\omega$使得$\omega^n=1$

（上述公式由于复数运算性质与指数运算性质大致相同，可以用泰勒公式证明，这里挖一个坑

## 一些定理

$$(\omega_n^k)^2=\omega_{n/2}^k$$

$$(\omega_n^{k+n/2})^2=\omega_n^{2k+n}=\omega_n^{2k}\cdot\omega_n^n=\omega_{n/2}^k$$

## DFT

$$y=DFT_n(a)$$

## FFT

$$A^{[0]}(x)=a_0+a_2x+a_4x^2+...+a_{n-2}x^{n/2-1}$$

$$A^{[1]}(x)=a_1+a_3x+a_5x^2+...+a_{n-1}x^{n/2-1}$$

$$A(x)=A^{[0]}(x^2)+xA^{[1]}(x^2)$$

$$y_k=y^{[0]}_k+\omega_n^ky^{[1]}_k$$

$$y_{k+n/2}=y^{[0]}_k-\omega_n^ky^{[1]}_k$$

$y_k=y^{[0]}_k+\omega_n^ky^{[1]}_k$

$=A^{[0]}(\omega_{n/2}^k)+\omega_n^kA^{[1]}(\omega_{n/2}^k)$

$=A^{[0]}((\omega_n^k)^2)+\omega_n^kA^{[1]}((\omega_n^k)^2)$

$=A(\omega_n^k)$

$y_{k+n/2}=y^{[0]}_k-\omega_n^ky^{[1]}_k$

$=A^{[0]}(\omega_{n/2}^k)+\omega_n^{k+n/2}A^{[1]}(\omega_{n/2}^k)$

$=A^{[0]}((\omega_n^{k+n/2})^2)+\omega_n^{k+n/2}A^{[1]}((\omega_n^{k+n/2})^2)$

$=A(\omega_n^{k+n/2})$

## 插值

$$\left[\begin{array}{aaa}y_0 \\ y_1 \\ y_2 \\ y_3 \\ \vdots \\ y_{n-1}\end{array}\right]=\left[\begin{array}{ccc}1 & 1 & 1 & 1 \ldots & 1 \\1 & \omega_n & \omega_n^2 & \omega_n^3 & \ldots & \omega_n^{n-1} \\1 & \omega_n^2 & \omega_n^4 & \omega_n^6 & \ldots & \omega_n^{2(n-1)} \\1 & \omega_n^3 & \omega_n^6 & \omega_n^9 & \ldots & \omega_n^{3(n-1)} \\\vdots & \vdots & \vdots & \vdots & \ddots & \vdots \\1 & \omega_n^{n-1} & \omega_n^{2(n-1)} & \omega_n^{3(n-1)} & \ldots & \omega_n^{(n-1)(n-1)}\end{array}\right]\left[\begin{array}{ppp}a_0 \\ a_1 \\ a_2 \\ a_3 \\ \vdots \\ y_{n-1}\end{array}\right]$$

$$\left[V_n^{-1}V_n\right]_jj'=\sum\limits_{k=0}^{n-1}(\omega_n^{-kj}/n)(\omega_n^{kj'})=\sum\limits_{k=0}^{n-1}\omega_n^{k(j'-j)}/n$$

$$a_j=\dfrac{1}{n}\sum\limits_{k=0}^{n-1}y_k\cdot \omega_n^{-kj}$$

## 高效实现FFT

$${a_0,a_4,a_2,a_6,a_1,a_5,a_3,a_7}$$

## Code

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <cmath>
using namespace std;
const int Maxn = 400010;
struct cp {
double r, i;
cp (double r = 0.0, double i = 0.0) : r(r), i(i) {}
}A[Maxn], B[Maxn]; int An, Bn, N;
cp operator +(cp a, cp b) { return cp(a.r+b.r, a.i+b.i); }
cp operator -(cp a, cp b) { return cp(a.r-b.r, a.i-b.i); }
cp operator *(cp a, cp b) { return cp(a.r*b.r-a.i*b.i, a.r*b.i+a.i*b.r); }
const double PI = acos(-1);
void fft(cp *a, int op) {
int i, j, k;
j = 0;
for(i = 0; i < N; i++){
if(i < j) swap(a[i], a[j]);
k = N >> 1;
while(j & k){ j -= k; k >>= 1; }
j += k;
}
for(i = 2; i <= N; i <<= 1){
cp wn = cp(cos(2.0*PI/i), op*sin(2.0*PI/i));
for(j = 0; j < N; j += i){
cp w = cp(1, 0);
for(k = j; k < j+i/2; k++){
cp x = a[k], y = w*a[k+i/2];
a[k] = x+y; a[k+i/2] = x-y;
w = w*wn;
}
}
}
if(op == -1) for(i = 0; i < N; i++) A[i].r /= N;
}
int main() {
int i, j, k;
scanf("%d%d", &An, &Bn);
N = 1;
while(N-1 <= An+Bn) N <<= 1;
for(i = 0; i <= An; i++) scanf("%lf", &A[i].r);
for(i = 0; i <= Bn; i++) scanf("%lf", &B[i].r);
fft(A, 1); fft(B, 1);
for(i = 0; i < N; i++) A[i] = A[i]*B[i];
fft(A, -1);
for(i = 0; i <= An+Bn; i++) printf("%d ", (int)(A[i].r+0.5));
printf("\n");
return 0;
}


## 完结撒花..

posted @ 2017-03-10 14:01  Ra1nbow  阅读(...)  评论(... 编辑 收藏