第十三講并行算法設計(一)并行算法設計的評判典型的并行應用問題矩陣乘FOX算法二維FFT問題Laplace問題、Jacobi迭代法多體問題N-Body的Barnes-Hut算法一維FFT問題Gauss-Seidel迭代法求解線性方程組常用的三種并行算法設計方法理想并行(embarrassinglyparallelcomputation)分治并行(partitioninganddivide-and-conquerstrategies)流水并行(pipelinedcomputation):本講理想并行：矩陣乘FOX算法、二維FFT問題、Laplace問題、Jacobi迭代法下一講分治并行：一維FFT問題、多體問題N-Body的Barnes-Hut算法流水并行：Gauss-Seidel迭代法求解線性方程組并行算法設計的評判:什么樣的并行算法是好的什么樣的并行程序是好的?效率高:使用N個處理器時,加速比能夠接近N伸縮性好當N很大時,也能充分利用N個處理器的計算能力無論問題規模W多大,通過增加處理器的數量N,總是不需要修改程序就能讓問題可解、計算需要的時間可以控制如何提高效率?按照BSP模型每個超級計算步上,各個處理器承擔的負載要均衡要交換的數據量少超級計算步上不能有臨界區如果要有臨界區,就拆成多個超級計算步(例如SUM運算)如果不能用多個超級計算步來避免臨界區,應該考慮修改計算算法了(N-body計算)臨界區:概念上說,多個子任務要并發讀寫某個公共的數據結構,且每次只能有一個任務占有該數據結構的訪問權;從實現上看,訪問權的管理機制包括pthread模型:互斥鎖操作,lock/unlockCELL/BE模型:MAILBOX操作,由PPE維護公共數據結構的值MPI模型:SEND-RECV操作,由master維護公共數據結構的值BSP模型:臨界區與任務分配方法臨界區:子任務間的相關讀寫相關:按照BSP模型,這兩個子任務不在一個超級計算步上寫寫相關:這意味著子任務間寫的次序不影響結果,也就是規約操作(reduction,如SUM/MAX/MIN).如果不是規約操作引起的寫寫相關,計算結果就不確定了.引入臨時變量,保存每個子任務的局部結果多個超級計算步:第一個超級計算步計算各個子任務，后面N個超級計算步規約它們的結果,N是處理器數P的對數任務分配方法:并行程序實現階段考慮的事情,如何優化超級計算步的實現效率peer-to-peercollaborating(靜態任務分配):master-slave(動態任務分配):超級計算步上子任務的數量不能預先知道子任務的負載不一致使用的處理器性能不一致(以網絡環境為計算平臺時,會出現這種情況)求素數問題的任務分配實現:優化BSP的任務分配pthread使用臨界區:MPI使用master-slaveBSP:并行程序的效率高超級計算步上沒有臨界區REDUCTION操作引起的臨界區超級計算步上,各個處理器的負載要均衡:有足夠的子任務子任務負載不均衡、或者各個處理器的性能不一致時,通過master-slave分配任務,可以達到基本負載的基本平衡超級計算步上子任務復雜性遠大于它的數據通信開銷需要動態任務分配時,子任務的復雜性遠大于一次任務分配的開銷

