Learning：多项式（一）（FFT）

最新推荐文章于 2023-05-23 16:13:53 发布

ezoiHQM

最新推荐文章于 2023-05-23 16:13:53 发布

阅读量1.1k

点赞数 1

分类专栏： FFT Learning 文章标签： FFT

本文链接：https://blog.csdn.net/ezoiHQM/article/details/81191697

版权

Learning 同时被 2 个专栏收录

8 篇文章 0 订阅

订阅专栏

FFT

7 篇文章 0 订阅

订阅专栏

多项式的表示

1 系数表示法

设 $A(x)$ 是一个关于 $x$ 的 $n$ 次多项式，则： $A(x)=\sum_{j=0}^na_jx^j$

2 点值表达法

我们可以把 $n$ 次多项式 $A(x)$ 看作一个函数，那么它可以用平面直角坐标系上的 $n$ 个点 $(x_1,y_1),(x_2,y_2),...,(x_n,y_n)$ 来确定。我们把这 $n$ 个点代入 $A(x)$ ，我们就可以得到一个 $n$ 元一次方程组，然后通过高斯消元就可以确定这个多项式。

系数表达法的多项式乘法时间复杂度显然是 $O(n^2)$ 的，但是点值表达法的多项式的乘法的时间复杂度却是 $O(n)$ 的（两个多项式的每一个点的横坐标都相等）。那我们就会希望可以利用点值表达法的这一优势来快速地进行多项式乘法。但是我们平时使用的都是系数表示法，想要利用这个优势，就要快速地将系数表达法转化为点值表达法。

复数

复数的定义

我们把形如 $a+bi$ （ $a$ , $b$ 均为实数）的数称为复数，其中 $a$ 称为实部， $b$ 称为虚部， $i$ 称为虚数单位。当虚部等于零时，这个复数可以视为实数；当 $z$ 的虚部不等于零时，实部等于零时，常称 $z$ 为纯虚数。复数域是实数域的代数闭包，也即任何复系数多项式在复数域中总有根。复数是由意大利米兰学者卡当在十六世纪首次引入，经过达朗贝尔、棣莫弗、欧拉、高斯等人的工作，此概念逐渐为数学家所接受。（摘自百度百科）

复数的运算法则

1 复数的加法

设 $z_1=a+bi,z_2=c+di$ ，那么： $z_1\pm z_2=(a+c)+(b+d)i$

2 复数的乘法

设 $z_1=a+bi,z_2=c+di$ ，那么： $z_1\times z_2=(ac-bd)+(ad+bc)i$

3 (一个奇奇怪怪的东西)

$e^{i\theta}=cos\theta+isin \theta$

单位根

定义

如果 $\omega^n=1$ ，那么我们称 $\omega$ 为 $n$ 次单位根。

单位根的值

事实上，在这里： $\omega_k=e^\frac{2\pi ik}{n}=cos\frac{2\pi k}{n}+isin\frac{2\pi k}{n}$
由欧拉恒等式 $(e^{i\pi}+1=0)$ 得到。

消去引理

$\omega_{dn}^{dk}=\omega_n^k$
证明：

$\omega_{dn}^{dk}=(e^\frac{2\pi i}{dn})^{dk}=e^\frac{2\pi ik}{n}=\omega_n^k$
折半引理

如果 $n$ 是大于 $0$ 的一个偶数，那么前 $\frac{n}{2}$ 个 $n$ 次单位根的平方的集合等于后 $\frac{n}{2}$ 个 $n$ 次单位根的平方的集合。
证明：对于任意的 $k<\frac{n}{2}$ ，有：

$(\omega_n^{k+\frac{n}{2}})^2=\omega_n^{2k+n}=\omega_n^{2k}\omega_n^n=\omega_n^{2k}=(\omega_n^k)^2$
因此， $(\omega_n^{k+\frac{n}{2}})^2=(\omega_n^k)^2$

DFT 离散傅里叶变换

终于切入正题了。
DFT就是将一个多项式的系数表达法转化为点值表达法。
具体操作：
对于 $n-1(2|n)$ 次多项式 $A(x)$ ，我们有：
$A(x)=a_0+a_1x+a_2x^2+...+a_{n-1}x^{n-1}$
令 $A^{[1]}=a_0+a_2x^2+...+a_{n-2}x^{\frac{n}{2}-1}$
$A^{[2]}=a_1+a_3x+...+a_{n-1}x^{\frac{n}{2}-1}$
则 $A(x)=A^{[1]}(x^2)+xA^{[2]}(x^2)$
把 $\omega_n^k(k<\frac{n}{2})$ 代入可得到： $A(\omega_n^k)=A^{[1]}((\omega_n^k)^2)+\omega_n^kA^{[2]}((\omega_n^k)^2)$
由折半引理我们会发现 $(\omega_n^k)^2=(\omega_n^{k+\frac{n}{2}})^2$
易证 $\omega_n^k=-\omega_n^{k+\frac{n}{2}}$
所以 $A(\omega_n^{k+\frac{n}{2}})=A^{[1]}((\omega_n^k)^2)-\omega_n^kA^{[2]}((\omega_n^k)^2)$
我们会发现这两个式子只有常数项不同
那么当我们在求第一个式子的时候，我们可以直接在 $O(1)$ 的时间内同时得到第二个式子的值
然后就这样递归下去，就可以得到 $A(x)$ 的点值表达式了

