快速傅里叶算法（灌水）

定义和声明

多项式表示方法

1. 系数表示法： $A(x) = \sum_{j=0}^{n-1} a_jx^j$
2. 点值表示法： ${ (x_0,y_0), (x_1,y_1), ..., (x_{n-1},y_{n-1}) }$
   
   • 可证这n个不同点唯一确定一个次数界为n的多项式
   • 也可以用公式$A(x) = \sum_{k=0}^{n-1}y_k*\prod_{j\neq k} (x-x_j)/\prod_{j\neq k}(x_k-x_j)$在$O(n^2)$的复杂度下求解系数表达式

单位复数根

1. n次单位复根，满足条件$w^n=1$，可得$w=e^{2\pi i k/n}$
2. n次单位复数根满足一些基本性质：
   
   • 消去性质：$w^{dk}_{dn}=w^k_n$,由此可得$对于偶数n，w^{n/2}_n=w_2=-1$
   • 折半性质： $w^{k+n/2}_n=-w^k_n$
   • 求和性质：$\sum_{j=0}^{n-1}{(w^k_n)}^j=0$

对于多项式A(x),计算它在n个n次单位复数根处的值，即可求出其离散傅里叶变换(DFT)。所以我们求其在n个n次单位复数根处的点值表达式，也相当于对其进行一次离散傅里叶变换。而使用其离散傅里叶变换得到的点值表达式能很方便地在O(n)的时间内计算得到多项式乘法结果。
所以整个算法的过程即： FFT->点值乘法->逆FFT

算法描述

对于多项式$A(x)$，定义

$$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)$$
所以，求A(x)在 $w^0_n,w^1_n,...,w^{n-1}_n$的值的问题被转化为求次数界为n/2的 $A^{[1]}(x)$和 $A^{[0]}(x)$在点 ${(w^0_n)}^2,{(w^1_n)}^2,...{(w^{n-1}_n)}^2$，这里实际上仅有n/2个不同值，此时问题规模被缩小了一半。

而计算逆FFT时，只需根据

$$a_j=1/n\sum_{k=0}^{n-1}y_kw^{-kj}_n$$ 即可
这里使用一段引自别人博客的代码（#

const double eps(1e-8);
typedef long long lint;

const double PI = acos(-1.0);

struct Complex
{
    double real, image;
    Complex(double _real, double _image)
    {
        real = _real;
        image = _image;
    }
    Complex(){}
};

Complex operator + (const Complex &c1, const Complex &c2)
{
    return Complex(c1.real + c2.real, c1.image + c2.image);
}

Complex operator - (const Complex &c1, const Complex &c2)
{
    return Complex(c1.real - c2.real, c1.image - c2.image);
}

Complex operator * (const Complex &c1, const Complex &c2)
{
    return Complex(c1.real*c2.real - c1.image*c2.image, c1.real*c2.image + c1.image*c2.real);
}

int rev(int id, int len)
{
    int ret = 0;
    for(int i = 0; (1 << i) < len; i++)
    {
        ret <<= 1;
        if(id & (1 << i)) ret |= 1;
    }
    return ret;
}

//当DFT= 1时是DFT, DFT = -1则是逆DFT
Complex* IterativeFFT(Complex* a, int len, int DFT)//对长度为len(2的幂)的数组进行DFT变换
{
    Complex* A = new Complex[len];//用A数组存储数组a分组之后新的顺序
    for(int i = 0; i < len; i++)
        A[rev(i, len)] = a[i];
    for(int s = 1; (1 << s) <= len; s++)
    {
        int m = (1 << s);
        Complex wm = Complex(cos(DFT*2*PI/m), sin(DFT*2*PI/m));
        for(int k = 0; k < len; k += m)//这一层结点的包含数组元素个数都是(1 << s)
        {
            Complex w = Complex(1, 0);
            for(int j = 0; j < (m >> 1); j++)//折半引理, 根据两个子节点计算父亲节点
            {
                Complex t = w*A[k + j + (m >> 1)];
                Complex u = A[k + j];
                A[k + j] = u + t;
                A[k + j + (m >> 1)] = u - t;
                w = w*wm;
            }
        }
    }
    if(DFT == -1) for(int i = 0; i < len; i++) A[i].real /= len, A[i].image /= len;
    return A;
}