快速多極子方法的并行技術(shù).ppt

1、快速多極子方法的并行技術(shù),國家973項目高性能科學(xué)計算研究 大規(guī)模并行計算研究,綱要,FMM Data Structures Parallelization,綱要,FMM Data Structures Parallelization,FMM in Computational Electromagnetics,EFIE MFIE CFIE Green函數(shù),積分方程的離散,RWG矢量基函數(shù) MOM 離散,Rao-Wilson-Glisson,Method of Moments,Surface is Discretized into Patches (Basis Functions),Basis

2、Functions Interact through the Greens Function,Generates a Dense Method of Moments Matrix,線性系統(tǒng)：,Mx=s M是NN矩陣，x、s是N矢量 Direct solution (Gauss elimination,LU Decomposition,SVD,) 空間復(fù)雜度為O(N2) ，需要O(N3)次運算； Iterative methods,空間復(fù)雜度仍為O(N2)，如果K(k largest N =32,768 Finding:快速矩陣乘向量的算法(NlogN)； 并行實施。,Fast Multipol

3、e Methods(FMM),Introduced by Rokhlin 需要構(gòu)建的函數(shù)，如Parent(n), ChildrenAll(n), Children(X;n,l), NeighborsAll(n,l), Neighbors(X;n,l),Grouping,每個盒子在l=0,.,L層的指標(biāo)n=Number=0,1,2ld-1.,由于E2(n,l)和E3(n,l)是互補的，因此對任意的n,l,綱要,FMM Data Structures Parallelization,Data Structure,構(gòu)造樹 壓縮八叉樹 建立連接 morton鍵 排序,構(gòu)造樹,離散點的空間層次劃分,對應(yīng)

4、的四叉樹及其壓縮四叉樹,確定層數(shù)L 根據(jù)讀入模型的最大幾何尺寸確定計算區(qū)域D,根據(jù)問題的參數(shù)確定最小盒子 尺寸d ,樹結(jié)構(gòu)的層數(shù)為L=log2(D/d) 第l-1層立方體等分為八個子立方體,形成第l層更小的立方體,l-1是l層的父層,l層是l-1層的子層. 形成相鄰組、次相鄰組、遠場組 第l層的節(jié)點數(shù)不超過2dl個,構(gòu)造樹八叉樹(1),構(gòu)造樹八叉樹（2）,procedure Octree_Build Octree = empty For i = 1 to n . 對所有的點做循環(huán) Octree_Insert(i, root) . 將點i插入八叉樹相應(yīng)的位置 end For . 八叉樹中可能有很

5、多空的葉節(jié)點, 但它們的兄弟節(jié)點非空 Traverse the tree (via, say, breadth first search) . 寬度周游 Eliminating empty leaves . 去掉空的葉節(jié)點 Compress Octree . 壓縮八叉樹 . 如果某中間節(jié)點僅包含一個子節(jié)點，則被其替換 每個節(jié)點用兩個鍵值描述: 包含相同數(shù)據(jù)的最小單元和最大單元,構(gòu)造樹八叉樹（3）,包含N個葉節(jié)點的壓縮八叉樹節(jié)點總數(shù)不超過2N-1 因此可以采用數(shù)組而不是鏈表來存儲壓縮八叉樹 MLFMM基于后序周游的壓縮八叉樹數(shù)據(jù)結(jié)構(gòu) 從鍵值最小的葉節(jié)點開始周游 每個節(jié)點都存儲在其子節(jié)點之后，且緊

6、挨其子節(jié)點存儲 節(jié)點排序時，使用的是其所表示的最小單元鍵值 對于兄弟節(jié)點，按鍵值從小到大排序 各節(jié)點所表示的單元指的是它所包含的最小單元 后序周游存儲方式是實現(xiàn)與分布無關(guān)的自動負載平衡并行MLFMM的關(guān)鍵 分布無關(guān)自適應(yīng)樹(Distribution-Independent Adaptive Tree, DIAT) 構(gòu)造2d維DIAT的復(fù)雜度為O(dnlogn),樹 節(jié) 點 存 儲,Morton鍵,為什么不用指針？ 用指針指向Children的指針可以很方便的建立父子節(jié)點之間的關(guān)系，從而構(gòu)成樹結(jié)構(gòu)的拓撲結(jié)構(gòu)。在串行程序，指針可以在全局存儲空間中尋址，效率很高也很準(zhǔn)確。但在內(nèi)存分布式并行環(huán)境中，一

7、個計算節(jié)點不能直接訪問另一個計算節(jié)點上的存儲空間，因此用于聯(lián)系樹結(jié)構(gòu)拓撲結(jié)構(gòu)的指針只能在其所在的計算節(jié)點上才有意義，如果要讓指針?biāo)赶虻臉涔?jié)點能夠存儲在其他節(jié)點上，就必須小心處理指針的變換關(guān)系。以便將指針的指向正確地映射到其他計算節(jié)點上的內(nèi)存空間。這種轉(zhuǎn)換，使得基于指針的樹拓撲方法在分布式并行環(huán)境中實現(xiàn)起來不僅很復(fù)雜，而且效率也不高。 Morton鍵技術(shù)是實現(xiàn)并行多層快速多極子的關(guān)鍵技術(shù)之一！ 原理：將空間坐標(biāo)的二進制序列按位交叉，把空間中某一點在x、y、z方向的坐標(biāo)/序號映射為一個值，這個值就是morton鍵值。,Morton次序,位于第m層,在該層中編號為n的盒子, 一般采用Morton次