pthread模型:一次臨界區操作CELL/BE模型:一次MAILBOX操作MPI模型:一次SEND-RECV操作BSP與計算機的體系結構基本沒有關系分布存儲計算機:消息傳遞模型實現數據交換、超級計算步的同步共享存儲計算機基于cache的多核處理器、SMP：鎖、信號量實現超級計算步的同步CELL處理器：DMA實現數據交換、mailbox實現同步BSP：并行程序的伸縮性好當N很大時,也能充分利用N個處理器的計算能力超級計算步上，要有足夠多的子任務無論問題規模W多大,通過增加處理器的數量N,總是不需要修改程序就能讓問題可解、計算需要的時間可以控制(Gustafson定律)超級計算步上，W增加了，子任務的數量也能增加，才能利用新增的處理器子任務的規模與W無關，否則處理器的尋址空間是個問題子任務的數據通信量與W無關如果有臨界區操作,子任務的復雜性遠超過臨界區操作的開銷規約操作引起的臨界區操作動態任務分配引起的臨界區操作并行程序要有好的伸縮性，需要面向分布存儲的計算機共享存儲計算機中CPU的數量很有限、總的存儲空間有限BSP：好的并行算法超級計算步上，有足夠多的子任務超級計算步上，或者沒有臨界區、或者只REDUCTION操作引起的臨界區超級計算步上子任務的數量隨著數據規模增加而增加子任務的數據規模、通信規模與總的數據規模沒有關系子任務的計算開銷遠大于數據傳輸開銷和臨界區操作開銷之和面向分布存儲結構的計算機并不排斥使用共享存儲結構的并行計算機：犧牲對問題規模的伸縮性好的并行算法是好的并行程序的基礎，但不意味著并行程序一定好實現代碼面向共享存儲結構的并行計算機消息交換中會在問題規模大時出現死鎖問題規模大時，消息交換出現死鎖示例if(myRank<comm_size-1)MPI_Send(buf1,buf_size,MPI_INT,myRand+1,tag,comm);elseMPI_Send(buf1,buf_size,MPI_INT,0,tag,comm);if(myRank>0)MPI_Recv(buf2,buf_size,MPI_INT,myRand-1,tag,comm,&status);elseMPI_Recv(buf2,buf_size,MPI_INT,comm_size-1,tag,comm,&status);每個處理器都要將buf1中的數據交給MPIRuntime后，才執行MPI_Recv語句MPIRuntime需要緩存了comm_size-1個處理器中的buf1數據后，才能開始執行MPI_Recv語句處理器數量越大，意味著可并行執行的子任務數量越多，需要緩存的數據越多，容易出現緩存空間不夠在SMP服務器上運行時，這種情況更容易出現：原來一臺4路的SMP沒有問題，新買一臺8路的SMP后就運行不了了解決辦法用非阻塞通信修改算法矩陣乘參考《AUser’sGuidetoMPI》第5.1節矩陣乘AB=CA是nm的矩陣B是mk的矩陣C是nk的矩陣并行性c(i,j)可以并行計算a(i,l)b(l,j)可以并行計算數據局部性：在計算c(i,j)的處理器上，需要A的第i行和B的第j列如何獲得高性能？在數據劃分與通信之間尋找trade-offFOX算法把處理器組織成一個二維的網格結構：

ppA、B、C都按照(block,block)劃分FOX算法：循環從每一行處理器中，選擇一個其中一個處理器，向所在行廣播A的本地局部片段curA用廣播得到的A片段curA與本地的B片段locB相乘，locC=locC+curAlocB把本地的B片段locB發送給“上鄰居”，把從“下鄰居”接收的數據作為本地的B片段locBA00A01A02A03A10A11A12A13A20A21A22A23A30A31A32A33B00B01B02B10B11B12B20B21B22B30B31B32C00C01C02C10C11C12C20C21C22C30C31C32=B03B13B23B33C03C13C23C33A00A01A02A03A10A11A12A13A20A21A22A23A30A31A32A33B00B01B02B10B11B12B20B21B32B30B31B32C00C01C02C10C11C12C20C21C22C30C31C32=B03B13B23B33C03C13C23C33A00A01A02A03A10A11A12A13A20A21A23A30A31A32A33B10B11B12B20B21B22B30B31B00B01B02C00C01C02C10C11C12C20C21C22C30C31C32B13B23B33B03C03C13C23C33A22B22=A00A01A02A03A10A11A12A13A20A21A22A23A30A31A32A33B20B21B22B30B31B32B00B01B12B10B11B12C00C01C02C10C11C12C20C21C22C30C31C32=B23B33B03B13C03C13C23C33A00A01A02A03A10A11A12A13A20A21A23A30A31A32A33B30B31B32B00B01B02B10B11B20B21B22C00C01C02C10C11C12C20C21C22C30C31C32B33B03B13B23C03C13C23C33A22B02=FOX算法分析FOX矩陣乘：代表了一種計算類型。有兩個對象集合A和BA的每個成員與B的每個成員都有聯系對A和B中的對象分組，每個作為一個子集：任取B中的一個子集，該子集分別與A的每個子集聯系特點：可伸縮性(scalability)增加處理器的數量，可以解決更大規模的問題：數據劃分單個處理器上，數組的大小有限制：內存空間、尋址空間增加處理器的數量，可以提高計算的效率計算開銷：nm2kCmultiple+n(m-1)kCadd通信開銷：p2CBStartup+nmCBmove+p2CPStartup+mkCPmovep2個處理器CBStartup：Broadcast通信的啟動開銷CBmove

：Broadcast通信中單個元素的開銷CPStartup：“點到點”通信的啟動開銷CBmove

