天體物理成團研究中的非規則訪存優化*

郝 赫，司雨蒙，韋建文，文敏華，林新華,2+

1.上海交通大學 高性能計算中心，上海 200240

2.NVIDIA Technology Center Asia Pacific，Singapore 999002

天體物理成團研究中的非規則訪存優化*

郝 赫1，司雨蒙1，韋建文1，文敏華1，林新華1,2+

1.上海交通大學 高性能計算中心，上海 200240

2.NVIDIA Technology Center Asia Pacific，Singapore 999002

HAO He,SI Yumeng,WEI Jianwen,et al.Optimizing irregular memory access in astrophysical clustering studies.Journal of Frontiers of Computer Science and Technology,2017,11(1)：80-90.

HGGF（halo-based galaxy group finder）算法實現了基于暗物質暈的星系找群，在研究宇宙大尺度結構及宇宙的演化等領域中占有至關重要的地位。但由于數據規模的增長，急需對HGGF算法進行優化，以縮短運行時間。經分析，算法的熱點部分耗時受到非規則訪存的嚴重影響，因此針對算法的結構和非規則訪存模型，提出了數據預排序方法，并分析了該方法如何影響訪存過程。在此基礎上，利用數據對齊、循環分解進一步優化訪存效率，利用負載均衡和互斥變量私有化的方法提高了OpenMP的并行效率，最終將HGGF應用使用12線程加速11.6倍，同時取得了更好的可擴展性。主要有三點貢獻：（1）分析了HGGF算法的非規則訪存問題；（2）提出并分析了數據預排序方法；（3）使用數據對齊、循環分解、負載均衡、互斥變量私有化方法提高了HGGF應用的并行性能。

天體物理成團；非規則訪存優化；數據預排序；并行計算

1 引言

隨著天文觀測能力的不斷提高，天文數據急劇增加，其中星系的觀測總量已經達到109量級。面對來自大型數字巡天計劃的海量數據，如何從數據中迅速準確地提取所需要的內容，直接影響著天文學的發展和研究進程。其中，如何將可觀測的星系和理論中不可見的暗物質暈聯系起來，是天體物理研究中的重要難題之一[1]。天文學家利用計算機設計出多種Group-Finding算法來尋找同一暗物質暈中的星系。FOF（friends-of-friends）是Group-Finding算法中最早的算法之一[2]，對于理解宇宙的結構和演化具有重大意義。HGGF（halo-based galaxy group finder）算法由FOF算法演變而來，雖然同樣是對星系進行分類，但是HGGF算法是依賴暗物質暈的性質，能實現對某些孤立的星系和暗物質暈的聯接[3]。HGGF算法還將亮度、紅移距離、暗暈質量等性質共同作為篩選的依據，相比僅依靠距離的算法更有實際意義和準確性。

星系成團類算法基本步驟是預先指定星系中心，根據連接長度[4]逐個將星系樣本分配到某一星系群，再反復迭代直到收斂，使用的數據結構多為基于指針的圖或樹，如FOF算法即基于K-D樹[5]實現。雖然整體框架基本相同，但是HGGF算法使用數組替代指針作為主要數據結構，而在天文學或分子動力學的N體模擬問題中，粒子間的位置和連接關系往往非常復雜，使用數組的存儲效率更高[6]。在HGGF算法的核心部分，由于連接體積和索引數組[7]的存在，主數組的訪問順序被徹底打亂，因此將HGGF算法的訪存方式視為非規則訪存。非規則訪存算法通常表現出極少的數據局部性和數據重用，大量的動態非連續存儲訪問以及細粒度的同步或通訊[8]，在現代處理器基于緩存的架構下，緩存的利用效率會對程序的性能起決定性作用。而非規則訪存算法由于訪問地址不連續，甚至是隨機訪問，緩存內的數據不能得到高效利用，加上CPU和內存之間的速度差異，即使CPU計算能力很強，也不能獲得好的性能。

隨著異構計算的發展，在科學計算或是工程計算中加速卡的使用已經非常普遍。人們嘗試將HGGF應用移植到Intel公司的MIC（many integrated core）上，但是由于處理器的性能相比CPU差距較大[9]，HGGF算法沒有在MIC上取得理想的性能，需要進一步優化。

