FFT與多項式乘法

2021-07-13 07:02:01 字數 2861 閱讀 5865

鑑於網上關於fft作用的文章並不少,這裡盡量少說廢話,直接說如何理解和具體實現

關於多項式乘法,通常是用係數乘積的方式,這樣的時間複雜度是o(n^2) n為多項式項數。這可以滿足大多數的乘法需求,然而當位數大於1000時,此法用起來就顯得捉襟見肘

於是我們便有了偉大的fft來解決這個問題

係數表示法:即平時看到的多項式:σi=0~n-1 a[i]*(x^i)

點值表示法:乙個最高次為n-1次的多項式f(x),可以表示為n個其影象上點(x,y),例如2x^2+3x+1可以表示為(0,1) (1,6) (2,15)

兩種方法可以互相轉換,而點值表示法下兩多項式相乘是非常方便的,只要將同乙個x對應的兩式y值相乘,就是新多項式對應的點值,當然,所取點值數量需要與結果多項式相符,例如:

(x+1)*(2x+1) = 2x^2+3x+1

x+1 -> (0,1) (1,2) (2,3)

2x+1 -> (0,1) (1,3) (2,5)

* -> (0,1) (1,6) (2,15) -> 2x^2+3x+1

而兩種表示法的轉換就是dft(離散傅利葉變換,discrete fourier transform)和idft(逆離散傅利葉變換,inversed discrete fourier transform),但是基於實數的dft和idft效率仍為n^2,在此不做累述

既然實數不行,那麼我們就用性質特殊的虛數,沒學過虛數的要理解起來可就費點勁啦

對於虛數就說兩句:首先設a,b為實數,則任意虛數可以表示為a+bi (i為-1的平方根)

那麼任意虛數就可以表示為下圖複數平面上乙個點

1的2^n次復根有2^n個,分別為複數平面的2^n等分線上單位向量

如4次復根分別平面上點(1,0)(0,1)(-1,0)(0,-1) 對應 1,i,-1,-i  //帶回可發現這四個值的四次冪真的都是1

其中從x軸正半軸開始逆時針找到的第乙個稱作單位根,如四次單位根為i,二次單位根為-1

n次單位根的k次冪可以用e^(2kπi/n)表示

復根有一些特殊性質,所以可以加速dft和idft

為方便表達,我們稱n次單位根的k次冪為w(n,k)

1、w(n,2*k)=w(n/2,k) //n為偶數

2、w(n,k)=-w(n,k+n/2) //n為偶數

理解了如上性質以後,就可以開始學習fft了

快速傅利葉之所以快,就是因為它使用分治策略,將多項式分為奇次項和偶次項處理。

對於a(x)=a[0]+a[1]*x+a[2]*x^2+...+a[n-1]x^(n-1) //n為偶數

設o(x)=a[1]+a[3]*x+a[5]*x^2+...+a[n-1]x^(n/2-1)

e(x)=a[0]+a[2]*x+a[4]*x^2+...+a[n-2]x^(n/2-1)

則有a(x)=x*o(x^2)+e(x^2)

fft將n次復根作為點值的x,根據復根性質1,x^2=(w(n,k))^2=w(n/2,k) 於是,結合上面的式子,我們成功地將n項問題轉化為了n/2項問題

那麼k大於n/2時怎麼辦呢?由復根性質2,w(n,k-n/2)=-w(n,k) 所以w(n,k)^2=w(n,k-n/2)^2,只要計算出w(n,k-n/2)對應的o和e函式值,w(n,k)也就自然求出了,這就是所謂的蝴蝶操作

分治的終點是當n=1,y[i]=a[i]*w(1,0)=a[i],時間效率o(nlogn)

因為每次分治都需要n為偶數,n必須為2的冪,如果不足就補上係數為0的項

根據如上敘述,可以寫出fft的遞迴實現方法,雖然後面的迭代實現更為優秀,但是強烈建議先理解遞迴演算法,因為這是理解fft的基礎

recursive-fft(a)

n=a.length

if n==1

return a

e= o=(a[1],a[3],...,a[n-1]}

y_e=recursive-fft(e);

y_o=recursive-fft(o);

for k=0 to n/2-1

w=e^(2πki/n)

y[k]=y_e[k]+w*y_o[k]

y[k+n/2]=y_e[k]-w*y_o[k]

return y

其中a為係數向量,返回的y陣列即為n次復根對應的n個值

上述**花費了大量空間用於建立和維護陣列,而我們可以使用迭代法避免這一部分的空間使用

為了迭代,首先陣列被分治到的部分必須連續,觀察遞迴呼叫的過程可以發現,分治是按照二進位制的末位開始分類的,如n=8時,分治順序如下

0,1,2,3,4,5,6,7

(0,2,4,6)(1,3,5,7)

(0,4)(2,6)(1,5)(3,7)

寫成二進位制000,100,010,110,001,101,011,111

將二進位制翻轉000,001,010,011,100,101,110,111這不正是0~7的單增序列嗎?

原理是翻轉後末位變首位,而高位恰恰決定了數的大小

於是我們遞推預處理出rev陣列

void get_rev(int bit)

然後根據rev對a重新排序,就可以進行迭代fft了,我們已知含乙個元素的dft為係數本身,那麼將其按次序使用蝴蝶操作兩兩合併就能得到兩個元素的dft,然後再合併得到四元素dft,以此類推,就可以得到整個式子的dft

typedef complexcd;//c++ 自帶複數類,需要標頭檔案complex

void fft(cd *a,int n)

void fft(cd *a,int n,int dft)

FFT與多項式乘法

係數表示法 即平時看到的多項式 i 0 n 1 a i x i 點值表示法 乙個最高次為n 1次的多項式f x 可以表示為n個其影象上點 x,y 例如2x 2 3x 1可以表示為 0,1 1,6 2,15 兩種方法可以互相轉換,而點值表示法下兩多項式相乘是非常方便的,只要將同乙個x對應的兩式y值相乘...

UOJ 34 多項式乘法(FFT)

description 給你兩個多項式,請輸出乘起來後的多項式 input 第一行兩個整數 n 和 m,分別表示兩個多項式的次數 第二行n 1 個整數,分別表示第乙個多項式的 0 到 n次項前的係數 第三行m 1 個整數,分別表示第乙個多項式的 0 到 m次項前的係數 output一行n m 1 個...

fjut 3283 FFT多項式乘法

本題簡單剖析出來就是固定其中乙個串然後另乙個串旋轉到乙個位置值得答案值最大.我們發現每次求值的複雜度是o n 一共有n種旋轉方法,所以總複雜度是o n 2 如果把這n種都綜合起來的話不就是近似多項式乘法了嗎?那麼怎麼轉換呢?看下面的圖 以n 5為例,那麼5種旋轉方式的值分部在多項式的 呢?從圖看出就...