8、序編號為(n, m). 左圖為第三層的Morton次序 ， 右圖為每一層編號,前三層分別有1, 4, 16個點.,Morton編碼,小的灰盒子在第3層,本層編號為11, 于是其Morton編號為(11,3); (023)4=(11)10 采用四進制編號為(0,2,3); Morton編號(Num,l)， 在2d叉樹中可以得到某盒子對應(yīng)的 2d進制編號: (N1, N2,Nl)2d,再按照下面的公式計算其Morton編號：,算法,空間編碼,尺度轉(zhuǎn)換 二進制編碼 d維空間編碼 二進制交錯及其解交錯 涉及到的算法：查找Parent、Children、Neighbor、盒子中心坐標(biāo)、盒子編碼等,尺度轉(zhuǎn)

9、換,2d樹結(jié)構(gòu)每個參數(shù)都有自己的參數(shù)Dj映射原始的盒子的大小 x1,min, x1,max x1,min, x1,max ,xj,max-xj,min=Dj,j=1,d 把原始的盒子映射到單位盒子 上 實際的三維物理空間Dj=D， xmin=(x1,min,，xd,min), 稱x為真實的坐標(biāo)， 為該點標(biāo)準(zhǔn)化的坐標(biāo) 目的：實施二進制編碼轉(zhuǎn)換，方便查找！,二進制編碼,For example d=1: 用十進制表示為 對于 不僅可以表示為 也可以表示為 用二進制表示 為： 對于,位交錯,在d維單位盒子里, 點坐標(biāo) 其中第k個分量 也可以寫成交錯二進制(bit interleaving)的形式:,解