數據預取是通過優化局部性改善訪存性能的常見方法。數據預取包括硬件預取和軟件預取，硬件預取一般由編譯器自動執行，而軟件預取需要在代碼中插入明顯的預取指令[10]。我們使用了兩種軟件預取方法，但效果一般。經過進一步分析算法結構，發現主數組訪存順序同時受到連接體積和索引數組影響，因此采用了按轉換坐標對輸入數據預排序的方法，明顯改善了數據局部性，減少了緩存未命中的次數，在單核上加速超過2倍。同時，利用OpenMP在單節點上實現了HGGF算法中篩選星系過程的多線程并行，并且從負載均衡、同步開銷等方面對算法進一步優化，最終在12線程上獲得11.6倍的加速，在將數據規模擴大16倍后，獲得15倍加速。

綜上，本文主要有以下三點貢獻：

（1）分析了HGGF算法的非規則訪存問題；

（2）提出并分析了數據預排序方法；

（3）使用傳統方法提高了HGGF應用的并行性能。

本文組織結構如下：第2章介紹相關工作；第3章介紹HGGF算法及其訪存模型；第4章詳細介紹優化方法；第5章是實驗環境及結果分析；最后總結全文，并提出下一步工作。

2 相關工作

目前，星系成團研究中對并行計算的使用還比較有限，現有工作局限在對FOF算法[11]、HOP算法[12]的研究中。文獻[11]介紹了卡內基梅隆大學的Gibson等人利用Hadoop集群對超大天文數據集分解處理的DiscFinder方法。DiscFinder基于初始版本的FOF算法，通過Hadoop和MapReduce實現數據并行。文獻[12]介紹了西北大學的Liu等人為HOP算法設計的并行方案。他們通過對K-D樹進行域分解，利用MPI（message passing interface）實現了多節點多核的并行加速。HGGF算法與FOF算法和HOP算法相比計算精度更高，與實際觀測數據更接近，在星系成團研究中的地位已得到廣泛認可[13-14]，但計算過程更復雜，數據依賴更多，尚未有文獻對其模擬運算性能進行研究。HGGF算法利用OpenMP共享內存實現的多任務并行效率較低，并且算法中數據互相依賴造成外層循環難以拆分，本文并未如文獻[12]采用多節點MPI并行，而是在深入分析算法結構后，針對訪存問題和共享變量等問題給出多種優化方法，在單節點多核環境下獲得很好的效果。

由于訪存延遲常常嚴重影響程序性能，至今已有多種方法通過改善數據局部性優化訪存性能，如緩存分塊[15]、單指令多數據（single instruction multiple data，SIMD）[16]、預取技術等。本文所述方法既包括如數據對齊、循環分解等在文獻[9,17]中記載的方法，也有針對算法特定結構提出的數據預排序方法。主數組-索引數組結構在N體模擬問題中應用廣泛。文獻[7]介紹的非規則訪存模型中包含主數組-索引數組結構，其中提出的數據重排方法基于Inspector-Executor機制[18]，通過在程序內插入預處理函數動態調整數據順序。這種方法雖然更靈活，但由于大規模數值模擬計算中循環迭代次數太多，節省的時間不足以抵消Inspector函數的耗時。

本文所述數據預排序方法在計算前的預處理階段根據索引數組特點重排輸入數據，與文獻[7]不同，對迭代次數多的大規模計算加速更明顯，同時還區別于文獻[19-20]中記載的直接排序方法，可應用于結構更復雜的訪存模式。

3 HGGF算法及其非規則訪存模型

3.1HGGF算法

本節將對算法的核心部分，同時也是程序熱點部分（占總運行時間超過90%）進行詳細介紹。

3.1.1 連接體積

文獻[4]首次提出了連接體積的概念。每個星系在篩選時，本應遍歷空間中所有其他星系，但若兩星系間距超過特定閾值，根據經驗它們不可能屬于同一類，沒有互相計算的必要。圖1為星系分布的二維空間抽象圖，每個點代表一個星系。圖中有效半徑外與圓心位置的星系肯定不同群，因此在后面的篩選過程中，只需計算球體（擴展到正方體）內的樣本。

初始化連接體積的偽代碼如下：

