并行計(jì)算MPI程序設(shè)計(jì)

1、*實(shí)踐教學(xué)*蘭州理工大學(xué)理學(xué)院 2016年春季學(xué)期并行計(jì)算課程設(shè)計(jì)專業(yè)班級(jí)： 2013級(jí)信息與計(jì)算科學(xué)姓名： 學(xué)號(hào)：13550418 指導(dǎo)教師： 成績：摘要FFT，即為快速傅氏變換，是離散傅氏變換的快速算法，它是根據(jù)離散傅氏變換的奇、偶、虛、實(shí)等特性，對(duì)離散傅立葉變換的算法進(jìn)行改進(jìn)獲得的。它對(duì)傅氏變換的理論并沒有新的發(fā)現(xiàn)，但是對(duì)于在計(jì)算機(jī)系統(tǒng)或者說數(shù)字系統(tǒng)中應(yīng)用離散傅立葉變換，可以說是進(jìn)了一大步。   設(shè)x(n)為N項(xiàng)的復(fù)數(shù)序列，由DFT變換，任一X（m）的計(jì)算都需要N次復(fù)數(shù)乘法和N-1次復(fù)數(shù)加法，而一次復(fù)數(shù)乘法等于四次實(shí)數(shù)乘法和兩次實(shí)數(shù)加法，一次復(fù)數(shù)加法等于兩次

2、實(shí)數(shù)加法，即使把一次復(fù)數(shù)乘法和一次復(fù)數(shù)加法定義成一次“運(yùn)算”（四次實(shí)數(shù)乘法和四次實(shí)數(shù)加法），那么求出N項(xiàng)復(fù)數(shù)序列的X（m）,即N點(diǎn)DFT變換大約就需要N2次運(yùn)算。當(dāng)N=1024點(diǎn)甚至更多的時(shí)候，需要N2=1048576次運(yùn)算，在FFT中，利用WN的周期性和對(duì)稱性，把一個(gè)N項(xiàng)序列（設(shè)N=2k,k為正整數(shù)），分為兩個(gè)N/2項(xiàng)的子序列，每個(gè)N/2點(diǎn)DFT變換需要（N/2）2次運(yùn)算，再用N次運(yùn)算把兩個(gè)N/2點(diǎn)的DFT變換組合成一個(gè)N點(diǎn)的DFT變換。這樣變換以后，總的運(yùn)算次數(shù)就變成N+2（N/2）2=N+N2/2。繼續(xù)上面的例子，N=1024時(shí)，總的運(yùn)算次數(shù)就變成了525312次，節(jié)省了大約50%的運(yùn)算

3、量。而如果我們將這種“一分為二”的思想不斷進(jìn)行下去，直到分成兩兩一組的DFT運(yùn)算單元，那么N點(diǎn)的DFT變換就只需要Nlog(2)(N)次的運(yùn)算，N在1024點(diǎn)時(shí)，運(yùn)算量僅有10240次，是先前的直接算法的1%，點(diǎn)數(shù)越多，運(yùn)算量的節(jié)約就越大，這就是FFT的優(yōu)越性關(guān)鍵字：FFT 蝶式計(jì)算 傅里葉變換目錄摘要2目錄3一、題目及要求41.1題目41.2要求4二、算法設(shè)計(jì)與算法原理52.1 算法原理與設(shè)計(jì)52.2設(shè)計(jì)求解步驟6三、算法描述與算法流程73.1算法描述73.2 流程圖9四、源程序代碼與運(yùn)行結(jié)果104.1源程序104.2運(yùn)行結(jié)果16五、算法分析及其優(yōu)缺點(diǎn)175.1 算法分析175.2優(yōu)缺點(diǎn)18