10、交錯,d維盒子在第l層的Morton鍵值為 可以解交錯(bit deinterleaving)為d個數(shù)值代表該盒子的d維坐標(biāo),Number=7689310=,=1112=710,=1110102=5810,=10012=910,(9, 58, 7),Number3=,Number2=,Number1=,尋找給定點所在的盒子,在d維單位盒子里, 點坐標(biāo)的交錯2d進制形式: 可以利用2d進制移位取整尋找給定點在第l層所在的盒子: 或者寫成dl位的平移形式:,查找分為5層八叉樹結(jié)構(gòu)的點(0.7681, 0.0459,0.3912)在第3層盒子的號 (1)二進制轉(zhuǎn)化(0.11000,0.00001,0

11、.01100)2 (2）形成單個串0.1001010010000102 (3）3*3=9 (Number,3)=1001010012=297 3*5=15 (Number,5)=1001010010000102=19010,查找盒子中心,單位盒子第l層編號為Num的子盒子中心的二進制d維坐標(biāo)形式: 或者不依賴計數(shù)體系,寫為: 例如: 八叉樹第5層10進制編號為533的盒子，計算過程為: (1)將Morton編號轉(zhuǎn)為2進制: 53310=1, 000, 010, 1012 (2)通過解交錯得到該盒子的三維坐標(biāo): Num3=10012=910 , Num2=0102=210 , Num1=0012

12、=110 (3)計算盒子中心坐標(biāo): x3c(533, 5)=2-5(9+0.5)=0.296875; x2c(533, 5)=2-5(2+0.5)=0.078125; x1c(533, 5)=2-5(1+0.5)=0.04875; 因此xc=(0.04875, 0.078125, 0.296875),算法,父盒子編碼,計算某盒子的父盒子 父盒子編碼與層次無關(guān), 其計算方法: 除法的商取整, 比如11/4=2，于是 例如： 于是#5981的父盒子是#747 除以2d相當(dāng)于作d次2進制位移就可以了,算法,子盒子編碼,計算某盒子的子盒子 子盒子是個集合,與所在的層無關(guān): 其計算方法為: For ex

13、ample,查找鄰居盒子（1）,基于交錯二進制的鄰居查詢算法 Step1.解交錯(deinterleaving).十進制(或二進制)編號可轉(zhuǎn)為d維坐標(biāo): Step2.每個分坐標(biāo)的鄰居: 于是鄰居集為: 其中nk定義為,比如編號為#26的 2維盒子,其鄰居 查詢過程如下: 2610=(01 10 10)2 解交錯形式為 (011,100) 2=(3,4) 10 其相鄰盒子為 (2,3),(2,4), (2,5) ; (3,3),(3,5); (4,3),(4,4),(4,5); 8個點的二進制: (010,011) 2 , (010,100) 2 , (100,101) 2,相鄰盒子編號的交錯

14、(interleaving)二進制形式: 001101, 011000, 011001, 001111, 011011, 100101, 110000, 110001 十進制編號為:13,24,25, 15,27, 37,48,49.,查找鄰居盒子（2）,綱要,FMM Data Structures Parallelization,并行實施技術(shù)(parallel implementation),MLFMM具有很好的并行效率 大多數(shù)操作都在本地機上進行 例如，Parant to Children,box to its interaction list ,八叉樹結(jié)構(gòu)是很大的 需要生成和分布式存儲

15、每個處理器只保持本地的基本樹 大多數(shù)工作只在本地的基本樹上進行 數(shù)據(jù)需要同步 父節(jié)點需要子節(jié)點上的值，但這兩個節(jié)點在不同的處理器上 荷載平衡問題,但是，還存在：,分布式八叉樹 負載平衡 相互作用表列 相鄰結(jié)點的通信 次相鄰點的通信,需要解決,并行計算步驟,構(gòu)造壓縮八叉樹 近場矩陣計算 建立轉(zhuǎn)移節(jié)點列表 遠場矩陣計算 聚合 轉(zhuǎn)移 發(fā)散,樹結(jié)構(gòu)的并行劃分（1）,二維計算區(qū)域?qū)?yīng)的分布式四叉樹,樹結(jié)構(gòu)的并行劃分（2）,構(gòu)造分布式壓縮八叉樹(1),層數(shù)L=log2(D/d), D和d為計算區(qū)域劃分的最大和最小盒子尺寸 將n個基函數(shù)在p個處理機上按次序平均分配，再按照基函數(shù)的位置生成包含這些基函數(shù)的葉節(jié)

16、點 由于不同基函數(shù)可包含于同一葉節(jié)點，因此這樣的葉節(jié)點會同時存儲在不同處理機上 RWG基函數(shù)定義在各邊上，并包含一對完整的三角形,P1,P2,P3,A1,A2,A3,用邊的中點代表整個邊, 各點都有相應(yīng)的最底層 非空盒子, 即葉節(jié)點.將 葉節(jié)點的Morton鍵值賦給中點所在的邊,構(gòu)造分布式壓縮八叉樹(2),采用并行排序算法對所有處理機中的基函數(shù)和葉節(jié)點排序，使個處理機包含相同數(shù)量的基函數(shù) 對每個處理機里的N/p個鍵值，采用快速排序(quicksort) 全局并行排序采用取樣排序(samplesort), 它用到位排序(bitornic sort) 排序時用到的通信為MPI_Allgather和