前面說的有效半徑不是實際的空間距離，而是將三維空間轉換成以L為邊長的正方體（簡稱L空間）后按比例換算出的數值。式（1）中，x、y、z分別代表利用輸入數據中星系的極坐標求得的直角坐標；x′、y′、z′為直角坐標經過放縮平移到L空間（0＜x′,y′,z′＜L）內的值，L′為經驗常數。式（2）中，x″、y″、z″為x′、y′、z′的取整，函數Func的值為Volume數組的下標（后面統稱為轉換坐標），并且可以保證僅當x″、y″、z″都相等時，Func的值才會相等。由于x″、y″、z″為原浮點數坐標的取整，誤差可能導致不同星系轉換坐標相同，此時，Linklist數組會以類似鏈表的方式逐一記錄下轉換坐標相等的星系編號，若未與其他星系相同，由Volume數組記錄星系編號。

Fig.1 Valid radius ofLspace圖1 L空間內的有效半徑

經過以上步驟，Volume和Linklist兩個數組記錄了所有的星系編號。進入篩選過程前，需要先確定有效半徑內包含了哪些星系。本文將以當前星系為中心（多次迭代后為星系群的中心），二倍有效半徑為邊長（如圖1虛線框）的正方體稱為連接體積，則連接體積必然包含所有有效半徑內的星系。建立索引數組時，按遞增順序篩選連接體積內的所有整數點，若該點的Func函數值剛好等于某一星系的轉換坐標，則將Volume中的星系編號加入索引數組，偽代碼如下：

掃描結束后，連接體積內所有星系的編號都存入了索引數組。

3.1.2 遍歷篩選星系

建立索引數組后進入篩選星系的過程，篩選的對象是連接體積內的星系，即索引數組和對應Linklist中存儲的星系編號。判斷星系是否同類的條件為式（3），Dist1為紅移距離，Dist2為空間距離，p1、p2、p3、p4和r、x、y、z分別代表被篩選星系和中心星系的坐標。每次篩選過程需要兩次判斷，若紅移距離不符合要求，則跳轉到下一個星系，若滿足條件，再判斷空間距離，如果兩個條件同時滿足，則記錄下該星系編號。此外，在每次篩選結束后會查看對應Linklist數組中是否存有其他星系。

3.2 非規則訪存模型

HGGF算法的非規則訪存模型為二重循環，外層循環為當前分類數，隨著星系成團的過程，星團數目會不斷減少，兩個內層循環分別對應上一小節建立索引數組和遍歷篩選星系的過程。第一個內層循環建立的索引數組用作第二個內層循環中主數組的下標，每個外層循環索引數組根據當前星系編號對應更新，包含新的要篩選的星系。訪存模型偽代碼如下：

因為編譯時主數組下標即索引數組的值不能立即獲取，所以訪存的順序完全由索引數組決定，且連接體積內的星系離散分布在所有樣本星系中。以最小輸入數據為例，樣本數為639 359個星系，則索引數組index的值不規則分布于區間1～639 359。主數組讀取數據時，會按一個緩存行的大小（64 B）讀入，即16個相鄰的數組元素，但受索引數組的影響，訪問地址不連續，需要的數組元素很難位于同一數據塊中，這會嚴重影響緩存性能。因此，可將這種訪存模式歸為非規則訪存。

4 優化方法

優化之前，使用VTune[21]測試出熱點函數每行代碼的具體耗時以及硬件性能計數器的數據。結果顯示篩選星系過程耗時占熱點函數的95%以上，因此優化工作主要針對這部分展開。

對于非規則訪存問題，軟件預取是優化局部性最常見的方法。由于沒有明顯性能提升，只簡要介紹軟件預取方法：通過插入具體預取指令，編譯器可提前將數據讀入緩存，并利用流水線上其他任務覆蓋這部分開銷。本文嘗試了兩種：

（1）使用#pragma prefetch指示編譯器進行預取；

（2）使用intrinsic代碼_mm_prefetch將數組元素預取到指定緩存。

從profile結果看，數據被成功預取，但開銷不僅沒有被覆蓋，還和優化前一樣。分析原因是需要預取的計算任務簡單，流水線不能覆蓋時間更長的訪存開銷[22]。

4.1 非規則訪存優化

4.1.1 數據預排序