4、六、總結(jié)19七、參考文獻(xiàn)20一、題目及要求1.1題目對(duì)給定的=（1,2，4，3，8，6，7，2），利用串行FFT遞歸算法(蝶式遞歸計(jì)算原理)計(jì)算其傅里葉變換的結(jié)果1.2要求利用串行遞歸與蝶式遞歸原理，對(duì)給定的向量求解傅里葉變換的結(jié)果二、算法設(shè)計(jì)與算法原理2.1 算法原理與設(shè)計(jì) 令 為n/2次單位元根，則有 . 將b向量的偶數(shù)項(xiàng)和奇數(shù)項(xiàng)分別記為 和 注意推導(dǎo)中反復(fù)使用圖2.12.2設(shè)計(jì)求解步驟三、算法描述與算法流程3.1算法描述n=8的FFT蝶式計(jì)算圖圖3.1圖3.2 FFT遞歸計(jì)算流程圖3.2 流程圖開始計(jì)算出前size_x/2個(gè)exp(-j*2*k/size_x)個(gè)值，即W的值輸入序列對(duì)應(yīng)值

5、（例如5+j3，輸入5 3）輸入序列長度size_x飛級(jí)數(shù)i>=？級(jí)數(shù)i加1 是 輸出fft結(jié)果序列結(jié)束 否 該級(jí)該組起始下標(biāo)j>=？計(jì)算出該級(jí)需要的W的個(gè)數(shù)l 是 否組起始下標(biāo)加2*l該級(jí)該組元素序數(shù)k>=？K加1 是Xj+k Xj+klXj+k+l*W(size_x/2/l)*k Xj+k+l -1 否圖3.3四、源程序代碼與運(yùn)行結(jié)果4.1源程序/*FFT*/ /整個(gè)程序輸入和輸出利用同一個(gè)空間xN，節(jié)約空間 #include <stdio.h> #include <math.h> #include <stdlib.h> #define

6、 N 1000 /定義輸入或者輸出空間的最大長度typedefstruct double real;doubleimg; complex; /定義復(fù)數(shù)型變量的結(jié)構(gòu)體 void fft(); /快速傅里葉變換函數(shù)聲明 void initW(); /計(jì)算W(0)W(size_x-1)的值函數(shù)聲明 void change(); /碼元位置倒置函數(shù)函數(shù)聲明 void add(complex,complex,complex *); /*復(fù)數(shù)加法*/ void mul(complex,complex,complex *); /*復(fù)數(shù)乘法*/ void sub(complex,complex,complex

7、 *); /*復(fù)數(shù)減法*/ void divi(complex,complex,complex *); /*復(fù)數(shù)除法*/ void output(); /*輸出結(jié)果*/ complex xN,*W; /*輸出序列的值*/intsize_x=0; /*輸入序列的長度，只限2的N次方*/ double PI; /pi的值int main() inti;system("cls"); PI=atan(1)*4;printf("Please input the size of x:n"); /*輸入序列的長度*/scanf("%d",&

8、size_x);printf("Please input the data in xN:(such as:5 6)n"); /*輸入序列對(duì)應(yīng)的值*/for(i=0;i<size_x;i+)scanf("%lf %lf",&xi.real,&xi.img);initW(); /計(jì)算W(0)W(size_x-1)的值 fft(); /利用fft快速算法進(jìn)行DFT變化 output(); /順序輸出size_x個(gè)fft的結(jié)果return 0; /*進(jìn)行基-2 FFT運(yùn)算,蝶形算法。這個(gè)算法的思路就是， 先把計(jì)算過程分為log(size_x

9、)/log(2)-1級(jí)（用i控制級(jí)數(shù)）； 然后把每一級(jí)蝶形單元分組（用j控制組的第一個(gè)元素起始下標(biāo)）； 最后算出某一級(jí)某一組每一個(gè)蝶形單元（用k控制個(gè)數(shù)，共l個(gè)）。*/voidfft() inti=0,j=0,k=0,l=0; complexup,down,product; change(); /實(shí)現(xiàn)對(duì)碼位的倒置 for(i=0;i<log(size_x)/log(2);i+) /循環(huán)算出fft的結(jié)果 l=1<<i; for(j=0;j<size_x;j+=2*l) /算出第m=i級(jí)的結(jié)果【i從0到（log(size_x)/log(2)）-1】 for(k=0;k<

10、;l;k+) /算出第i級(jí)內(nèi)j組蝶形單元的結(jié)果 /算出j組中第k個(gè)蝶形單元mul(xj+k+l,W(size_x/2/l)*k,&product); /*size/2/l是該級(jí)W的相鄰上標(biāo)差，l是該級(jí)該組取的W總個(gè)數(shù)*/add(xj+k,product,&up); sub(xj+k,product,&down); xj+k=up; /up為蝶形單元右上方的值 xj+k+l=down; /down為蝶形單元右下方的值 void initW() /計(jì)算W的實(shí)現(xiàn)函數(shù) inti; W=(complex *)malloc(sizeof(complex) * size_x); /*

11、申請(qǐng)size_x個(gè)復(fù)數(shù)W的空間(這部申請(qǐng)的空間有點(diǎn)多，實(shí)際上只要申請(qǐng)size_x/2個(gè)即可)*/ for(i=0;i<(size_x/2);i+) /*預(yù)先計(jì)算出size_x/2個(gè)W的值，存放,由于蝶形算法只需要前size_x/2個(gè)值即可*/ Wi.real=cos(2*PI/size_x*i); /計(jì)算W的實(shí)部 Wi.img=-1*sin(2*PI/size_x*i); /計(jì)算W的虛部 void change() /輸入的碼組碼位倒置實(shí)現(xiàn)函數(shù) complex temp; unsigned short i=0,j=0,k=0; double t; for(i=0;i<size_x;

12、i+) k=i; j=0; t=(log(size_x)/log(2); while(t-)>0) j=j<<1; j|=(k & 1); k=k>>1; if(j>i) temp=xi; xi=xj; xj=temp; void output() /輸出結(jié)果實(shí)現(xiàn)函數(shù) inti; printf("The result are as followsn"); for(i=0;i<size_x;i+) printf("%.4f",xi.real); /輸出實(shí)部 if(xi.img>=0.0001) /如果

