計算方法研究報告

1、2016年計算方法研討課題目 Toeplitz矩陣與向量的乘積課題小組：組長：李靈萱成員：仇天樂、肖晟遠(yuǎn)、趙海寧、壽立夫、張連煒、李怡寧一個n階Toeplitz矩陣形如一 tottjt1to*aa+-t工- .t1t。1可以將T嵌入一個“循環(huán)矩陣”，然后利用FFT來快速計算矩陣-向量積Tx,其計算復(fù)雜度為 O(nlog n)。1.理論基礎(chǔ)1)一個n階循環(huán)矩陣C具有如下結(jié)構(gòu)Co2)C =circc0, g , 。=c。,cnA一個n階離散Fourier矩陣Fn為一1FnHC=Fn AFn，A =-001= diag%,%；：%上入nV12 7i ,八(nJ)e n包(n4)(n)e n這里的i是

2、虛數(shù)單位。3)數(shù)學(xué)上業(yè)已證明：一個 n階循環(huán)矩陣C能被Fn對角化，即F 2T：jk且鼠=標(biāo)(C) = cj e n (k =0,1，n 1)是C的特征值?？梢杂肐FFT計算，即 j =0IA cici2-2r-1c2c2:=VnFn : =IFFT(:)3 n 1fn 1jcn 14) 利用FFT和IFFT ,則矩陣-向量積Cx可以表為Cx = FnHAFnx=、n FFT 上IFFT(x) =FFT(A，IFFT(x).5)對于n階Toeplitz矩B$ T ,首先將T嵌入一個2n階循環(huán)矩陣C ,即TIT，其中：T一0titn0ti 1.tn0將n維列向量x嵌入一個2n維列向量y =rI這時

3、有_oT Cy =M；.T叫(PlT JLoJ bTx_根據(jù)4)的事實，可用2n階的FFT和IFFT計算矩陣-向量積Cy ,于是間接獲得矩陣-向量積Tx。這樣，矩2、陣-向重積Tx的計算受雜度從 O(n )降低到了 O(n log n)。6) 公式的數(shù)學(xué)證明我們令矩陣A=FnCFnH,則需要證明A為對角化的矩陣。根據(jù)矩陣的結(jié)合律，令 B=FnC一Tij巖j記列向重 c(i) =c0,G,，cn. ,i =0,1,.,n1,且Wn =en ,并且記c(i j) =c(i j)nR(n),即c(i j)以n為周期延拓后序列的主值序列，也就是循環(huán)矩陣 C 的第 j 個列向量，得到 C =c(i),c

4、(i -1),c(i -2),., c(i -n-1)on 1-1記 y(k)= JnMFnc(i)= c(i)W1k =/IDFT(c(i), i =0n1jk則由DFT的性質(zhì)可以得到 Fnc(i j) =Ty(k)Wnj ,所以 、-nB =FnC =Fnc(i),c(i -1),c(i -2),.,c(i-n-1)=y(k),y(k)W：,., y(kW(n)k .n而 A = BFnH ,1 n1(k1)(j1)aij = - ；(i，j) =bikWn n yy(i) i=j=0 i = jn(J)k=-y(i)Wnn7綜上，命題得證。.編程實現(xiàn)1)存儲方案設(shè)計由于Toeplitz矩

5、陣與循環(huán)矩陣中有大量重復(fù)數(shù)據(jù)，因此，根據(jù)兩者的特點，我們用長度為2n-1的向量aT來儲存Toeplitz矩陣，即為=口1,.,3梟/1,.,3;用長度為n的向量出來儲存循環(huán)2一矩陣，即ac =Q,G,.,a。將儲存仝間從O(n )降低為O(n)。2)計算步驟設(shè)計.循環(huán)矩陣直接計算因為我們將循環(huán)矩陣存儲在長度為n的向量a。中，所以可以利用n次向量積求得答案，在計算時只需要注意在每次向量積時調(diào)整出中元素的順序即可。利用fft求解首先利用ifft函數(shù)求得循環(huán)矩陣的特征值， 注意到matlab中ifft函數(shù)并未做歸一化處理， 因此需 要在ifft得到的向量上乘以系數(shù) n；之后對向量X取ifft,并與特