通過分析連接體積算法的原理可以發現，循環內主數組元素的訪問順序由索引數組決定，而索引數組作為Volume數組的子集，存儲了當前星系有效半徑內的所有星系編號（除轉換坐標相等的星系由Linklist數組存儲）。因為遍歷連接體積時是按轉換坐標從小到大的順序，所以索引數組中星系編號的存儲位置同樣應該按轉換坐標從小到大排列。設想一下，如果能讓輸入數據在程序開始前按同樣順序排好，則通過遍歷連接體積內整數點得到的索引數組，其中的星系編號很可能已經是連續的，即1～639 359中的任意一個連續區間，但不包含轉換坐標相等的星系。

原輸入數據中星系的性質與轉換坐標并沒有直接關系，也沒有與轉換坐標順序一致或相似的一列數據，因此只有通過對輸入數據進行預處理，才能達到使輸入數據按轉換坐標有序的目的，這個預處理過程稱之為數據預排序。從優化的思路來看，數據預排序與Inspector-Executor有一定相似之處，兩種方法都是通過改變數據的分布來優化訪存性能，但前者是在程序運行前通過預處理輸入數據改變數據局部性，而Inspector-Executor是通過在程序中插入數據重排函數，動態改變數據的分布。

下面介紹數據預排序的具體操作方法。如圖2，首先利用與原程序相同的算法求出所有星系的轉換坐標，然后將得到的值按原輸入數據順序添加到文件中（原輸入數據中沒有轉換坐標那一列，相當于給每個星系增加一個新的性質），再按轉換坐標從小到大對輸入文件按行排序，并作為新的輸入文件，整個過程并不需要修改源程序的代碼。

Fig.2 Procedure of data pre-sorting圖2 數據預排序步驟

圖3是數據預排序前后索引數組分配的抽象圖，可看出預排序前，各星系編號被隨機分配到索引數組中，這樣讀取數據時，由于星系編號差別很大，導致每次讀入的一行數據中可能只有一個有效數據。預排序后，由于星系編號已經按轉換坐標排序，會按順序被分配到索引數組中，使得每次讀入緩存的數據塊中的全部數組元素都是本次循環或下個循環要被使用的數據。圖4是數據預排序前后索引數組（陰影部分）及Linklist數組（突出部分）存儲的星系編號，因為遍歷的順序類似于深度優先搜索，所以預排序后數組內星系編號可視為完全連續。雖然主數組-索引數組的結構未變，但是數據局部性已經徹底得到改善。

Fig.3 Distribution of index array圖3 索引數組分配

Fig.4 Change of index and Linklist array圖4 索引數組及Linklist數組變化

相比訪存延遲的巨大額外開銷，數據預排序本身的開銷可以被抵消掉，而且僅通過一次預排序，可以為之后多次模擬實驗節省時間，顯然是可以接受的。

4.1.2 數據對齊

篩選星系過程的主數組是二維數組，記錄了星系的全部16個性質，并且聲明在內存的相鄰位置，但是計算過程中只訪問到16個元素中的前4個（三維坐標、紅移距離），這會影響數據讀取的效率。如圖5（a）是數據未對齊時最好的情況，即一次將4個數組元素全部讀入，但是浪費了后面的12個單元。圖5（b）則是最差的情況，4個元素由于前端未對齊被分開，此時需要訪問兩次內存才能取出全部4個元素。圖5（c）是使用數據對齊，并且將數組分開聲明的情況，再配合數據預排序，一次可同時讀入4組有效數據，比之前最差的情況效率提高8倍。

Fig.5 Data alignment圖5 數據對齊

數據對齊方法：首先確保數組的起始地址對齊，使用可以將普通的內存分配和釋放改寫為對齊的intrinsic原語，即將malloc()改為_mm_malloc()，free()改為_mm_free()，由于剛好被16整除，不需要考慮后端對齊。

4.1.3 循環分解

由于CPU緩存空間有限，循環中無關數據過多可能使原本應保留在緩存中的數據被擠出，而換入的數據若再被換出，會形成惡性循環，導致有效數據總是不能保留在緩存中。如果將沒有依賴關系的計算拆開，分成兩個或多個循環，不僅可以提高各部分對緩存的利用效率，還可以分別展開進一步優化。

4.2 OpenMP并行效率優化

本文使用OpenMP對篩選星系過程實現并行，但在12線程上只加速不到4倍，因此使用了以下方法進一步優化算法。

4.2.1 負載均衡

