正文:
文章寫的有點急。有錯誤的地方望指出
我學習 FFT 是一個比較慢的過程。 期間反反復復。 我寫這篇博文只是一個非常淺顯的理解。同時也可以幫助初學者在學習FFT的時候。有所偏重。避免太多思維上的負擔。
直接正題吧:
首先:本身并不負責多項式之間的乘法。
只是一種變換。
則是DFT的快速算法。(分治 提高效率)
利用 。我們快速的將多項式變換為利于計算的形式
用這種方便計算的形式計算出來兩個多項式的乘積。
這時候我們雖然已經得到目標多項式。但其形式并不是我們想要的
所以 之后利用 的 逆運算 又快速的變換回去
我們記:的逆運算為 1
對于一個次多項A的表示。最常見的形式(系數(shù)表達):
小學生手算 系數(shù)形式 的 多項式乘法 的 復雜度是
如果我們知道了一個次多項式 的曲線上的個不同的的點。
我們是可以計算出來這個多項式的
把系數(shù)看做未知數(shù)。列出來n個方程組。解出來這n個系數(shù)。
也就是說,給出曲線上的n個點:
我們可以確定其系數(shù)形式;
我們稱這種表達一個多項式的方法叫:點值表達
他有很多優(yōu)點。比如可以時間內計算一個多項式乘法
再給出一個多項式:
設與在取相同x的點值表達為:
這里 ,
則為:
非常方便。順其自然;
之所以與多取一些點。是因為的次數(shù)界增加。
取點不足會導次數(shù)界達不到
對于一個次多項式。隨機求其n個不同的點的樸素方法復雜度是
假設 n為偶數(shù)。那么我們把。重組為兩個多項式:
設其系數(shù)的下標為偶數(shù)組成的多項式為:
設其系數(shù)的下標為奇數(shù)組成的多項式為:
所以有:
好像并不會優(yōu)化運算。
因為
依然有可能會組成n個不同的數(shù)字。
那么我們要計算的規(guī)模不會減半。
如果我們恰當選取一些x。是否會優(yōu)化運算呢。
例如:
各個數(shù)字平方后。得到的集合大小減半:
因為
那么
但是這種關系不會傳遞下去。平方后得到的數(shù)字全是不小于0的數(shù)字。
再次遞歸又回到了原來的形式。很尷尬。
不過這啟迪我們。取相反數(shù)。然后平放??梢园褑栴}規(guī)模減半。但再一次遞歸就失效了。是否存在不失效的取值呢。
如果對于一個有偶數(shù)個元素的數(shù)字集合。每個元素平方后。得到的新集合。去除重復元素后。集合大小能夠減半。并且得到的新集合如果為偶數(shù)。新集合依然滿足上面的性質。我們稱這個集合有 折半 的性質。
如果我們可以快速的找到一個滿足折半性質的自變量的取值集合.
分治就是可行的。
有一種東西??梢詽M足我們的需求。
那就是 —— n次單位復數(shù)根:
(使用復數(shù)確實讓人有顧慮。這也是我學習FFT時最大 思想負擔)
我們定義
那么形如
的復數(shù)有諸多良好的性質。
例如:
上面的性質可以用數(shù)學歸納法得到(較為簡單。這里不做證明
我們記:
我會告訴你 多項式在x取下面的值時。有助于我們進行DFT
上面的n個復數(shù)稱為。n次單位復數(shù)根。(n次單位復數(shù)根是指這n個數(shù))
如果我們x取時。根據(jù)剛才的分解:
而
因為三角函數(shù)的周期性:
我們有:
所以:
則:
驚不驚喜,意不意外。
的取值集合取單位復數(shù)根。不但滿足折半的性質。而且還有一定的規(guī)律性。與原集合保持一致。
這意味著我們只需要計算取:
的
問題范圍減半。并且如果n為偶數(shù)。再次平方依然集合大小減半。
記:
那么,當我們得到:
這意味著我們得到了所有的
則:
即得到了:
其中
當然。FFT的計算 n要取2的整數(shù)次冪 。這是因為每次減半。我們可以把不足的系數(shù)用0填充。
上面過程可以看作 ,取單位復數(shù)根計算目標y向量,如下圖:
本著盡可能簡單的原則。我們不在特別的說這個矩陣。
因為之前說關于系數(shù)的n個方程必然可解。所以上面的矩陣:
必然存在逆矩陣5。(有線性代數(shù)基礎的,或許不陌生)
可能你還是有點迷惑。沒關系:
如果我們把FFT看作計算上面特別的矩陣相乘的算法呢。
利用FFT我們可以快速得到:
同時 。我可以直接告訴你:
的逆矩陣為:
因為我直接給你了答案。。。所以懷疑的話。我們可以計算一下看看。
我們記Y=E*D
上面可以看作等比數(shù)列求和。
當時:
當時,令:
因為:
所以,時:
所以。Y為單位矩陣。E,D互為逆矩陣得證。
所以:
所以。其實就是的變?yōu)榈玫健?/h4>
并且所得目標向量每個元素除掉n就可以了6
FFT的高效實現(xiàn):
通常。我們希望FFT的實現(xiàn)盡可能快。所以FFT算法也都使用迭代結構而不是遞歸結構。
下面介紹一種常用的去遞歸方法。
首先。上面介紹的FFT是只能處理2的整數(shù)次冪的次數(shù)界的多項式
所以不特別說明。都有(系數(shù)不足的用0補足)
對于遞歸的第一步:
輸入數(shù)組
被分解為,偶數(shù)和奇數(shù)兩個部分數(shù)
下標二進制形式右起第一個字符為 0 被分配到左集合
下標二進制形式右起第一個字符為 1 被分配到右集合
更一般的。對于第k次遞歸:
下標右起第k個字符為 0 被分配到它所在問題的左集合
下標右起第k個字符為 1 被分配到它所在問題的右集合
那么對于一個2進制形式的數(shù)字 B。將其表示為(注意:):
根據(jù)之前分析。這個數(shù)字遞歸后將會被分配到的位置為:
所以。對于一個數(shù)字的二進制形式的前k位對稱過來。就是遞歸后的位置。
例如:,
我們調用FFT前對數(shù)組進行一次上述重拍。
便可 自底向上迭代 實現(xiàn)FFT
總結:
由于遞歸結構的FFT較好理解。而理解遞歸結構的FFT后。
不難寫出非遞歸結構
所以在這里我們對遞歸結構的FFT基本框架做一個總結。
因為使用到了復數(shù)。所以不可避免的 。我們以下計算都在復數(shù)范圍的
方便起見。我們用符號:來表示復數(shù)
那么多項式:
令:
當n=1時:
n>1時:
令: ,
計算與
計算完成后,令:
返回
————————————————————————
我對上述思想解釋一下
A在處的值保存在Y中。
并返回給上一層。所以對于當前調用:
當 時:
模版
給出一個迭代結構的FFT模版:
定義復數(shù):
struct Complex
{
double x,y;
Complex(double x1=0.0 ,double y1=0.0)
{
x=x1;
y=y1;
}
Complex operator -(const Complex &b)const
{
return Complex(x-b.x,y-b.y);
}
Complex operator +(const Complex &b)const
{
return Complex(x+b.x,y+b.y);
}
Complex operator *(const Complex &b)const
{
return Complex(x*b.x-y*b.y,x*b.y+y*b.x);
}
};
- 1
- 2
- 3
- 4
- 5
- 6
- 7
- 8
- 9
- 10
- 11
- 12
- 13
- 14
- 15
- 16
- 17
- 18
- 19
- 20
- 21
去遞歸:
void change(Complex y[],int len)
{
int i,j,k;
for(i=1,j=len/2;i<len-1;i++)
{
if(i<j)swap(y[i],y[j]);
k=len/2;
while(j>=k)
{
j-=k;
k/=2;
}
j+=k;
}
}
迭代結構的FFT
on== 1 時:
on==-1 時:
void FFT(Complex y[],int len,int on)
{
change(y,len);
for(int h=2;h<=len;h<<=1)
{
Complex wn(cos(on*2*pi/h),sin(on*2*pi/h));
for(int j=0;j<=len;j+=h)
{
Complex w(1,0);
for(int k=j;k<j+h/2;k++)
{
Complex u=y[k];
Complex t=w*y[k+h/2];
y[k]=u+t;
y[k+h/2]=u-t;
w=w*wn;
}
}
}
if(on==-1)
for(int i=0;i<len;i++)
y[i].x/=len;
}
- 1
- 2
- 3
- 4
- 5
- 6
- 7
- 8
- 9
- 10
- 11
- 12
- 13
- 14
- 15
- 16
- 17
- 18
- 19
- 20
- 21
- 22
- 23
|