13、虛部的值大于0.0001，輸出+jx.img的形式printf("+j%.4fn",xi.img); else if(fabs(xi.img)<0.0001) /如果虛部的絕對(duì)值小于0.0001，無需輸出 printf("n"); elseprintf("-j%.4fn",fabs(xi.img); /如果虛部的值小于-0.0001，輸出-jx.img的形式 void add(complex a,complexb,complex *c) /復(fù)數(shù)加法實(shí)現(xiàn)函數(shù) c->real = a.real + b.real; /復(fù)數(shù)實(shí)部相

14、加 c->img = a.img + b.img; /復(fù)數(shù)虛部相加 void mul(complex a,complexb,complex *c) /復(fù)數(shù)乘法實(shí)現(xiàn)函數(shù) c->real = a.real*b.real - a.img*b.img; /獲取相乘結(jié)果的實(shí)部 c->img = a.real*b.img + a.img*b.real; /獲取相乘結(jié)果的虛部 void sub(complex a,complexb,complex *c) /復(fù)數(shù)減法實(shí)現(xiàn)函數(shù) c->real = a.real - b.real; /復(fù)數(shù)實(shí)部相減 c->img = a.img -

15、b.img; /復(fù)數(shù)虛部相減 void divi(complex a,complexb,complex *c) /復(fù)數(shù)除法實(shí)現(xiàn)函數(shù) c->real = (a.real*b.real + a.img*b.img) / (b.real*b.real+b.img*b.img); /獲取相除結(jié)果的實(shí)部 c->img = (a.img*b.real - a.real*b.img) / (b.real*b.real+b.img*b.img); /獲取相除結(jié)果的虛部 4.2運(yùn)行結(jié)果（1）處理器p=4圖4.2.1（2）處理器p=2圖4.2.2五、算法分析及其優(yōu)缺點(diǎn)5.1 算法分析（1）FFT算法的

16、基本原理是把長序列的DFT逐次分解為較短序列的DFT。按照抽取方式的不同可分為DIT-FFT（按時(shí)間抽?。┖虳IF-FFT（按頻率抽?。┧惴?。按照蝶形運(yùn)算的構(gòu)成不同可分為基2、基4、基8以及任意因子（2n,n為大于1的整數(shù)），基2、基4算法較為常用。（2）總體結(jié)構(gòu)說明 輸入數(shù)據(jù)為串行的數(shù)據(jù)流，故在第一級(jí)蝶形運(yùn)算模塊前加入串并轉(zhuǎn)換模塊，將串行數(shù)據(jù)流轉(zhuǎn)換為并行的兩列數(shù)據(jù)流以適應(yīng)基2蝶形運(yùn)算模塊的輸入信號(hào)要求。 由于每級(jí)蝶形運(yùn)算一次處理的兩個(gè)輸入數(shù)據(jù)不能直接由前一級(jí)蝶形運(yùn)算一次性輸出，故在兩個(gè)蝶形運(yùn)算單元之間插入延時(shí)對(duì)齊模塊，將前一級(jí)蝶形運(yùn)算的結(jié)果（兩列并行的數(shù)據(jù)流）作適當(dāng)?shù)难訒r(shí)

