首页 > 其他分享 >FFT学习笔记

FFT学习笔记

时间:2022-08-15 21:56:39浏览次数:49  
标签:frac cdot 多项式 sum FFT 笔记 学习 cdots omega

-1. 前置知识

基础的复数知识。

0. 什么是多项式乘法

众所周知,多项式本质是一种特殊的函数,可以表示为自变量的若干次幂之和,即

\[F(x)=\sum_{i=0}c_i\cdot x^i \]

其中 \(c_i\) 被称为 \(x^i\) 的系数。

已知 \(F,G\) 是两个多项式函数,考虑定义一个新的函数 \(H(x)=F(x)G(x)\)。我们称 \(H\) 为 \(F\) 与 \(G\) 的乘积,记作 \(H=F\cdot G\),这种通过已知的两个多项式函数,以乘积形式生成另一个函数的运算,称为多项式乘法。

本文讨论的就是计算多项式乘法,即求出 \(H\) 解析式的过程。

根据定义,有

\[\begin{aligned} H(x)&=F(x)\cdot G(x) \\ &= (\sum_{i=0} a_ix^i)(\sum_{j=0}b_jx^j) \\ &= \sum_{i=0}(a_ix^i\cdot \sum_{j=0}b_jx^j)\\ &= \sum_{i=0}\sum_{j=0} a_ib_jx^{i+j}\\ \end{aligned} \]

我们发现,两个多项式的乘积仍是多项式。

由于多项式的特性,我们只需知道该多项式的每一项系数即可,即使用一个\(n\) 维向量去描述一个\(n-1\)次多项式。

显然,有一种暴力的方法 \(H[i]=\sum_{j=0}^iF[j]* G[i-j]\)。我们称这种形如\(c_n=\sum_{i\otimes j=n}a_i\cdot b_j\)运算为卷积。(其实只是卷积的一种,即加法卷积)

显然,这样做复杂度是 \(\mathcal{O}(n^2)\)的,实在太慢了。

\(\frac{1}{2}\). 多项式的点值表示法

既然多项式是一个函数,那我么可以画出它的函数图象,并在上面取几个点。

显然,由待定系数法可知,给定 \(n+1\) 个不同的点,可以唯一确定一个次数不超过 \(n\) 的多项式。

我们可以考虑用点值反推解析式,即如果已知函数 \(F,G\) 的\(n\) 个对应点,由定义式 \(H(x)=F(x)G(x)\) 可以直接 \(\mathcal O(n)\) 算出 \(H\) 对应的点!

但如果随意带入点,用待定系数法高斯消元暴力反推解析式,复杂度 \(\mathcal O(n^3)\),直接爆炸。

当然我们可以用拉格朗日插值公式做到 \(\mathcal O(n^2)\),然并卵,还不如暴力乘法……

\(\frac{9}{10}\). 单位根的引入及其性质

既然带一些普通的东西进去没啥用,那我们直接躺平带一些奇怪的东西进去,比如复数。考虑带入 \(n\) 次单位根。

定义:方程 \(x^n-1=0\) 在复数域的全部解称为 \(n\) 次单位根。由代数基本定理知,\(n\) 次单位根有 \(n\) 个,分别记为 \(\omega_n^0 \cdots \omega_n^{n-1}\) 。由复数乘法模相乘,幅角相加可知,\(n\) 次单位根全部在单位圆上,且幅角为 \(\frac{2\pi}{n}\) 的整数倍。显然,\(\omega_n^k=\cos \frac{2\pi k}{n}+i\sin \frac{2\pi k}{n}\)。当然,其它单位根是 \(\omega_n^1\) 的整数次幂,因此我们称其为主\(n\)次单位根,简记为 \(\omega_n\)。

根据欧拉公式 \(e^{ix}=\cos x+i\sin x\),有\(\omega_n^k=e^{\frac{2\pi k}{n}i}\)

一些性质:

  1. \(\omega_n^k=(\omega_n^1)^k\)

  2. \(\omega_n^i\cdot\omega_n^j=\omega_n^{i+j}\)

  3. \(\omega_{dn}^{dk}=\omega_n^k\) (消去引理)

证明:\(\large\omega_{dn}^{dk}=e^{\frac{2\pi dk}{dn}i}=e^{\frac{2\pi k}{n}i}=\omega_n^k\)

  1. 若 \(2|n\),\(\omega_n^{k+\frac{n}{2}}=-\omega_n^k\)

证明:\(\omega_n^{k+\frac{n}{2}}=\omega_n^{k}\omega_n^{\frac{n}{2}}=\omega^1_2\omega_n^k=-\omega_n^k\)

推论:\((\omega_n^k)^2=(\omega_n^{k+\frac{n}{2}})^2=\omega_{n/2}^k\)(折半引理)

  1. 对于任意不为 \(n\) 倍数的 \(k\),有:

\[\sum_{i=0}^{n-1}(\omega_n^k)^i=0 \quad \text{(求和引理)} \]

证明:

\[\sum_{i=0}^{n-1}(\omega_n^k)^i=\frac{(\omega_n^k)^n-1}{\omega_n^k-1}=\frac{1^k-1}{\omega_{n}^k-1}=0 \]

知道这些结论,就可以学习 FFT 了!

1. DFT

考虑求出一个多项式在 \(\omega_n^0\cdots\omega_n^{n-1}\) 的值,我们称这种运算为 DFT。

首先,我们将该多项式的次数扩充为 \(2\) 的幂次。

将 \(A(x)\) 写成系数向量形式 \(a=[a_0,a_1,a_2,\cdots,a_{n-1}]\)

考虑将其按奇数项和偶数项分为

\[a^{[0]}=[a_0,a_2,a_4,\cdots,a_{n-2}] \]

\[a^{[1]}=[a_1,a_3,a_5,\cdots,a_{n-1}] \]

分别记其对应的多项式为 \(A^{[0]}(x)\) 和 \(A^{[1]}(x)\)。

\[A(x)=a_0+a_1x+a_2x^2+\cdots+a_{n-1}x^{n-1} \]

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

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

换成 \(x^2\):

\[A(x)=a_0+a_1x+a_2x^2+\cdots+a_{n-1}x^{n-1} \]

\[A^{[0]}(x^2)=a_0+a_2x^2+a_4x^4+\cdots+a_{n-2}x^{n-2} \]

\[A^{[1]}(x^2)=a_1+a_3x^2+a_5x^4+\cdots+a_{n-1}x^{n-2} \]

\(A^{[1]}\) 乘上 \(x\):

\[A(x)=a_0+a_1x+a_2x^2+\cdots+a_{n-1}x^{n-1} \]

\[A^{[0]}(x^2)=a_0+a_2x^2+a_4x^4+\cdots+a_{n-2}x^{n-2} \]

\[xA^{[1]}(x^2)=a_1x+a_3x^3+a_5x^5+\cdots+a_{n-1}x^{n-1} \]

由此得到 \(A(x)=A^{[0]}(x^2)+xA^{[1]}(x^2)\)

说人话:如果知道 \(A^{[0]}(x^2)\) 和 \(A^{[1]}(x^2)\) ,就可以求出 \(A(x)\)。

问题转化为求出所有 \(A^{[0]}((\omega_n^k)^2)\) 和 \(A^{[1]}((\omega_n^k)^2)\)。

由折半引理知,求出所有 \(A^{[0]}((\omega_n^k)^2)\) 和 \(A^{[1]}((\omega_n^k)^2)\) 就是求出所有 \(A^{[0]}(\omega_{n/2}^k)\) 和 \(A^{[1]}(\omega_{n/2}^k)\)。

注意到这是原问题的子问题,可以直接分治处理。递归边界为 \(n=1\),此时对应的多项式为常量,直接返回即可。

考虑合并。

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

当然,如果 \(k\geq \frac{n}{2}\),我们让 \(k\bmod \frac{n}{2}\) 即可。

显然,合并复杂度为 \(\mathcal O(n)\)。

则 \(FFT\) 的复杂度 \(T(n)=2T(\frac{n}{2})+\mathcal O(n)\),故 \(T(n)=\mathcal O(n\log n)\)

全剧终。

我们发现,我们得到的只是一些点值罢了,根本不是系数表示。

2. IDFT

即 DFT 的逆运算。

如何将点值转换为系数?注意到上文提到 DFT 时,说的是将 \(A(x)\) 看成系数向量,即一维矩阵。

如果我们将得到的点值记为一维向量 \(b\),则有:\(a=IDFT(b)\),即 \(b=DFT(a)\)。

我们现在知道 \(b\) 求 \(a\)。

由定义知 \(b[k]=\sum_{i=0}^{n-1}(\omega_n^k)^i\cdot a[i]\)

如果你会单位根反演,你可能能够推出 \(a\) 关于 \(b\) 的式子。这里给出单位根反演的原理。

还记得求和引理吗?

对于任意不为 \(n\) 倍数的 \(k\),有:

\[\sum_{i=0}^{n-1}(\omega_n^k)^i=0 \]

此外,若 \(k\) 为 \(n\) 的倍数,显然上个式子的值为 \(n\)。

\[\frac{1}{n}\sum_{i=0}^{n-1}(\omega_n^k)^i=[n|k] \]

考虑如下式子:

\[a_k=\sum_{i=0}^{n-1} a_i\cdot[k=i] \]

感觉是句废话……

考虑把 \([k=i]\) 往 \([n|k]\) 上靠。

我们发现,\(i,k\) 的取值范围是 \([0,n-1]\),即 \(k=k\bmod n,i=i\bmod n\)