為避免各線程任務分配不均，可利用OpenMP制導語句選項nowait或schedule實現任務的動態分配。加入選項后，編譯器會根據各線程的執行情況動態分配任務，先執行完的線程會被分配新的任務。兩種方法對HGGF算法的效果基本相同。

4.2.2 互斥變量私有化

因為程序中的部分操作和變量在迭代過程中只能串行執行或由多線程共享，所以需要使用鎖機制避免該部分程序因為多線程并行出現讀寫沖突。OpenMP中實現互斥最常見的方法是使用指導選項#ompcritical或#ompatomic。添加指導選項后，同一時間只能有一個線程對互斥區進行訪問。需要特別指出，atomic比critical雖然在性能上會有小幅提升，但由于atomic選項只能作用于簡單語句，當互斥區內指令比較復雜時，可能造成結果錯誤。

雖然指導選項使用便捷，但線程間只能互相等待，嚴重影響程序并行效率。HGGF算法中互斥區內指令的作用是計數及記錄星系編號，線程執行的先后順序不影響最終結果。利用這個特點，將線程序號與數組下標結合，為每個線程單獨聲明一組變量，使得原來需要互斥訪問的共享變量可以被各線程私有和同步讀寫，實現了線程間的完全并行。偽代碼如下：

4.3 MIC移植介紹

MIC采用傳統x86架構，擁有數十個精簡核心，提供高度并行的計算能力，支持現有程序的平滑移植。MIC有靈活的編程模式，既可以用作協處理器，也可以作為獨立的節點，Native和Offload是兩種比較常見的模式[23]。但是，MIC并不總是能提供好的性能，其CPU核心不支持亂序執行，復雜算法在MIC上的表現一般。與傳統多核CPU相比，MIC的SIMD和緩存對算法更敏感，訪存密集或訪存不規則的程序受到的影響更加嚴重。

為了進一步驗證文中優化方法的有效性，以及探索未來程序向眾核平臺移植的可能性，在MIC的Natvie模式下測試了HGGF應用，僅比較并行部分的耗時，即使開啟全部60個核上的60/120/180/240個線程，并行性能仍落后于CPU。

要在MIC上獲得好的性能，往往需要從多方面進行優化。受到算法限制，緩存分塊、單指令多數據等技術并不適用于HGGF算法。但考慮到同樣是x86架構，前面介紹的優化方法應該也適用于MIC。因此本文將同樣的優化方法應用在MIC的程序中，取得了與CPU上接近的加速比。

5 實驗環境及結果分析

5.1 實驗環境

5.1.1 數據集

原始HGGF算法已經在單線程機器上驗證并使用過。其使用的輸入數據是基于P3M算法得到的仿真數據[24]，輸出數據正確性可由2度視場星系紅移巡天計劃（2-degree Field Galaxy Redshift Survey，2dFGRS）和斯隆巡天計劃（Sloan Digital Sky Survey，SDSS）數據驗證。在本次并行化的升級算法中，將和原來的單線程結果相互比較驗證。

5.1.2 軟硬件環境

實驗的硬件環境如表1所示，軟件環境如表2所示。

Table 1 Hardware表1 硬件環境

Table 2 Software表2 軟件環境

5.2 結果分析

5.2.1 非規則訪存優化

為了驗證數據預排序對訪存性能的影響，利用VTune測出了程序單核運行時的各級緩存丟失率。從表3中可以看出，使用該方法后，各級緩存丟失率均有明顯下降。這說明數據預排序已經改善了局部性，使得CPU每次都能讀取更多的有效數據，與本文的分析基本一致。

Table 3 Cache miss rate during pre-sorting表3 預排序前后緩存丟失對比 %

確認了預排序方法的作用后，又在多線程上測試了程序的運行時間。如圖6，單獨使用預排序方法時（紅色線），單線程上的加速效果最明顯，有1.9倍加速，而使用全部12個線程后，僅有1.2倍加速。原因是使用多核并行后，每個核有獨立的一、二級緩存，等效于每次循環可將更多數據讀入緩存，削弱了數據預排序的效果。此外，數據對齊和緩存分解方法是已有記載的優化方法，不僅針對單一訪存模型，無需在此詳細分析。使用全部3種非規則訪存優化方法（藍色線）的運行時間，相比僅使用預排序方法有1.2倍加速，在不同核數上加速效果基本相同。因為核數增多后，數據預排序加速效果減弱，所以設計了使用更大數據集的實驗。如圖7，數據規模從1到16倍，預排序的加速比從1.2倍增加到1.6倍。原因是隨著數據規模擴大，緩存容量再次不足，數據預排序的效果逐漸回升，證明其弱可擴展性較好。

