大数乘法(快速傅立叶变换)上
博客搬家:最爱午后红茶
在谈论大数乘法前,先来看看多项式乘法,假设
A = 7x^3 + 5x^2 + 3x + 4
B = 4x^2 + 6
C = A * B
那传统的做法就是这样:
7x^3 + 5x^2 + 3x + 4
4x^2 + 6
---------------------------------------
+42x^3 + 30x^2 + 18x + 24
28x^5 + 20x^4 +12x^3 + 16x^2
-----------------------------------------------------------
28x^5 + 20x^4 + 54x^3 + 46x^2 + 18x + 24
若A 和 B 都有 n 个系数,很明显这样做的时间复杂度为O(n^2),当n很大时效果相当不理想
像 B = 4x^2 + 6 这样的式子叫系数表示法;如果我们把 B 一般化,即 B = b0 + b1x + b2x^2;则我们至少需要 3 对 (x, y) 值才能求出 b0, b1, b2。这里 b0 = 6, b1 = 0, b2 = 4;例如 B 在x0 = 1, x1 = 2, x2 = 3 的情况下的解分别为 y0 = 10, y1 = 22, y2 = 42。那我们也可以用 {(1, 10),(2, 22),(3, 42)},即 {(x0, y0),(x1, y1),(x2, y2)} 来表示唯一的多项式B,即假如我们在不知道b0, b1, b2的情况下可以用{(1, 10),(2, 22),(3, 42)} 来唯一表示多项式B,即唯一确定b0, b1, b2。当然 x0, x1, x2 也可以取其他值,但对于这个例子一定要取够3个,因为有3个未知数 b0, b1, b2;对应的 y0, y1, y2 也会不同,这就叫点值表示法。
为什么要这样做?还是上面那个例子,我们用点值表示法算一遍C = A * B,(A 和 B 的xi要对应相同)还是假设取x0 = 1, x1 = 2, x2 = 3,我们发现 A 有4个未知数,所以要多一个x3,假设x3 = 4;为了方便计算,我们把 B 的 (x3, y3) 也算出来。
A = {(1, 19),(2, 86),(3, 247),(4, 544)}
B = {(1, 10),(2, 22),(3, 42),(4, 70)}
xi 不变,yi对应相乘得到 Ci = (xi, yai * ybi)
C = {(1, 190),(2, 1892),(3, 10374),(4, 38080)}
到这里我们发现对于一般情况的等长度的两个多项式相乘,结果的长度会接近加倍,比如
A = a0 + a1x + a2x^2 + a3x^3
B = b0 + b1x + b2x^2 + b3x^3
C = A * B = c0 + c1x + c2x^2 + c3x^3 + c4x^4 + c5x^5 + c6x^6
因此有n个系数的两个多项式相乘(比如A * B)时仅用其n个点值相乘是不能正确表示结果的(比如C),所以我们要使用足够多的点值,比如这里 A 和 B 都应该至少使用7个点值。如下(数有点大,偷下懒。。。)
A = {(x0, y0),(x1, y1),(x2, y2),(x3, y3),(x4, y4),(x5, y5),(x6, y6)}
B = {(x0, y0'),(x1, y1'),(x2, y2'),(x3, y3'),(x4, y4'),(x5, y5'),(x6, y6')}
C = {(x0, y0*y0'),(x1, y1*y1'),(x2, y2*y2'),(x3, y3*y3'),(x4, y4*y4'),(x5, y5*y5'),(x6, y6*y6')}
这样C的点值才能表示唯一的6次多项式(7个未知系数有6次幂),最后,用某方法(重点下面讲)把C的点值表示法转换成系数表示法就达到了目的,这个过程叫插值
所以,先理下思路
1)分别计算A、B足够多的点值(也叫求值),使其能唯一表示C(如A、B的系数都有n个,则通常扩展到2n个点值)
2)点乘:即A、B点值对应相乘得到C
3)插值:把C的点值表示法转换为系数表示法
其中点乘的复杂度显然为O(n),我们记N = 2n(记住了下面要用的~),所以我们如能找到一种方法能够快速对A、B进行求值,插值,那目的就达到了,事实上也的确有这样的方法,这也是多项式乘法的重点,不然都没法往下讲。。。
----------------------------------------------------------------------------------------------------------------------------------------------------
先看求值,y = a0 + a1x + a2x^2 + ... + anx^n 给出一个 x 要求 y 如果按正常顺序求需要做 (1 + n) * n / 2 次乘法和 n 次加法;这种方法肯定是不可取的,因为算一个值都要 O(n^2),那算 N 个值就是 O(n^3) 了,比系数表示法直接乘还慢~
再来看另一种方法,国产方法,没错就是它,秦九韶算法(在国外也叫霍纳(Horner)法则);为了理解方便,我们举个短的例子 y = a0 + a1x + a2x^2 + a3x^3;那它就可以化为 y = a0 + x(a1 + x(a2 + a3x)),这样对某个x求y就需要做 n 次乘法和 n 次加法,算 N 个值的复杂度为 O(n^2)。比上面的好了点,不过跟系数表示法相比并没有什么优势~
所以我们需要第三种方法,这种方法的特点其实是在 xi 的选取上很特殊;它选的 xi 不是从 x0 = 1, x1 = 2, ... , xn = n 或者其他一些相异的整数或实数;它选的是复数!当然也不是乱选的,因此我们要先在这里停一下,吸口气,去打下复数的基础。。。
---------------------------------------------------------------------------------------------------------------------------------------------------------
z = x + yi 这是常见的复数表示形式大家都知道,|z| = sqrt(x^2 + y^2) 表示复数的模
用建立了笛卡尔直角坐标系的平面表示复数的平面叫做复平面
下面必须盗图了(来自《复变函数与积分变换》)不难的,请坚持看完!
复数的指数式也称为欧拉公式,它更像是一种定义。欧拉公式是非常美妙的式子,为什么要引入它呢?主要是为了简化计算!因为 z = x + yi 可以表示成 z = r*(cosθ + i*sinθ),r = |z|。看个例子:
所以用指数式表示复数是很简洁的,我们最后的求值与插值的方法也是用复数的指数式表示。
这里有一点要提一下,有人说 θ = 0 时,e^(θi) = e^(0*i) = 1;结果没有错,但不是这样算的。e^(θi) 不能理解为 e^(θ*i);那个不是乘号,只是一种记号。正确解法应该把它展开算 e^(θi) = cosθ + isinθ = cos0 + isin0 = 1
设 z 为复数,n 为正整数,若存在复数 w 满足方程 w^n = z;
重要的事情说三遍:
在几何上 z^(1/n) 的 n 个值就是以原点为中心,r^(1/n) 为半径的圆的内接正 n 边形的 n 个顶点。
在几何上 z^(1/n) 的 n 个值就是以原点为中心,r^(1/n) 为半径的圆的内接正 n 边形的 n 个顶点。
在几何上 z^(1/n) 的 n 个值就是以原点为中心,r^(1/n) 为半径的圆的内接正 n 边形的 n 个顶点。
还有上面的公式 z = e^(iθ) != 0,wk = (r^(1/n)) * e^(i(θ+2*k*pi) / n)。θ 为 z 在复平面中向量表示的幅角。
做个题目理解下:
再举个幅角不是0的例子,根据公式不难算的:
这个就不画了,想加深理解可以自己画一下~
--------------------------------------------------------------------------------------------------------------------------------------------------------
接下来才是高能部分:
我们前面讲了,在对 A = a0 + a1x + a2x^2 + ... + a(N-1)x^(N-1) (我们假设A的系数有 N(前面讲过N = 2n) 个,即a0 ~ a(N-1))求值时我们选的 N 个 xi (0 <= i < N) 用的是复数,那现在告诉你这 N 个 xi 就是 w^N = 1 的 N 个解,w0 ~ wN 叫N 次单位复根。即 x0 ~ xN 分别等于 w0 ~ wN;
这里的 N 也就是一种记号,如果看不顺,理解成 n 也是没问题的~
再上一遍公式
对于 N 次单位复根的求解,θ 值为 0 ,r^(1/N) = 1。其中解 w1 = e^(i(2*pi/N)) 叫主单位复根;其实也就是把一个圆分成 N 份,w1 表示每份的大小,只不过是在复平面上表示。
下面我们要换下记号了
记 N 为 n
则 A = a0 + a1x + a2x^2 + ... + a(n-1)x^(n-1)
也可记为
原来的各个单位复根 w0, w1, ... , w(n-1) 记为 ,下标的 n 就是对单位复数 1 开 n 次方的意思,很简单大家都能明白。这样是因为下面的计算会改变下标。
那用 代入 xk 来求解 A(xk) 就可以表示为 (0 <= k <= n-1)
如果把 A 的系数用向量表示 a = {a0, a1, a2, ..., a(n-1)},把求出的所有 y 值也用向量表示 y = {y0, y1, y2, ... , y(n-1)},那我们可以称向量 y 为 向量 a 的离散傅立叶变换(Discrete Fourier Transform, DFT)。我们的目的之一就是把向量 a 的离散傅立叶变换(即向量 y )快速地算出来,我们使用的是一种叫快速傅立叶变换(FFT)的算法。这里我们可以猜到,FFT 算法起码要比系数表示法直接相乘的 O(n^2) 的复杂度更低才行,幸运的是,FFT 的时间复杂度是 O(n*logn);下面讲的是 FFT 是怎样快速地计算向量 y 的。
在讲FFT之前需要讲清楚一件事,就是上面所讲的 n 到底是多大?
举个例子吧,比如:
A = 3 + 2X^2 + 4X^3 + 6x^4
B = 2x + X^2 + 5X^3
以 A 举例(B类似),它可以写成标准式:A = a0 + a1x + a2x^2 + a3x^3 + a4x^4;这里有 5 个系数,根据前面讲的为了能正确表示 A 与 B相乘的结果 C,我们把 A 放大到 A = a0 + a1x + ... + a9x^9,使它变成10个系数(当然 a5 ~ a9 的值都是 0,B 类似),同时用 A 和 B 的 10 个点值计算C。不过对于 FFT 算法,这样扩大是不行的;根据算法要求(下面会讲清楚,现在不用管),我们要确保 n 为 2 的幂,所以 A 有 5 个系数,不是 2 的幂,我们应该把它先扩大到比5大且离 5 最近的那个 2 的幂,即 2^3 = 8。然后再把 8 扩大两倍变成 16;这才是 A B 两个多项式相乘的最终的 n。当然原来系数个数就是 2 的幂的话就只需要把它扩大两倍就行了。当系数个数为 513 的话,没错,要先把它扩大成1024个系数,然后再放大两倍变成 2048个系数;就是这么狠!想想这样有点耗空间~
所以最终的 A = a0 + a1x + a2x^2 + a3x^3 + ... + a(n-1)x^(n-1) 中的系数个数,即 n 必须是 2 的幂,我们总能够通过添加 ai = 0 满足需要。
FFT 算法运用了分治策略,这下有点理解为什么 n 要是 2 的幂了吧~,什么是分治策略就不展开讲了,不懂的就搞清楚了再往下看。
次数界就是多项式最高次幂是多少次方。
莫慌,我们看看 A(k) = a0 + a1k + a2k^2 + a3k^3 + ... + a(n-1)k^(n-1) 怎样算
把 k^2 代入 和 有
不难看出 A(k) = + k *
这样,求次数界为 n 的多项式 A(x) 在 处的值就转化为:
求次数界为 n / 2 的多项式 和 在 处的值。
形象地举个 n = 4 的例子就是:
在这里我们引入定理:
相消引理:对任何整数 n >= 0,k >= 0,d > 0;
也不难证明:
可以得出推论:对于任意偶数 n > 0, 有 (举个例子上下标不断除以2嘛)
结合上面的 n 个点 ,其实只有 n/2 个不同的值。
上面 n = 4 的例子中
可以知道 ,因此
所以, 前面 n/2 项跟后面 n/2 项是重复出现的;如上面 和 重复出现了2次。
因此,之前的式子可以改为:
这样,求次数界为 n 的多项式 A(x) 在 处的值就完全转化为:
求次数界为 n / 2 的两个多项式 和 在 处的值;规模神奇地缩小了一半,之后我们递归处理和 就可以把向量 y 整个儿求出来!伪代码应运而生:
如果伪代码还需要解释的话就说明之前的知识没有理解透彻,请往回多看几次(主要是高能部分)!
这就是多项式系数表示法转换为点值表示法(也称为求值过程)的 O(n*logn) 复杂度所用到的FFT算法的整个过程。
所以我们可以将多项式 A 和 B 按照上面的方法分别转化为点值表示,然后 y 值对应相乘就得到了多项式 C 的点值表示,接下来把 C 的点值表示转换为系数表示就完成了多项式乘法的运算。所以,下面讲介绍插值过程的算法。
插值过程的算法: 大数乘法(快速傅立叶变换)下
大数乘法(快速傅立叶变换)上相关推荐
- 大数乘法(快速傅立叶变换)下
博客搬家:最爱午后红茶 大数乘法(快速傅立叶变换)上 上篇已经已经讲了多项式乘法由系数表示法转化为点值表示法(即求值)的FFT算法的过程:接下来讲插值算法,它需不需要用新的算法写一遍呢?并不用这么麻烦 ...
- 一步一步的无障碍理解快速傅立叶变换
/// 作者:tt2767 声明:本文遵循以下协议自由转载-非商用-非衍生-保持署名|Creative Commons BY-NC-ND 3.0 查看本文更新与讨论请点击:http://blog.cs ...
- 快速傅立叶变换_FFT
快速傅立叶变换FFT 下面的排版比较乱,可以直接下word版的:http://minus.com/dbdrNFtIEL89Eb.docx: 原理:利用多项式点值求积的便捷性. 对于一个次数界为n的多项 ...
- 快速傅立叶变换的C语言实现方法
转自:http://www.beamsky.com/fft-c-language/ 傅立叶变换的重要性不用我说,想必大家也很清楚,有了傅立叶变换,我们就可以从信号的频域特征去分析信号.尤其在无线通信系 ...
- 分治算法--快速傅立叶变换
快速傅立叶变换 FFT算法 快速傅立叶变换 FFT算法 快速傅立叶变换 FFT算法 卷积和多项式相乘介绍 卷积和多项式相乘计算 蛮力算法 n阶多项式相乘运算 1.多项式的值表示法 2.利用值表达式法计 ...
- 线性代数28——复矩阵和快速傅立叶变换
原文 | https://mp.weixin.qq.com/s/YzPoPnRb-gEm_EiV9et0TA 实矩阵也可能碰到复特征值,因此无可避免地在矩阵运算中碰到复数. 矩阵当然也有可能包含复数, ...
- 快速傅立叶变换(FFT)
FFT 作用:快速求两个多项式的乘积/卷积 文章目录 FFT 前置知识 复数(Complex) 单位根 离散傅立叶变换(Discrete Fourier Transform , DFT) 快速傅立叶变 ...
- 快速傅立叶变换算法详解(小白也能看懂的FFT教程)
写在前面 我为什么要写这篇博客? \quad 如果你随便拉几个OI党,问他们最难理解的几个算法,FFT一定榜上有名.我从开始尝试理解FFT算法到真正理解FFT算法之间间隔了一年之久,但是当我真正理解了 ...
- 快速傅立叶变换的意义及应用
快速傅立叶变换的意义及应用 1.为什么要进行傅里叶变换,其物理意义是什么? 傅立叶变换是数字信号处理领域一种很重要的算法.要知道傅立叶变换算法的意义,首先要了解傅立叶原理的意义.傅立叶原理表明:任何连 ...
最新文章
- python3字典详解_python3中字典详解
- leecode第一百四十八题(排序链表)
- iOS开发之实现毛玻璃效果及图片模糊效果
- Nokia7610彩信设置
- Win32-Application的窗口和对话框
- git项目根据不同需求进行独立开发
- PHP 正在干掉 Python
- 深度学习-Tensorflow1.x-CNN中的padding参数
- 推荐系统笔记三、基于近邻的推荐系统进阶篇
- java 聊天机器人_java实现自动回复聊天机器人
- 使用Mono.cecil修改Unity游戏内存
- 从零了解进程(操作系统定位,进程的概念,特征,虚拟地址)
- 机器学习(周志华)-支持向量机课后习题:
- windows安装spacemacs
- 计算机专业名称bios翻译,电脑BIOS中英文对照翻译有哪些?
- JavaWeb笔试题
- 值得学习的C语言开源项目
- 恋爱话术小程序源码,土味情话,恋爱导师支持多种流量主模式
- Amazon.com 美国亚马逊 直邮中国 手把手教程(转)
- [KO机器学习] Day 7 模型评估:评估指标的局限性
热门文章
- 9006 - ProxySQL Error: connection is locked to hostgroup 2 but trying to reach hostgroup 1
- Shiro权限控制+整合shiro
- mysql 重新初始化
- 【面向对象程序设计】侩子手游戏(Java、JavaFX)
- 净水行业首家,安吉尔新水效国标检测能力获CNAS认可
- 实验吧-密码学-Fair-Play(Playfair解密)
- 【Python】如何利用python对c程序源码进行协助解读(学会事半功倍)
- Freeswitch 安装
- 游戏建模师具体干什么!30岁后进入行业算不算晚?
- Twitter开发者账号及开发者APPs的创建 2019.05