：“點到點”通信中單個元素的開銷如何實現FOX算法？分布存儲結構的并行計算機處理器組需要組織成一個二維的拓撲結構需要把處理器劃分成子組，而且要有兩種劃分辦法二維拓撲結構的每一行作為一個子組：廣播A的數組片段二維拓撲結構的每一列作為一個子組：循環移動B的數組片段用“點到點”通信完全能夠滿足上述要求，只是編程比較麻煩MPI提供了一組通信子管理函數，這樣編程起來比較方便進程組拓撲結構通信子的創建進一步理解MPI的進程組、通信域、通信子進程組：參與計算的進程，無拓撲結構。在代碼設計階段，度量伸縮性通信域：參與計算的處理器、何時參與、及其拓撲結構。代碼實現階段，表達對運行平臺的要求通信子：進程組+通信子。代碼運行階段，滿足程序的要求不同階段，要不同的通信子多個并行程序的并發執行二維FFT什么是一維DiscreteFourierTransform：給定一個序列(a0,a1,…,an-1)，希望得到另一個序列(b0,b1,…,bn-1)，其中什么是二維DiscreteFourierTransform：給定一個nm的陣列(as,t)，希望得到另一個nm的陣列(cs,t)，其中

其中ei=cos+isine-i=cos-isin行變換列變換jkaj,kcj,kbj,k二維FFT分析計算特征：對nm的陣列(as,t)中每一行進行FourierTransform，得到一個新的nm陣列(bs,t)對nm的陣列(bs,t)中每一列進行FourierTransform，得到一個最終的nm陣列(cs,t)并行性從(as,t)nm變換到(bs,t)nm，每一行都是可以并行的數據局部性：計算bs,t時，涉及(as,t)nm的第s行全部元素從(bs,t)nm變換到(cs,t)nm，每一列都是可以并行的數據局部性：計算cs,t時，涉及(bs,t)nm的第t列全部元素同步：在開始計算(cs,t)nm前，必須完成(bs,t)nm的計算二維FFT的并行顯然，(按照BSP)二維FFT問題的并行策略進程的拓撲結構是一個線性序列計算任務劃分成兩個步驟步驟1：(as,t)nm按照(block,*)的方式劃分，每個進程負責局部(as,t)nm的行變換步驟2：(cs,t)nm按照(*,block)的方式劃分，每個進程負責局部(cs,t)nm的列計算通信：步驟1結束后，(bs,t)nm按照(block,*)的方式劃分步驟2要求(bs,t)nm按照(*,block)的方式劃分需要對(bs,t)nm進行一次轉置矩陣轉置A0B0C0D0E0F0A1B1C1D1E1F1A2B2C2D2E2F2A3B3C3D3E3F3A4B4C4D4E4F4A5B5C5D5E5F5A0A1A2A3A4A5B0B1B2B3B4B5C0C1C2C3C4C5D0D1D2D3D4D5E0E1E2E3E4E5F0F1F2F3F4F5進程的拓撲結構是一個線性序列原分布：A[N,M]按照(block,*)的方式劃分；各個進程中，局部的A片段按“列優先”方式存儲轉置在每個進程上，對自己的局部的A片段執行MPI_Scatter對每次MPI_Scatter所接收到的數組(Ai、Bi、Ci、Di、Ei、Fi)執行轉置二維FFT的特征總結符合BSP模型整個計算問題從時間軸上可以劃分成一組super-step，總的super-step數量是有限的在每個super-step上，有一組并行執行的子任務在相鄰的兩個super-step之間，數據訪問的局部性不同，需要通信實際上，任何一個應用問題都具備下列特征從時間軸上，可以把整個計算問題的執行劃分成有限個順序執行的階段(super-step)。例如在Laplace計算中，可把整個問題劃分成三個步驟分發初始數據執行迭代計算收集計算結果每個階段分別完成自己的任務，相對獨立一個階段的計算任務由一個進程完成：Laplace計算中的第一個階段計算任務(讀文件數據、組織成數組結構)和第二個階段計算任務(寫文件數據)都是由一個進程完成一個階段的計算任務由多個進程并行完成：Laplace計算中的第二個階段在相鄰的兩個super-step之間，數據訪問的局部性不同，需要通信理想并行什么是理想并行？(embarrassinglyparallelcomputation)計算任務可以劃分成一組子任務子任務之間是完全獨立的——不需要交換數據每個子任務處理一份數據——每個子任務處理的數據可以相同，也可以不同子任務之間不需要交換計算的結果/中間結果對子任務的結果不需要額外的處理(如分發distribute、收集collect、合并combine等)，就可以得到整個任務的結果舉例：二維FFT整個計算被劃分成兩個super-step，每個super-step上的計算任務都是理想并行的一個super-step可以劃分成一組子任務，每個子任務執行一行(列)數據的一維FFT在單個super-step上，每個一維FFT的子任務處理完全不同的數據，沒有任何的數據交換FOX算法也是理想并行FOX算法：循環從每一行處理器中，選擇一個其中一個處理器，向所在行廣播A的本地局部片段curA用廣播得到的A片段curA與本地的B片段locB相乘，locC=locC+curAlocB把本地的B片段locB發送給“上鄰居”，把從“下鄰居”接收的數據作為本地的B片段locB在每個super-step上,要改變一次數據的劃分輸入數據計算結果進程Laplace問題典型的“場”問題：連續的數據區域（datadomain），每一點都有若干感興趣的、隨時間變化的物理量（速度，溫度，壓力等），且通常它們在時間和空間上都是連續的。“標準的”處理方法：在空間和時間上做離散化處理。