Fig.6 Strong scalability of memory optimization圖6 訪存優化強可擴展性

Fig.7 Weak scalability of data pre-sorting圖7 數據預排序弱可擴展性

Fig.8 Strong scalability of OpenMP optimization圖8 OpenMP并行優化強可擴展性

5.2.2 OpenMP并行效率優化

如圖8，分別測試了負載均衡（紅色線）和互斥變量私有化（藍色線）兩種方法的運行時間?？梢钥闯鲐撦d均衡對性能的提升并不明顯，只有5%左右，說明原算法線程任務分配已經比較均勻。而互斥變量私有化隨著使用核數的增加加速效果越來越明顯，12線程直接加速2倍，這也說明原算法中互斥區的存在嚴重影響了線程的并行。

最后來看一下綜合所有優化方法的加速效果及可擴展性。如表4，優化后12線程的加速比為11.6倍，相比優化前僅有3.9倍的加速有了大幅提高。而優化后整個程序（不僅Kernel部分）的加速比也達到8倍。圖9顯示了程序經整體優化后弱可擴展性良好。

5.2.3 MIC優化結果

本文分別測試了使用數據預排序方法和全部優化方法在MIC上的運行時間。如圖10，僅使用數據預排序時，最多加速1.4倍，但不同線程數的加速效果相差不大，主要是因為在MIC上每個核最多可開啟4個線程，雖然線程數不同，但是使用的核數都是60個，所以并未改變緩存的容量，就不會對數據預排序的效果產生影響。使用全部優化方法后，最多在240個線程上加速2.2倍，主要原因是線程越多，互斥變量私有化方法的效果越明顯。

Table 4 Speedup of all optimizations表4 所有優化方法加速情況

Fig.9 Weak scalability of final optimization圖9 最終優化弱可擴展性

Fig.10 Result on MIC圖10 MIC測試結果

6 總結與下一步工作

本文通過對天文應用Halo-based Galaxy Group Finder進行分析與優化，歸納出算法核心部分的非規則訪存模型，并對該類型問題使用數據預排序方法優化訪存性能，取得了不錯的優化結果。還分析了數據預排序方法對該非規則訪存模型產生影響的原因，為類似的訪存問題提供了新的優化思路和方法。在數據預排序的基礎上，使用了幾種傳統優化方法：（1）使用數據對齊、循環分解優化訪存性能；（2）使用負載均衡、互斥變量私有化提高算法的并行效率。綜合上述方法，在單節點12個線程上對HGGF應用加速11.6倍，并且隨著數據規模的擴大，加速效果更加明顯。將HGGF應用以本地模式移植到MIC卡上，并測試了性能，驗證了本文優化方法的有效性。綜上，本文的優化成果對于HGGF應用涉及的天文學問題的研究將有很大幫助，也為相似訪存類型應用提供了更多優化途徑。

在下一步工作中，將繼續通過優化算法改善其數據局部性和可擴展性，并利用多節點（MPI）或異構平臺（CPU+GPU，CPU+MIC）進行加速。

致謝感謝上海交通大學高性能計算中心提供的實驗環境。感謝上海交通大學物理與天文系楊小虎教授在程序優化過程中的悉心指導，以及對本文的修改意見。

[1]Yang Xiaohu,Mo H J,van den Bosch F C,et al.A halobased galaxy group finder:calibration and application to the 2dFGRS[J].Monthly Notices of the Royal Astronomical Society,2005,356(4):1293-1307.

[2]Huchra J P,Geller M J.Groups of galaxies:I-nearby groups [J].TheAstrophysical Journal,1982,257:423-437.

[3]Yang Xiaohu,Mo H J,van den Bosch F C,et al.Galaxy groups in the SDSS DR4.I.the catalog and basic properties [J].TheAstrophysical Journal,2007,671(1):153.