17、MPI_Alltoall 每個處理機還包含下一處理機的第一個葉節(jié)點, 并根據(jù)這些葉節(jié)點建立本地壓縮八叉樹，通過后序周游存儲 各處理機將本地壓縮樹中位于從下一處理機得到的葉節(jié)點之后的節(jié)點發(fā)送到下一處理機，并按后序周游插到對應(yīng)的位置 共享葉節(jié)點個數(shù)不超過L, 每個處理機接收的非本地節(jié)點不超過7L 各處理機從下而上構(gòu)造本地樹的復(fù)雜度為O( (N/p)log(N/p),近場計算,Morton次序保證兄弟節(jié)點的鍵值是相鄰的，但并不是只有兄弟節(jié)點相鄰, 因此需要調(diào)用位交錯和解交錯函數(shù)去查找鄰居節(jié)點 近場矩陣Anear只需要考慮最細層(第L層)每個節(jié)點及其相鄰節(jié)點所包含的基函數(shù)相互作用; 必須按照MOM計算

18、, 并在迭代前存儲 每個葉節(jié)點最多有26個相鄰節(jié)點。若最細層每個盒子至多包含c個基函數(shù), 則每個葉節(jié)點上的計算量為27c2 如果同一棵子樹上的相鄰節(jié)點位于同一處理機, 則無需通信 如果相鄰節(jié)點位于不同的處理機上, 則 使用MPI_Allgather獲得每一處理機的第一個和最后一個葉節(jié)點的鍵值.,并存儲在長為2p的數(shù)組中, 通過該數(shù)組的檢索得到相鄰葉節(jié)點 調(diào)用MPI_Alltoall實現(xiàn)數(shù)據(jù)(葉結(jié)點及其包含的基函數(shù))的分發(fā); 整個近場計算復(fù)雜度為O( (N/p)log(N/p),遠場計算,局部的不變項(如最細層的D和A)只需要分配到它對應(yīng)的子樹所在處理機上，全局不變項(如常數(shù))則分配到所有處理機

19、上；由于下一層的計算依賴于上一層的信息，因此完成每一層的計算時，各個處理機都需要進行一次同步 存儲時采用內(nèi)存循環(huán)策略，它依賴于數(shù)據(jù)的相關(guān)性。聚集項S在層層上聚時分配內(nèi)存，每當(dāng)層層下推時某層的發(fā)散項B計算完畢，就將該層的S的內(nèi)存釋放掉；發(fā)散項B在層層下推時分配內(nèi)存, 每當(dāng)處理機上某層的所有的B計算完畢，就將其父層的B釋放掉 歸并各處理機得到的結(jié)果, 即為遠場矩陣向量乘 通過歸約得到整個系數(shù)矩陣與向量的乘積 通過并行的迭代法得到計算結(jié)果(等效電流), 每次迭代的矩陣向量乘法計算完成時，需要進行一次同步 完成計算結(jié)果的后處理，比如RCS的計算。,樹結(jié)構(gòu)代碼,上聚 A M2M 內(nèi)插值 下推 C M2L

20、 B L2L 外插值 二叉樹的例子,建立相互作用表列,相互作用表列(interaction list)包含每一層、每個節(jié)點的次相鄰節(jié)點 次相鄰節(jié)點指它們本身不相鄰，但它們的父節(jié)點相鄰 因此每個節(jié)點最多有63-33=189個次相鄰節(jié)點，遠多于相鄰節(jié)點（26） 需要在表列中注明次相鄰點是否位于其它處理機 每層都要建立次相鄰表列，但次相鄰點可能不是物理意義的同一層,empty,M2L,相互作用表列 在迭代前存儲, 每次遠場作用 的轉(zhuǎn)移項都會 調(diào)用該表列,聚集,對于每一層的每個盒子, 聚集相當(dāng)于將子層組(t)中心的平面波移置到父層組(Pt )的中心(Ct ), 并通過內(nèi)插值得到大數(shù)目的平面波: 父節(jié)點(或部分子節(jié)點)在其它處理機上的節(jié)點構(gòu)成剩余八叉樹,它至多包含8pL個子節(jié)點, 因此數(shù)據(jù)交換的量級為O(logp+logL) 對壓縮八叉樹的所有節(jié)點(包括剩余節(jié)點)應(yīng)用并行聚合算法,因此每個處理機的計算