1. 程式人生 > 其它 >快速傅立葉變換FFT修訂版

快速傅立葉變換FFT修訂版

目錄

以前我也學習過FFT
但是你也看見了,內容很冗雜,而且當時什麼鬼都不懂。
所以寫的很馬虎
現在我學廢了一點點
就回來寫了這篇文章。
所以在這片文章裡不會有太多的前置知識,需要的知識會使用超連結告訴你。
主要集中在原理方面
最重要一點:
請各位讀者有任何問題務必在評論區提出,也好讓我知道在這方面我有什麼不足。

多項式乘法

快速傅立葉變換最簡單的用法就是用在多項式乘法
現在給你兩個函式\(f,g\)

他們的次數分別是\(n_1,n_2\)
那麼現在要求\(h=f\cdot g\),他的次數就應該是\(n=n_1+n_2\)下文的\(n\)就是在這裡定義的
那麼我們知道要確定一個\(n\)次的多項式,如果我們把它看作一個函式的話,我們就需要知道在函式圖象上的至少\(n\)個點
FFT就是通過構造這\(n\)個點來求出\(h\)的解析式的。
所以就涉及到係數轉點值點值轉系數兩個關鍵操作。
係數是什麼呢?我們知道我們有一個多項式,用\(f\)舉例吧,這是一個\(n_1\)次多項式,所以就可以表示為\(n_1+1\)個單項式相加,他們就是

\[f(x)=\displaystyle\sum_{i=0}^{n_1}a_ix^i \]

其中\(a_i\)

就被稱為\(i\)次項係數我們就可以用\(n_1+1\)個係數來表示一個\(n_1\)次多項式。
同樣,我們要知道\(h\)的解析式,我們就需要確定\(n+1\)個係數
然後初中的待定係數法告訴我們,我們有\(n+1\)個未知的係數,我們就需要至少獲得\(n+1\)個點
所以我們就取\(n+1\)個互不相同的\(x_i\)代入到\(h(x_i)=f(x_i)\cdot g(x_i)\),就可以求解出\(h\)的解析式。
理論上。
實際實操起來還是比較艱難的。
首先就是怎麼選取\(x_i\)
雖然說這個\(x_i\)是隨便取的,但是這樣我們計算\(n+1\)個點,每個點的時間複雜度都是\(n_1\)
,這樣時間複雜度很不划算,還不如直接\(n_1n_2\)搞暴力。
我們就要取一些\(x_i\)的值,使得這個時間複雜度可以蜜汁迷之降低。
這個我們後面再討論

快速傅立葉變換

我們先看看我們可以怎麼轉化掉這個\(f\),把它變成點值。
看好了

\[f(x)=a_0x^0+a_1x^1+a_2x^2+\dots \]

這是係數表達方法
然後我們按照次數的奇偶性分類

\[f_0(x)=a_0x^0+a_2x^1+a_4x^2+\dots\\ f_1(x)=a_1x^0+a_3x^1+a_5x^2+\dots \]

這樣分類之後就出來了兩個次數為\(n_1\div 2+1\)的多項式
然後我們是有

\[f(x)=f_0(x^2)+x\cdot f_1(x^2) \]

是不是很神奇?我也不知道FFT怎麼想出來的但是反正成立就行了。
我們就可以抓住這個不放,想一下在這個式子上怎麼優化。
我們根據高中時候做偶函式題目的直覺
我們發現後面傳進去的都是\(x^2\),唯獨前面的係數是一個\(x\)
就是如果沒有這個\(x\),我們就會有\(f(x)=f(-x)\)
但是他有
怎麼辦?
其實也是一樣的。你在求\(f(x)\)的時候就要求\(f_0(x^2),f_1(x^2)\)的值,然後在合起來。
那麼\(f(x)\)\(f(-x)\)只是一個係數上的不同,因為這時候你的自變數已經確定了,所以我們可以暫且認為,這兩個式子只是在一個係數上不一樣。
所以我們在求\(f(x)\)的時候求出來的\(f_0(x^2),f_1(x^2)\)的值在求\(f(-x)\)的時候一樣也需要用,這樣我們相當於用\(O(1)\)求出了\(f(-x)\)
所以我們只要讓我們取的\(x_i\)滿足

\[\tag1 \forall i,j\in[0,n+1], i\neq j, x_i\neq x_j\\ \forall i\in[0,n+1], \exist ! j \in [0,n+1], j \neq i, x_i+x_j=0 \]

這樣我們就可以把求解的範圍縮小一半。
但是這一半貌似不太夠。
所以我們觀察\(f_0,f_1\)
發現他們也是多項式!
他們轉點值也可以這樣二分。
前提是滿足\(n_1 \div 2 +1\)同樣是奇數。
然後他們又二分,又二分。。。
一路二分下去,最後的時間複雜度就降到了\(O(n \log_2n)\)
需要滿足的條件就是當前多項式的次數必須是奇數,除非你遞迴到了邊界:0次
所以我們就要求對於所有的次數必須都是偶數,也就是一開始多項式的次數必須是2的整數次冪。
所以前面我們提到的\(n\)\(\leftarrow\)點一下)最小取值就應該是

