辽宁省沈阳市第十五中学2020年高中数学论文 图形计算器应用能力测试活动学生 用快速傅立叶变换实现多项式乘法

辽宁省沈阳市第十五中学2020年高中数学论文 图形计算器应用能力测试活动学生 用快速傅立叶变换实现多项式乘法宁波效实中学 1502班 周杨Input:Expand(x+2)(x+3)Output:x2+5x+6在CLASSPAD330中，main模块中Interactive栏目Transformation中的expand一直是广大高中生喜闻乐见的好功能。用法简洁易懂。但是，当你需要算两个很长很长的多项式的乘积的时候，我相信你一定会被冗长的表达式输入问题伤透脑筋（我*，怎么要输这么多的“xm”各种残念）。今天，我就带大家认识一种高贵冷艳的计算多项式乘法的方法。在进入正式介绍之前，我们先要普及一些前提知识。那就是FFT。FFT。俗称“法法塔”。学名快速傅立叶变换。FFT有什么用呢？它可以在O(nlogn)的时间复杂度内完成一个函数对n个单位虚根的求值，比朴素的O(n2)要厉害得多。下面的几段是关于FFT具体操作过程的介绍，如果没有兴趣可以跳过，因为对我们的实际操作没有太大影响,CLASSPAD330将会为我们完成它的具体工作的。单位根定义n次单位根为满足wn=1的复数w。恰好有n个单位根e2ik/n(k=0, 1, .,n-1)，其中i是虚数单位(i2=-1)，而eiu = cos u + i sin u定义wn=e2ik/n，称为n次主单位根(principal nth root of unity)，则所有其他单位根都是它的幂。这样，n个单位根可以写成wn0,wn1,.,wnn-1 ，它们在乘法下构成群，且结构和加法群Zn一样。容易看出，wnk=-wnk+n/2。设多项式为，那么我们将奇数项和偶数项分开处理。A0(x)=a0+a2x+.+an-2xn/2-1A0(x)=a1+a3x+.+an-1xn/2-1那么A(x) = A0(x2) + xA1(x2)。前面提到过wnk=-wnk+n/2，因此(wnk)2=(wnk+n/2)2这样，我们可以同时计算A(x)在这两个点的值。这样，计算A(x)在n个点上值的任务被转化成了计算A0(x)和A1(x)在n/2个点上值的任务了。上式的操作实际上是用y0k和y1k 算yk和yk+n/2，这样的操作称为蝴蝶操作(butterfly operation) ，如下图：这样，只需要花O(n)的时间就可以把序列y0和y1合并成序列y了。递归方程为T(n)=2T(n/2)+O(n)，因此解为T(n)=O(nlogn)。远低于O(n2)。实际操作的时候可以用迭代来完成，这样就可以避开傲娇卖萌极慢无比的系统栈了。唯一可惜的是CLASSPAD330上的program无法使用数组，所以暂时无法在CLASSPAD330上手动实现一个FFT（包括后面要讲的IFFT）。读者可以自己用C+或其他高级语言在电脑上实现一下。难度应该不大。为了能够更好的解决问题我们还需要IFFT，即FFT的逆操作。因为是逆操作，所以我们可以不加思索地推出只需要把上面计算过程中的w换为-w即可。接下来让我们切入正题。也许你会觉得这个高贵冷艳的FFT和我们要用的多项式乘法八竿子打不着。但是它的确可以在实践中加速多项式乘法的操作，从O(n2)到O(nlogn)。对于一个多项式，我们有两种表示法。系数表示法和点值表示法。系数表示法就是记录系数向量a=(a0,a1,.,an-1)。注意要用列向量。这样若两个多项式的系数向量为a,b。那么多项式加法就相当于a+b（这应该是非常直观的），而多项式乘法就对应于向量卷积ab。点值表示法是一个点-值对的集合(x0,y0),(x1,y1),.,(xn-1,yn-1)，其中所有xi两两不相等且yk=A(xk)。这就意味着一个多项式有无穷多种点值表示法。但是如果限定多项式的次数小于n，那么点值表示法和多项式就一定是一一对应的了。那么我们可以这样进行多项式乘法的计算：1. 将两个多项式的系数表示法高位补0，得到两个次数为2n的多项式。2. 将两个多项式的系数表示法转换为点值表示法。3. 把两个点值对应的对应点的数值即yk相乘得到C即结果的点值表示法。4. 再把C的点值表示法转换为系数表示法。这也是应用FFT进行多项式计算的原理。为了应用FFT，我们不妨将点都取单位虚根，这样第2步就变成了FFT，第4步也就编程了IFFT。那么我们重新表述一下上述过程1. 将两个多项式的系数表示法高位补0，得到两个次数为2n的多项式。2. 把2n次单位根作为求值点，用FFT计算两个多项式的点值表示法。3. 把两个点值对应的对应点的数值即yk相乘得到C即结果的点值表示法。4. 用IFFT把C转换为系数表示法。接下来我们用CLASSPAD330来实现一下上述过程。首先我们引入一个具体问题。A(x)=1+2x+6x2B(x)=3+5x+9x2+7x3求两个多项式的乘积。则系数向量a=(1,2,6),b=(3,5,9,7)。第一步添0，则a=(1,2,6,3,0,0,0,0,0)，b=(3,5,9,7,0,0,0,0,0)。则Input:a:=1,2,6,0,0,0,0,01,2,6,0,0,0,0,0Input:b:=3,5,9,7,0,0,0,03,5,9,7,0,0,0,0第二步计算FFT。不过在计算前要记得把模式调为Cplx（虚数）模式。Input:x:=FFT(a)Input:y:=FFT(b)第三步求对应点乘。Input:c:=x*y第四步求IFFTInput:z:=IFFT(c)Output: 3,11,37,55-1.25e-12i,68,42,0,-1.25e-11+1.25e-12i这个答案看上去似乎很奇怪。但仔细观察你会发现后面跟的那个虚数几乎可以忽略不计。其实这是在FFT和IFFT计算过程中产生的精度问题。一般情况下我们可以直接忽略后面的虚数。那么得到的就是乘积的系数表示，即c=(3,11,37,55,68,42).但实际上，这个精度的问题也是可以避免的。只需要将第二步到第四步写在一个表达式里面就可以了。如下所示：Input:IFFT(FFT(a)*FFT(b)Output: 3,11,37,55,68,42,0,0最后，附上一张模拟器运行图。最后我们验证一下。记得先用Edit栏里的Clear All Variables清空一下变量。Input:expand