Laplace問題特征并行性：在時間軸上，每個時間步串行執行在空間軸上，數據區域上的每一點都是可以并行的數據局部性：在局域范圍內，只與“鄰居”有數據交換關系按照BSP，我們將每個時間步看作一個super-step每個super-step上的通信模式相同每個super-step上的計算劃分相同每個super-step上的數據劃分相同Jacobi迭代法一種線性方程組的數值解法什么是Jacobi迭代法：AX=B表示一個方程組A是mn的矩陣，且

A(i,i)0X是長度n的未知向量B是長度為m的已知向量令X(0)=(c0,c1,…,cn-1)當max(|X[i](t)-X[i](t-1)|)<時，X=X(t)Jacobi迭代法的理解在工程計算和科學分析領域，經常需要對系統的穩定性進行分析：在一定的約束條件下，我們關系的物理量在整個問題空間中的分布態勢。例如：在考慮三峽大壩時，設計者必須關心的一件事情是三峽壩區的應力場分布，即這個區域中各個地質斷層板塊之間的穩定性大壩和蓄水都會增加有關地質斷層板塊上的壓力在哪里建壩、建多大的壩、蓄水位最大可到多高時，這種新的壓力不至于引起原有地質斷層板塊的穩定性我們不知道地質斷層板塊之間的作用力，但我們知道有哪些地質斷層板塊地質斷層板塊之間的拓撲結構山體的高度、各個地質斷層板塊的剛度、……人類已有的認識可以讓我們確定一組約束條件：為了保持地質斷層板塊的穩定性，它們之間的各種作用力(壓力、摩擦力等)與地質斷層板塊的特征(大小、剛度、位置、密度等)一定滿足某個關系我們需要根據這些約束關系推算作用力在整個壩區的分布，從數學模型上考慮，可以把每個作用力看作是該系統中的一個自由度X的每個分量表示一個自由度A和B表示約束條件Jacobi迭代法的特征在每個時間步上，執行的計算任務完全相同并行性X[0](t)、X[1](t)、…、X[m-1](t)的計算可以并行進行數據局部性：計算X[i](t)

，涉及A的第i行和B[i]、X(t-1)數據相關性：每個X[i](t)都需要X(t-1)并行算法進程的拓撲結構是一個線性序列A按照(block,*)方式劃分、B按照“block”方式劃分對X(t-1)進行數據復制X(t)按照“block”方式劃分通信：維護X(t-1)按照BSP模型分析把每個時間步上的計算看作一個super-step在單個super-step上時間的scalability：每個處理器有一次通信，執行n/p個X[i](t)的計算規模的scalability：A和B都被劃分，只有X被復制運行的速度取決于問題的收斂速度——循環迭代的次數近似理想并行什么是近似理想并行？(nearlyembarrassinglyparallelcomputation)計算任務可以劃分成一組子任務子任務之間是完全獨立的——不需要交換數據每個子任務處理一份數據——每個子任務處理的數據可以相同，也可以不同子任務之間不需要交換計算的結果/中間結果經過對子任務結果額外的處理(如分發distribute、收集collect、合并combine等)，得到整個任務的結果舉例：Jacobi迭代(求解方程組AX=B)，在每個時間步上的計算任務都是近似理想并行的每個X[i](t)的計算只涉及A的第