17、并通過轉(zhuǎn)接器對(duì)齊，形成后一級(jí)蝶形運(yùn)算模塊所需要的2列輸入序列。 在最后一級(jí)蝶形運(yùn)算后加入串并轉(zhuǎn)換模塊，將2列并行的數(shù)據(jù)流合成為1列。最后加入倒序模塊將DIF-FFT得到的倒序輸出序列整理為順序輸出。 旋轉(zhuǎn)因子產(chǎn)生模塊產(chǎn)生各級(jí)基2蝶形運(yùn)算所需的旋轉(zhuǎn)因子。由運(yùn)算流圖可以看出最后一級(jí)的旋轉(zhuǎn)因子其實(shí)是1，故可省略最后一級(jí)蝶形運(yùn)算單元中的旋轉(zhuǎn)因子乘法器。因此用一個(gè)雙口ROM將兩組數(shù)據(jù)分別輸出到第一級(jí)和第二級(jí)的蝶形運(yùn)算單元即可。 基2蝶形運(yùn)算模塊由兩個(gè)復(fù)數(shù)加法器和一個(gè)復(fù)數(shù)乘法器構(gòu)成。旋轉(zhuǎn)因子由ROM產(chǎn)生后，作為復(fù)數(shù)乘法器的輸入之一，與前面復(fù)數(shù)加法器得到的結(jié)果相乘完成一次蝶形運(yùn)

18、算。為提高系統(tǒng)的運(yùn)行速度可在蝶形運(yùn)算單元中插入流水線寄存中間結(jié)果（3）蝶形運(yùn)算單元如下所示：5.2優(yōu)缺點(diǎn)（1）優(yōu)點(diǎn)：結(jié)構(gòu)清晰，可讀性強(qiáng)，而且容易用數(shù)學(xué)歸納法來證明算法的正確性，因此它為設(shè)計(jì)算法、調(diào)試程序帶來很大方便。 （2）缺點(diǎn)：遞歸算法的運(yùn)行效率較低，無論是耗費(fèi)的計(jì)算時(shí)間還是占用的存儲(chǔ)空間都比非遞歸算法要多。六、總結(jié)通過這次課程設(shè)計(jì)，讓我更加深刻了解課本知識(shí)，和以往對(duì)知識(shí)的疏忽得以補(bǔ)充。雖然這次課程是那么短暫的1周時(shí)間，我感覺到這些天我的所學(xué)勝過我這一學(xué)期所學(xué)，這次任務(wù)原則上是設(shè)計(jì)，其實(shí)就是一次大的作業(yè)，是讓我對(duì)課本知識(shí)的鞏固和對(duì)基本公式的熟悉和應(yīng)用課程設(shè)計(jì)是一個(gè)重要的教學(xué)環(huán)節(jié)，通過課程設(shè)計(jì)使我們了解到一些實(shí)際與理論之間的差異。通過課程設(shè)計(jì)不僅可以鞏固專業(yè)知識(shí)，為以后的工作打下了堅(jiān)實(shí)的基礎(chǔ)，而其還可以培養(yǎng)和熟練使用資料，運(yùn)用工具書的能力，把我們所學(xué)的課本知識(shí)與實(shí)踐結(jié)合起來，起到溫故而知新的作用。課程設(shè)計(jì)誠然是一門專業(yè)課，給我很多專業(yè)知識(shí)以及專業(yè)技能上的提升，同時(shí)又是一門講道課，一門