6、征值做 matlab中的點乘操作，最后 對所得向量做fft運算即可得到結(jié)果。. Toeplitz 矩陣直接計算Toeplitz矩陣與向量乘法的步驟與循環(huán)矩陣類似，同樣需要注意所存儲向量在實際運算中元素的順序；利用fft求解Toeplitz矩陣與向量乘法的關(guān)鍵是將Toeplitz矩陣轉(zhuǎn)化為循環(huán)矩陣，只需要在存儲n階Toeplitz矩陣的向量中添加一個0,并調(diào)整元素順序即可得到2n階循環(huán)矩陣，之后的計算步驟與循環(huán)矩陣完全類似。3)編制Matlab程序Matlab程序已包含在文件夾中3.測試算例設(shè)計測試方法。進(jìn)行測試時由隨機(jī)數(shù)生成器生成矩陣為保證實驗的準(zhǔn)確性，我們要進(jìn)行多階多次實驗分別取平均時間以減

7、小實驗誤差，因此我們對1-1000階矩陣每階分別進(jìn)行 1000次實驗并作出平均用時-階數(shù)圖像來直觀地進(jìn)行對比。測試步驟：1、確定實驗矩陣的最小階數(shù)與最大階數(shù)2、循環(huán)矩陣C根據(jù)階數(shù)n生成n個隨機(jī)數(shù)，Toeplitz矩陣則要生成2n-1個隨機(jī)數(shù)，并根據(jù) 5）變?yōu)?循環(huán)矩陣3、對生成的循環(huán)矩陣分別直接計算向量積和利用FFT計算向量積，分別統(tǒng)計時間，并取時間平均值4、分別作出用時-階數(shù)圖像，進(jìn)行直觀對比1）給定循環(huán)矩陣C不變，對隨機(jī)生成的多個向量X連續(xù)計算矩陣-向量積Cx,并記錄總時間。對比直接計算與借助FFT計算的速度。519 U一 一 eO 1 5 O o O ono00直接計算實額循環(huán)矩陣計算時

8、間一矩陣階數(shù)圖像20030040050060。7003009001000矩陣階數(shù)nFFT計算實數(shù)循環(huán)矩陣計算時間一距陣階數(shù)圖像13-12003004005006007008009001000矩庠階數(shù)fi直接計算復(fù)數(shù)循環(huán)矩陣計篁時間一矩陣階數(shù)圖像2003004006006007003009001000矩陣階數(shù)nFFT計算復(fù)數(shù)循環(huán)矩陣計算時間一矩陣階數(shù)圖像02 O.5 19O1 5 o aH.OO600700S0090010009 O J20 0矩陣階數(shù)n可以看到在計算實數(shù)循環(huán)矩陣時，利用fft的方法相比于直接計算大約塊50100倍；而在計算復(fù)數(shù)矩陣時加速的效果則要差一些。2）給定Toeplitz

9、矩B$ T不變，對隨機(jī)生成的多個向量X連續(xù)計算矩陣-向量積Tx ,并記錄總時間。對比直接計算與借助 FFT計算的速度。x 10 1直接計算實繳T-pim福環(huán)矩陣 討算時間一走眸階數(shù) 圖像二IIIII001002003004005006007008009001000炬陣階數(shù)n1.5尸F(xiàn)T計算實數(shù)同汩匕循環(huán)矩陣計算時間一矩陣階數(shù)圖俅x 1 口.001004007Doo矩陣階弟n直接L算復(fù)數(shù)Tgpli泛循環(huán)矩陣計算時間一粗障階數(shù)圖像FFT計算復(fù)數(shù)Tqm心循環(huán)矩陣計算時間一矩陣階數(shù)圖像1002003004005006。07008009001000矩陣階數(shù)n利用fft計算Toeplitz的計算時間在512階矩陣左右處有一個跳變，這是因為在