[4]Eke V R,Baugh C M,Cole S,et al.Galaxy groups in the2dFGRS:the group-finding algorithm and the 2PIGG catalogue[J].Monthly Notices of the Royal Astronomical Society,2004,348(3):866-878.

[5]Shevtov M,Soupikov A,Kapustin A.Highly parallel fast KD-tree construction for interactive ray tracing of dynamic scenes[J].Computer Graphics Forum,2007,26(3):395-404.

[6]Badawy A H,Aggarwal A,Yeung D,et al.The efficacy of software prefetching and locality optimizations on future memory systems[J].Journal of Instruction-Level Parallelism,2004,6(7).

[7]Lin Yuan,Padua D.Analysis of irregular single-indexed array accesses and its applications in compiler optimizations [C]//LNCS 1781:Proceedings of the 2000 Joint European Conferences on Theory and Practice of Software,Berlin, Mar 25-Apr 2,2000.Berlin,Heidelberg:Springer,2000: 202-218.

[8]Hennessy J L,Patterson D A.Computer architecture:a quantitative approach[M].San Francisco,USA:Morgan Kaufmann Publishers Inc,2011.

[9]Lin Xinhua,Li Shuo,Zhao Jiaming,et al.Node-level memory access optimization on Intel Knights Corner[J].Computer Science,2015,42(11):37-42.

[10]Gornish E H,Granston E D,Veidenbaum A V.Compilerdirected data prefetching in multiprocessors with memory hierarchies[C]//Proceedings of the 25th International Conference on Supercomputing Anniversary Volume,Tucson, USA,May 31-Jun 4,2011.New York:ACM,2014:128-142.

[11]Fu Bin,Ren Kai,Lopez J,et al.DiscFinder:a data-intensive scalable cluster finder for astrophysics[C]//Proceedings of the 19th International Conference on High Performance Distributed Computing,Chicago,USA,Jun 20-25,2010. New York:ACM,2010:348-351.

[12]Liu Ying,Liao Weikeng,Choudhary A.Design and evaluation of a parallel HOP clustering algorithm for cosmological simulation[C]//Proceedings of the 17th International Parallel and Distributed Processing Symposium,Nice,France,Apr 22-26,2003.Piscataway,USA:IEEE,2003:82.

[13]Faltenbacher A,Li Cheng,Mao Shude,et al.Three different types of galaxy alignment within dark matter halos[J]. TheAstrophysical Journal Letters,2007,662(2):L71-L74.

[14]Yang Xiaohu,Mo H J,van den Bosch F C.Galaxy groups in the SDSS DR4.II.halo occupation statistics[J].The Astrophysical Journal,2008,676(1):248-261.

[15]Nishtala R,Vuduc R W,Demmel J W,et al.When cache blocking of sparse matrix vector multiply works and why [J].Applicable Algebra in Engineering,Communication and Computing,2007,18(3):297-311.

[16]Knobe K,Lukas J D,Steele G L.Data optimization:allocation of arrays to reduce communication on SIMD machines [J].Journal of Parallel and Distributed Computing,1990,8 (2):102-118.

[17]Hager G,Wellein G.Introduction to high performance computing for scientists and engineers[M].Boca Raton,USA: CRC Press,Inc,2010.

[18]Arenaz M,Tourino J,Doallo R.An inspector-executor algorithm for irregular assignment parallelization[C]//LNCS 3358:Proceedings of the 2nd International Symposium on Parallel and Distributed Processing and Applications,Hong Kong,China,Dec 13-15,2004.Berlin,Heidelberg:Springer, 2005:4-15.

[19]Oh G H,Kim J M,Kang W H,et al.Reducing cache misses in hash join probing phase by pre-sorting strategy[C]//Proceedings of the 2012 SIGMOD International Conference on Management of Data,Scottsdale,USA,May 20-24,2012. New York:ACM,2012:864-864.

[20]Abouzied A,Abadi D J,Silberschatz A.Invisible loading: access-driven data transfer from raw files into database systems[C]//Proceedings of the 16th International Conference on Extending Database Technology,Genoa,Italy,Mar 18-22,2013.New York:ACM,2013:1-10.