\[\displaystyle n=2^{\lceil\log_2(n_1+n_2)\rceil} \]

然後我們回到\(x_i\)的取值問題。
我們設當前處理的多項式次數為\(m\),很明顯這個東西是2的整數次冪。然後\(hf=m\div 2\)
那麼我們分出來的兩個子函式\(f_0,f_1\)的次數都應該是\(hf\)
然後由\(f_1\)組成\(f\)中次數為奇數的,\(f_0\)組成是偶數的。
那麼我們選取這個\(x_i\)就很講究了,應該要滿足

\[\forall i \in [0,hf) x_i+x_{i+hf}=0 \]

也就是我加上一半的下標之後的數和我是相反數。
那很簡單,我們只需要一個數組\(x\),在裡面隨便構造。
然後你就發現了一個問題。
當你要二分下去的時候,你發現你的子函式。。。他所需要的下標不知道怎麼取了吧。。。因為你的子函式除了係數下標奇偶性要考慮之外你還要讓他相同的係數代入相同的\(x_i\)
很難取吧。。。

如果你覺得這都有可能的話可以學完單位根的方法之後再把這個方法寫出來。
你寫的出來,我當場就
把它加入到這個部落格裡面以表示你的大聰明。

這個時候聰明的不知道誰(FLY?)提出了使用單位根
首先證明充分性,為什麼代入單位根可以呢?因為單位根是複數,所以。。。我上面也沒說用複數求解不行啊。。。複數求解是可以的,這個理由很充分。
然後我們知道當前處理的多項式的次數\(m\)是2的整數次冪,所以我們的\(m\)一定是偶數,然後我們構造\(m\)次單位根,設他們之中輻角最小的為\(\omega_m\),那麼所有的單位根的集合就可以表示為\(\{\omega_m^i|i\in[0,m)\}\)。因為\(\omega_m^0=\omega_m^m=1\)(參考單位根定義)所以這裡是左閉右開區間。
就會滿足\(hf=m\div 2,\omega_m^i=-\omega_m^{i+hf},i\in[0,hf)\)
所以用單位根也滿足了我們的式子\((1)\)(往上翻一翻,在頁面右邊有標記某個數學公式的編號為\((1)\)
然後證明必要性,我們為什麼有使用單位根的必要呢?
這裡我們已經通過感性想象知道了使用陣列提前搞我是搞不定的。
然後我們看看單位根啊。
按照奇偶性分好子函式之後,他們的次數就是\(hf\)了,他們執行的時候分別構造\(hf\)次單位根,可以保證

\[\{\omega_{hf}^i|i\in[0,hf)\} \subsetneq \{\omega_m^i|i\in[0,m)\} \]

我們粗略證明一下,\(hf\)次單位根必定都是\(m\)次單位根。
\((\omega_{hf}^i)^m=(\omega_{hf}^{m=2hf})^i=(\omega_{hf}^{hf})^{2i}=1^{2i}=1\)他的\(m\)次冪為\(1\)他就是\(m\)次單位根。
而且根據一些感性的理解(在單位圓裡面隔一個單位根取一個)我們可以得知,這個\(\omega_hf^i\)和子函式獲得的係數在原函式裡是一一對應的(因為係數也是隔一個取一個才能滿足下標的奇偶性相同)
這就完美解決了我們搞不定的問題。
所以我們每次就二分,求解,可以通過\(O(n\log_2n)\)求出\(h\)的點值表達。

快速傅立葉逆變換

上一句話我們說,這個是點值表達。
但是我們需要的是係數表達。
那麼我們要給他轉換回來。
我們設\(h\)的係數分別是\({a_0,a_1,a_2,\dots,a_{n-1}}\)
然後我們求出來的點值分別是\({y_0,y_1,y_2,\dots,y_{n-1}}\)
那麼就應該滿足有全稱命題\(\displaystyle\forall i\in[0,n),h(\omega_n^i)=y_i\)
\(h(\omega_n^i)\)暴力拆分一下(貌似這裡以前那個還推錯了。。。)

\[\displaystyle h(\omega_n^i)=\sum_{j=0}^na_j\cdot(\omega_n^i)^j=y_i \]

看上去還是很難的。
然後我們考慮一下DFT的某些過程
這裡是DFT不是FFT
DFT就是離散傅立葉變換,FFT可以認為是求解DFT的一種快速的方法。
求解DFT的正常方法很簡單,就是我在開頭說的,直接隨便取\(n\)個值帶入進去完事,不需要什麼單位根甚至不需要滿足相反數。
這樣的做法是很高時間複雜度的。
但是他給我們提供了一個思路
他求點值的方法就是直接暴力。
如果我們代入單位根,但是我們同樣暴力。
那麼這個過程就變成了一個矩陣乘向量的過程(注意不是矩陣乘法

\[\begin{bmatrix} h(\omega_n^0)=y_0\\ h(\omega_n^1)=y_1\\ h(\omega_n^2)=y_2\\ \vdots\\ h(\omega_n^3)=y_n \end{bmatrix} = \begin{bmatrix} (\omega_n^0)^0&(\omega_n^0)^1&(\omega_n^0)^2&\cdots&(\omega_n^0)^{n-1}\\ (\omega_n^1)^0&(\omega_n^1)^1&(\omega_n^1)^2&\cdots&(\omega_n^1)^{n-1}\\ (\omega_n^2)^0&(\omega_n^2)^1&(\omega_n^2)^2&\cdots&(\omega_n^2)^{n-1}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ (\omega_n^{n-1})^0&(\omega_n^{n-1})^1&(\omega_n^{n-1})^2&\cdots&(\omega_n^{n-1})^{n-1}\\ \end{bmatrix} \begin{bmatrix} a_0\\a_1\\a_2\\\vdots\\a_{n-1} \end{bmatrix} \]

然後現在倒過來,我知道了等號左邊的向量要求右邊的向量。
爾曹可知逆元?
算了就是逆矩陣。
我們就有

\[\begin{bmatrix} a_0\\a_1\\a_2\\\vdots\\a_{n-1} \end{bmatrix} = \begin{bmatrix} (\omega_n^0)^0&(\omega_n^0)^1&(\omega_n^0)^2&\cdots&(\omega_n^0)^{n-1}\\ (\omega_n^1)^0&(\omega_n^1)^1&(\omega_n^1)^2&\cdots&(\omega_n^1)^{n-1}\\ (\omega_n^2)^0&(\omega_n^2)^1&(\omega_n^2)^2&\cdots&(\omega_n^2)^{n-1}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ (\omega_n^{n-1})^0&(\omega_n^{n-1})^1&(\omega_n^{n-1})^2&\cdots&(\omega_n^{n-1})^{n-1}\\ \end{bmatrix}^{-1} \begin{bmatrix} y_0\\y_1\\y_2\\\vdots\\y_n \end{bmatrix} \]

根據大量的草稿紙和代數運算我們求出中間最大的那個矩陣的逆元

\[\begin{bmatrix} (\omega_n^0)^0&(\omega_n^0)^1&(\omega_n^0)^2&\cdots&(\omega_n^0)^{n-1}\\ (\omega_n^1)^0&(\omega_n^1)^1&(\omega_n^1)^2&\cdots&(\omega_n^1)^{n-1}\\ (\omega_n^2)^0&(\omega_n^2)^1&(\omega_n^2)^2&\cdots&(\omega_n^2)^{n-1}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ (\omega_n^{n-1})^0&(\omega_n^{n-1})^1&(\omega_n^{n-1})^2&\cdots&(\omega_n^{n-1})^{n-1}\\ \end{bmatrix}^{-1}= \frac1n \begin{bmatrix} (\omega_n^0)^{-0}&(\omega_n^0)^{-1}&(\omega_n^0)^{-2}&\cdots&(\omega_n^0)^{-(n-1)}\\ (\omega_n^1)^{-0}&(\omega_n^1)^{-1}&(\omega_n^1)^{-2}&\cdots&(\omega_n^1)^{-(n-1)}\\ (\omega_n^2)^{-0}&(\omega_n^2)^{-1}&(\omega_n^2)^{-2}&\cdots&(\omega_n^2)^{-(n-1)}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ (\omega_n^{n-1})^{-0}&(\omega_n^{n-1})^{-1}&(\omega_n^{n-1})^{-2}&\cdots&(\omega_n^{n-1})^{-(n-1)}\\ \end{bmatrix} \]

然後我們就發現,我們求\({a_0,a_1,\dots}\)的過程和我們求\({y_0,y_1,\dots}\)的過程是很像的。
只是我們的單位根是倒著取的。這個小小的變動很容易實現,只需要更改一下迴圈變數的更新就可以了。
然後最後輸出的時候不要忘了逆矩陣前面還有一個係數\(\frac1n\)
此處的參考資料:快速傅立葉變換(FFT)——有史以來最巧妙的演算法?

Code

建設中。。。




寫在後面 本文來自部落格園,作者:IdanSuce,轉載請註明原文連結:https://www.cnblogs.com/LASS/p/15661018.html 覺得寫的闊以的什麼都不用留下,哪裡寫的不行的歡迎在評論區指出或聯絡作者IdanSuce
作者喜歡的網站:https://oi-wiki.org/




作者喜歡的圖片:。。。也不知道為什麼沒有設定成背景。。。可能要致敬紙嫁衣吧