快速求n阶多项式乘积

设有两个n阶多项式

A(n)=an-1x^n-1+an-2x^n-2+...+a0

B(n)=bn-1x^n-1+bn-2x^n-2+...+b0

则如何求A(n)与B(n)的乘积？

通常的方法是

C(n)的表达形式是

C(n)=c(2n-2)x^(2n-2)+c(2n-1)x^(2n-1)+...+c0

易知C(n)的第k项系数ck=sigma(ai*b(k-i))(0<=i<=k)

很明显这个算法是Θ(n^2)

 

所以我们从多项式的基本表达上去研究

一个(n)阶多项式的表达方式有几种？

最常见的方法就是上文的解析式，也就是系数表达式。f(x)=an-1x^n-1+an-2x^n-2+...+a0

但是这种方法在乘法时的时间复杂度是n^2太慢了！

另一种方法就是点值表达式，这种方法我们初中就学到过——在求解一次函数多项式的解析式（系数表达时），我们只需要找2个点就确定了其结构，在求解二次函数时也需三个点就能求出解析式，那么这个结论能否向高次方程扩张呢？

答案是可以的

猜想：对于n阶多项式，我们知道n个数对(xi,yi)且任意xi不相等，就能确定一个n阶多项式经过这些点

证明：这是与线性代数相关的命题

{1 x1 x1^2 x1^3... x1^n-1}    a0      y1

{1 x2 x2^2 x2^3... x2^n-1}    a1      y2

......................................... *       =  

.........................................

{1 xn xn^2 xn^3... xn^n-1}  an        yn

其行列式不为0（详见算导） 所以有且仅有一解；

所以我们可以由一组n个元素的数对(xi,yi)求解出n阶多项式的系数表达。这称为插值。

点对表达式的有点在于若两个n阶多项式相乘，如A（n）*B（n）

设A（n）的点对表达为n个（Axi，Ayi）

B（n）的点对表达为n个（Bxi，Byi）

那么A（）与B（）的乘积的点值表示就是

n个点（Axi*Bxi,Ayi*byi）

这里确定一个2n的多项式仅用n个点是不够的

所以上文中把A（）和B（）都得取2n个点对

现在还有一个插值的逆操作——求值

及知道一个n阶多项式的系数表达，如何在θ（n log n）的时间内求出n个元素的数对(xi,yi)？

这里需要用到单位复数根的性质

另单位复数(模长=1)e=cosα+isina;

则e^2=cosα*cosα-sinα*sinα+(2cosα*sinα)i=cos2α+sin2αi

及e的模长不变，角度转过一倍。

则方程x^n=1的n解是什么

由上述知解的模长必为1，切设其与x轴夹角为α，nα=2kπ。

xi=sin(2πi/n)+cos(2πi/n)i

x1..xn成为n阶的单位复数根

于是我们可以取x1,x2..xn这几个单位复数根，再求出这几个点时的y值

这个方法可以用递归进行，就是传说中的DFT啦

f(xi)=an-1xi^n-1+an-2xi^n-2+...+a0

定义F'(xi)=a0+a2xi+a4xi^2+...+an-2 x^(n/2-1)

和 F"(xi)=a1+a3xi+a5xi^2+...+an-1 x(n/2-1)

f(xi)=F'(xi^2)+xi*F"(xi^2)

我们还发现其实对于x1..xn这个数在复数系中的位置的个数

x1^2,x2^2..xn^2更少（再转了α后半圈xi与前半全重合）

甚至有对于even(n),前一个的规模恰好是后一个的两倍

现在我们要计算x1,x2..xn这n个点的平方在当前系数下的值（就变成了x2,x4..xn这n/2个点，数据规模缩小了一半）

所以T(n)=2T(n/2)+θ（n）（对每个xi利用F'和F“”的值组合出f）

根据主定理可知时间复杂度为n log n

//见算导P533 这里求的是w0..wn-1,与上文中的1..n做出一点调整（wn=w0）还有不用记录当前传进去的几个单位复数根，因为单位复数根的个数与系数向量的模长相等，直接生成即可。
FFT(a)// a is a vector including n elements
{
　　n=a.length();
　　if(n==1) return a;
　　wn=e(2πi/n);//就是最小转一次需要乘的值
　　w=1;//方便线扫求出wi
　　A0=(a0,a1,...,an-2);
　　A1=(a1,a3,...,an-1);
　　y0=FFT(A0);
　　y1=FFT(A1);
　　// y0,y1 is a vector including n/2 elements就是求出的F',F"啦
　　for k=0 to n/2-1
　　　　y[k]=y0[k]+wy1[k];
　　　　y[k+(n/2)]=yk[0]-wy1[k]；
　　　　w=w*wn;
　　return y;
}