[21]Wang P.How to use Vtune amplifier XE 2013 on Intel Xeon Phi coprocessor[R/OL].Intel Corporation(2012)[2015-10-19].https://software.intel.com/zh-cn/blogs/2012/11/05/vtuneamplifier-xe-2013-intel-xeon-phi-coprocessor?language=it. [22]Lee J,Kim H,Vuduc R W.When prefetching works,when it doesn’t,and why[J].ACM Transactions on Architecture and Code Optimization,2012,9(1):2.

[23]Jeffers J,Reinders J.Intel Xeon Phi coprocessor high performance programming[M].Beijing:Posts and Telecom Press,2014.

[24]Jing Y P,Suto Y.Triaxial modeling of halo density profiles with high-resolutionN-body simulations[J].The Astrophysical Journal,2002,574(2):538-553.

附中文參考文獻：

[9]林新華,李碩,趙嘉明,等.Intel Knights Corner的結點級內存訪問優化[J].計算機科學,2015,42(11):37-42.

[23]Jeffers J,Reinders J.Intel Xeon Phi協處理器高性能編程指南[M].北京:人民郵電出版社,2014.

HAO He was born in 1989.He is an M.S.candidate at Shanghai Jiao Tong University.His research interest is high performance computing.

郝赫（1989—），男，遼寧錦州人，上海交通大學碩士研究生，主要研究領域為高性能計算。

SI Yumeng was born in 1992.He is an M.S.candidate at Shanghai Jiao Tong University.His research interest is high performance computing.

司雨蒙（1992—），男，山西晉城人，上海交通大學碩士研究生，主要研究領域為高性能計算。

WEI Jianwen was born in 1986.He is an assistant engineer at the Center for High Performance Computing of Shanghai Jiao Tong University.His research interest is high performance computing.

韋建文（1986—），男，廣西百色人，上海交通大學高性能計算中心助理工程師，主要研究領域為高性能計算。

WEN Minhua was born in 1988.He is an assistant engineer at the Center for High Performance Computing of Shanghai Jiao Tong University.His research interest is high performance computing.

文敏華（1988—），男，江西會昌人，上海交通大學高性能計算中心助理工程師，主要研究領域為高性能計算。

LIN Xinhua was born in 1979.He is the vice director at the Center for High Performance Computing of Shanghai Jiao Tong University.His research interests include computer architecture and code optimization.

林新華（1979—），男，浙江紹興人，上海交通大學高性能計算中心副主任，主要研究領域為體系結構，代碼優化。

Optimizing Irregular MemoryAccess inAstrophysical Clustering Studies*

HAO He1,SI Yumeng1,WEI Jianwen1,WEN Minhua1,LIN Xinhua1,2+
1.Center for High Performance Computing,Shanghai Jiao Tong University,Shanghai 200240,China
2.NVIDIATechnology CenterAsia Pacific,Singapore 999002
+Corresponding author:E-mail:james@sjtu.edu.cn

Halo-based galaxy group finder(HGGF)tries to find galaxies in the same dark matter halo which is not directly visible.It plays a very important role in the research of large-scale structure of the universe.However,because of the growth of data scale,it’s extremely necessary to increase the running speed by optimizing the group finder coding algorithm.After a thorough investigation on the original HGGF code,it is found that the kernel part of the algorithm is seriously affected by the irregular memory access.This paper proposes a specific data pre-sorting approach and analyzes how it affects the process of memory access according to the structure of the algorithm and the irregular memory access pattern.Moreover,this paper uses data alignment and loop fission to optimize the memory access as well as improving the efficiency of OpenMP with load balance and mutex privatization.Eventually the HGGF application gets 11.6 times speedup on 12 threads,and gets better weak scalability.The following is the original contributions:(1)Analyze the irregular memory access of the HGGF application;(2)Propose and analyze the data pre-sorting;(3)Improve the parallel performance of HGGF application with another four approaches including data alignment,loop fission,load balance and mutex privatization.

A

：TP391

10.3778/j.issn.1673-9418.1512078

*The National High Technology Research and Development Program of China under Grant No.2014AA01A302(國家高技術研究發展計劃(863計劃));the Project of Japan Society for the Promotion of Science RONPAKU Fellowship(日本學術振興會項目).

Received 2015-11,Accepted 2016-01.

CNKI網絡優先出版:2016-01-07,http://www.cnki.net/kcms/detail/11.5602.TP.20160107.1540.008.html

Key words:astrophysical clustering studies;optimizing irregular memory access;data pre-sorting;parallel computing