那么,\([k=i]\) 可以写成 $[k\bmod n=i\bmod n] $,即 $ [k\equiv i\pmod n]$ 或 \([n|(k-i)]\)

带回原式试试!

\[\begin{aligned} a_k &=\sum_{i=0}^{n-1}a_i\cdot[n|(k-i)]\\ &=\sum_{i=0}^{n-1}a_i\cdot\frac{1}{n}\sum_{j=0}^{n-1}(\omega_n^{k-i})^j \\ &=\frac{1}{n}\sum_{i=0}^{n-1}\sum_{j=0}^{n-1} \omega_n^{ik}\omega_n^{-ij}a_i \\ &=\frac{1}{n}\sum_{i=0}^{n-1}(\omega_n^{ki})\cdot \sum_{j=0}^{n-1}\omega_n^{-ij}a_i \end{aligned} \]

完蛋了……推不出来了

但 \([n|(k-i)]\) 与 \([n|(i-k)]\) 是等价的。

再试一次。

\[\begin{aligned} a_k &=\sum_{i=0}^{n-1}a_i\cdot[n|(i-k)]\\ &=\sum_{i=0}^{n-1}a_i\cdot\frac{1}{n}\sum_{j=0}^{n-1}(\omega_n^{i-k})^j \\ &=\frac{1}{n}\sum_{i=0}^{n-1}\sum_{j=0}^{n-1} \omega_n^{-ik}\omega_n^{ij}a_i \\ &=\frac{1}{n}\sum_{i=0}^{n-1}(\omega_n^{-ki})\cdot \sum_{j=0}^{n-1}\omega_n^{ij}a_i \\ &=\frac{1}{n}\sum_{i=0}^{n-1}(\omega_n^{-k})^i\cdot b_i \end{aligned} \]

这次成功了。变形一下:

\[n\cdot a[k]=\sum_{i=0}^{n-1}(\omega_n^{-k})^ib[i] \]

我们发现右边和 DFT 过程几乎完全一致,只是将 \(\omega_n^k\) 改为 \(\omega_n^{-k}\)。

到这里,FFT 的理论部分就讲完了。

3. 代码实现

标签:frac,cdot,多项式,sum,FFT,笔记,学习,cdots,omega
From: https://www.cnblogs.com/pref-ctrl27/p/16589789.html

相关文章

  • 长春大学 第二组刘禹彤 学习笔记
    打卡32天###学习内容Mysql数据库数据库数据库【按照数据结构来组织,存储和管理数据的仓库】,是一个长期存储在计算机内的,有组织的,可共享的,统一管理大量的数据的集合数......
  • 2022-8-15MySQL的学习
    MySQL数据库数据库数据库【按照数据结构来组织来存储和管理数据的仓库】。是一个长期存储在计算机内的有组织的可共享的,统一管理的大量数据的集合。数据对于公司来......
  • 2022-8-15 第一组 (≥▽≤) 学习笔记
    目录1.MySQL数据库1.1数据库1.2MySQL1.3基本操作1.4表2.SQL语言2.1SQL语句的分类2.1.1DCL(数据控制语言)给用户授权撤销授权查看指定用户授权删除用户2.1.2DDL(数据......
  • 笔记day03
    上周内容回顾循环结构之while循环1.while语法结构while条件:  条件成立之后执行的循环体代码2.while+break3.while+continue4.while+else5.死循环、循环......
  • Linux 0.11学习
    学习链接:sunym1993/flash-linux0.11-talk:你管这破玩意叫操作系统源码—像小说一样品读Linux0.11核心代码(github.com)1.从开机到运行main.c的过程在主板上写死......
  • 操作系统学习笔记2 | 操作系统接口
    这部分将讲解上层应用软件如何与操作系统交互,理解操作系统到底发生了什么事情,理解操作系统工作原理,为以后扩充操作系统、设计操作系统铺垫。参考资料:课程:哈工大操作系......
  • 2022-08-15 第四小组 王星苹 学习笔记
    学习心得       开始MySQL数据库的学习,创建库,再创建表,在表中保存多条数据,一列就是一个字段,也可以增删改查心情 又新学个新东西,新起点,出发加油掌握情况:背......
  • RecyclerView 的学习记录
    官方文档**RecyclerView样式与适配器等解耦**:通过设置不同的LayoutManager,就可以实现不同的布局展示样式;通过设置不同的ItemDecoration,可以实现不同的间......
  • Hadoop学习第一天
    学习课程是B站上的黑马程序员第一阶段主要是基础的概念,数据、大数据;大数据特点;数据分析的基本流程、方向;分布式、集群;操作系统,虚拟机。基本上就是这些基本概念的学习。第......
  • 三更 SpringSecurity笔记
    SpringSecurity从入门到精通课程介绍0.简介​SpringSecurity是Spring家族中的一个安全管理框架。相比与另外一个安全框架Shiro,它提供了更丰富的功能,社区资源......