IDFT 离散傅里叶逆变换

我们用点值表达法求出了两个多项式的乘积，但是我们还需要将这个点值表达法转换回系数表达法。那这个又该怎么操作呢？
首先，我们可以吧从系数表达法转换为点值表达法的过程用矩阵来表示：
$\begin{bmatrix}ans_0\\ans_1\\ans_2\\ans_3\\...\\ans_{n-1}\end{bmatrix}\quad=\begin{bmatrix}\omega_n^0&\omega_n^0&\omega_n^0&...&\omega_n^0\\\omega_n^0&\omega_n^1&\omega_n^2&...&\omega_n^{n-1}\\\omega_n^0&\omega_n^2&\omega_n^4&...&\omega_n^{2(n-1)}\\...&...&...&...&...\\\omega_n^0&\omega_n^{n-1}&\omega_n^{2(n-1)}&...&\omega_n^{(n-1)^2}\end{bmatrix}\quad\begin{bmatrix}a_0\\a_1\\a_2\\a_3\\...\\a_{n-1}\end{bmatrix}\quad$
事实上，我们会发现，把中间那个矩阵的所有元素全部取一个倒数，再乘上一个 $\frac{1}{n}$ ，就可以得到原本的矩阵。
这样，我们就可以把点值表达法转换回系数表达法。

递归版FFT

从上面我们可以很容易得到FFT的递归写法。
代码实现：

void init(){
    int limit=1;
    while(limit<=n+m)
        limit<<=1;
}
void FFT(int limit,cmplx *a,int type){
    if(limit==1)
        return;
    cmplx a1[limit>>1],a2[limit>>1];
    for(int i=0;i<=limit;i+=2){
        a1[i>>1]=a[i];
        a2[i>>1]=a[i+1];
    }
    FFT(limit>>1,a1,type);
    FFT(limit>>1,a2,type);
    const double x=2.0*pai/(const int)limit;
    cmplx wn=cmplx(cos(x),type*sin(x)),w=cmplx(1,0);
    for(int i=0;i<(limit>>1);i++,w*=wn){
        const cmplx y=w*a2[i];
        a[i]=a1[i]+y;
        a[i+(limit>>1)]=a1[i]-y;
    }
    return;
}

但是递归版的FFT的常数非常大，一不小心就会被卡飞掉。所以我们需要一个更加高效的做法。

迭代版FFT

在递归版FFT中，我们要将奇数项和偶数项分离。
以 $n=8$ 时为例，我们可以观察递归版的FFT每一次进行分离之后的序列：
0 1 2 3 4 5 6 7

0 2 4 6 1 3 5 7

0 4 2 6 1 5 3 7

0 4 2 6 1 5 3 7
再把它们转化为二进制看看：
000 100  010  110  001 101 011  111
再把它们的二进制反转后：
000 001 010 011 100 101 110 111
神奇的事情发生了！得到的这个序列的二进制反转后与原序列是一样的。
所以我们就可以通过这个很轻松地求出这个序列。
然后我们就可以迭代实现FFT了！
代码如下：

void init(){
    int limit=1,l=0;
    while(limit<=n+m){
        limit<<=1;
        l++;
    }
    for(int i=0;i<limit;i++)
        r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
    return;
}
void FFT(int limit,cmplx *a,int type){
    for(int i=0;i<limit;i++)
        if(i<r[i])
            swap(a[i],a[r[i]]);
    for(int mid=1;mid<limit;mid<<=1){
        const double x=pai/(const int)mid;
        cmplx wn=cmplx(cos(x),type*sin(x));
        for(int R=mid<<1,j=0;j<limit;j+=R){
            cmplx w=cmplx(1,0);
            for(int k=0;k<mid;k++,w*=wn){
                cmplx y=a[j+k],z=w*a[j+k+mid];
                a[j+k]=y+z;
                a[j+k+mid]=y-z;
            }
        }
    }
    return;
}

给大家几道题：
【BZOJ3160】万径人踪灭
【BZOJ3771】Triple
【BZOJ4827】礼物

如果有误在评论区吼一声哦！

ezoiHQM

关注

1
点赞
踩
4

收藏

觉得还不错? 一键收藏
0
评论
Learning：多项式（一）（FFT）

多项式的表示1 系数表示法设是一个关于的次多项式，则：2 点值表达法我们可以把次多项式看作一个函数，那么它可以用平面直角坐标系上的个点来确定。我们把这个点代入，我们就可以得到一个元一次方程组，然后通过高斯消元就可以确定这个多项式。系数表达法的多项式乘法时间复杂度显然是的，但是点值表达法的多项式的乘法的时间复杂度却是的（两个多项式的每一个点的横坐标都相等）。那我们就会希望可以...
复制链接

扫一扫