田英齊,畢玉江,賀雨晴,馬運(yùn)恒2,,劉朝峰2,,徐 順

1(中國(guó)科學(xué)院 計(jì)算機(jī)網(wǎng)絡(luò)信息中心,北京 100190)

2(中國(guó)科學(xué)院大學(xué),北京 100049)

3(中國(guó)科學(xué)院 高能物理研究所,北京 100049)

格點(diǎn)量子色動(dòng)力學(xué)(Lattice Quantum Chromo Dynamics,LQCD) 是描述夸克和膠子之間的強(qiáng)相互作用的理論,它將連續(xù)時(shí)空離散成有限大小的四維超立方晶格,將QCD 理論中的夸克場(chǎng)放在超立方晶格的格點(diǎn)上,晶格點(diǎn)之間的連線上放置膠子場(chǎng).可測(cè)量的物理量是這兩種場(chǎng)的函數(shù),夸克場(chǎng)和膠子場(chǎng)的分布情況反映了QCD 的動(dòng)力學(xué)性質(zhì),所以只要在格子上給定了這兩種場(chǎng)的分布,原則上就能計(jì)算所有物理量(的期望值),因此格點(diǎn)QCD 模擬[1]是一個(gè)重要的基礎(chǔ)方法,具體來說是目前已知能系統(tǒng)研究夸克及膠子間低能強(qiáng)相互作用的不可缺少的非微擾計(jì)算方法.

格點(diǎn)QCD 的一個(gè)重要研究課題是膠球性質(zhì)的研究.膠球是傳統(tǒng)夸克模型中沒有的強(qiáng)子態(tài),由能發(fā)生自相互作用的膠子組成.我國(guó)的BESIII 合作組一直在尋找膠球,并發(fā)現(xiàn)了一些可能的候選者.但要確認(rèn)它們,需要將理論計(jì)算與實(shí)驗(yàn)測(cè)量結(jié)合起來,因此膠球性質(zhì)的理論研究是一個(gè)重要的前沿課題,膠球的理論計(jì)算工作有助于基礎(chǔ)物理科學(xué)的新發(fā)現(xiàn).

格點(diǎn)QCD 是一個(gè)典型的計(jì)算密集型同時(shí)也是通信密集型的計(jì)算工作.常用馬可夫鏈實(shí)現(xiàn)蒙特卡洛重點(diǎn)抽樣,并且采用各種高級(jí)算法如pseudo-heat-bath 算法來加速組態(tài)生成過程的收斂；高精度大規(guī)模組態(tài)抽樣的文件輸出需要高效的計(jì)算IO；格點(diǎn)化的場(chǎng)近程作用是四維stencil 計(jì)算的一種具體形式,其計(jì)算效率不僅受處理器頻率影響也受內(nèi)存緩存大小影響；另外格距越小,計(jì)算物理觀測(cè)結(jié)果越精確,但帶來的計(jì)算量指數(shù)增加,造成了很多物理研究中采用近似算法.總之,格點(diǎn)QCD 對(duì)并行計(jì)算優(yōu)化提出了很高的需求.

針對(duì)格點(diǎn)QCD 的并行計(jì)算優(yōu)化受到高度重視,超級(jí)計(jì)算機(jī)上格點(diǎn)QCD 計(jì)算優(yōu)化及作業(yè)調(diào)度優(yōu)化多次獲得戈登貝爾(Gordon Bell) 獎(jiǎng)(分別是1987年、1995年、1998年和2006年,以及2018年入圍決賽).2005年前后,IBM、RBC/UKQCD和RIKEN 等設(shè)計(jì)了專門用于格點(diǎn)QCD 計(jì)算的超級(jí)計(jì)算機(jī)QCDOC(QCD on Chips),取代了其更早的QCDSP (QCD on DSPs)機(jī)器,IBM 后來還打造了更強(qiáng)大的Blue Gene/L平臺(tái),其Pavlos Vranas 團(tuán)隊(duì)基于該平臺(tái)實(shí)現(xiàn)了3.2 萬核12 Tflops 的格點(diǎn)QCD 應(yīng)用高性能計(jì)算,從而獲得了2006年的戈登貝爾獎(jiǎng)；2018年P(guān)avlos Vranas 再次采用GPU 加速格點(diǎn)QCD 計(jì)算,在美國(guó)Sierra和Summit超算平臺(tái)上以高效的計(jì)算性能晉級(jí)了當(dāng)年的戈登貝爾決賽[2].格點(diǎn)QCD 已逐漸成為美國(guó)、日本及歐洲超算平臺(tái)上的重要應(yīng)用[3].

當(dāng)前以美國(guó)USQCD 組織研發(fā)的格點(diǎn)QCD 軟件套件為領(lǐng)域內(nèi)主流軟件工具,其為四層軟件框架,最上層應(yīng)用層包括Chroma、CPS和MILC 等,下面三層為支撐庫(kù)層,分別實(shí)現(xiàn)消息通訊、數(shù)據(jù)結(jié)構(gòu)定義和基本算法實(shí)現(xiàn).MILC 軟件工具也成為SPEC CPU2006 基準(zhǔn)測(cè)試程序[4].QUDA 是基于CUDA 的格點(diǎn)QCD 計(jì)算庫(kù)[5],實(shí)現(xiàn)了較底層的計(jì)算加速,可以和上層的Chroma和MILC 等集成.Mikhail Smelyanskiy 等人基于多核處理器的Cache 架構(gòu)進(jìn)行了格點(diǎn)QCD 的線程化MPI 優(yōu)化,相關(guān)工作發(fā)表在2012年Supercomputing Conference會(huì)議上[6],Issaku Kanamori 基于Intel KNL和Oakforest-PACS 超算系統(tǒng)做了格點(diǎn)QCD 設(shè)計(jì)實(shí)現(xiàn)的討論,涉及SIMD和MPI 優(yōu)化方面[7].日本高能加速器研究機(jī)構(gòu)(KEK) 也于2009年開始了用于格點(diǎn)QCD模擬的Bridge ++軟件工具研發(fā)[8].我國(guó)存在一個(gè)“中國(guó)格點(diǎn)QCD”的合作組,當(dāng)前正在努力推出國(guó)際上較有影響的格點(diǎn)QCD 計(jì)算軟件.

本文研究工作探討格點(diǎn)QCD 計(jì)算的瓶頸問題,并提出有效的并行優(yōu)化計(jì)算方案.相關(guān)程序原始代碼由中科院高能物理研究所研究員、“中國(guó)格點(diǎn)QCD”成員陳瑩提供.為了實(shí)現(xiàn)其大規(guī)模并行計(jì)算、本文進(jìn)行多方面的并行計(jì)算分析和性能優(yōu)化研究.

1 格點(diǎn)QCD 模擬的基本流程

本節(jié)給出在淬火近似下格點(diǎn)QCD 生成規(guī)范場(chǎng)組態(tài)的基本流程.一個(gè)具體的膠子場(chǎng)分布叫做一個(gè)組態(tài),產(chǎn)生上述的組態(tài)需要涉及名為馬可夫鏈的隨機(jī)行走過程.整個(gè)鏈條通常從一個(gè)最簡(jiǎn)單的或隨機(jī)的膠子場(chǎng)分布開始起步.每一次步進(jìn)更新就意味著由當(dāng)前的場(chǎng)分布更新到另一種物理上允許的分布.

組態(tài)生成與測(cè)量流程可以分為如下幾步：

1) 從一個(gè)給定的場(chǎng)分布開始,計(jì)算每個(gè)link 的周遭環(huán)境,使用pseudo-heat-bath 算法計(jì)算link 更新的幾率,據(jù)此決定是否將其更新,對(duì)更新造成的link 形變加以修正；

2) 重復(fù)以上過程將整個(gè)格子上的所有l(wèi)ink 更新若干次,直到達(dá)到熱平衡,熱平衡之后繼續(xù)不斷更新格子上的link,使之成為新的組態(tài),將新的組態(tài)導(dǎo)出保存；

3) 在生成好的組態(tài)上測(cè)量膠球關(guān)聯(lián)函數(shù),繼續(xù)生成組態(tài),保存和測(cè)量,直到測(cè)量了足夠多的組態(tài)樣本.

更詳細(xì)的程序流程見圖1.

圖1 膠子場(chǎng)組態(tài)生成與測(cè)量的流程圖

在我們的程序中使用的pseudo-heat-bath 算法用于加快收斂.當(dāng)我們利用馬可夫鏈抽樣,每個(gè)樣本點(diǎn)對(duì)應(yīng)的膠子場(chǎng)分布都是我們所要的組態(tài).在有夸克場(chǎng)時(shí),膠球與通常介子態(tài)之間的混合會(huì)給研究增加復(fù)雜度,需要更大的計(jì)算量.因此作為第一步,忽略夸克場(chǎng)的所謂淬火近似下的模擬計(jì)算具有重要意義.對(duì)淬火組態(tài)產(chǎn)生程序的優(yōu)化,目的在于給出更多的組態(tài),從而得到更大的統(tǒng)計(jì)量和更好的信號(hào).淬火近似下膠球譜的研究工作和J/psi 粒子輻射衰變到標(biāo)量膠球或張量膠球的淬火近似研究工作[9-11],其結(jié)果都發(fā)表在物理學(xué)一流學(xué)術(shù)期刊上.這些計(jì)算中用到組態(tài)產(chǎn)生程序正是本研究?jī)?yōu)化的對(duì)象.

程序?qū)崿F(xiàn)中主要涉及的功能模塊有：

1) 數(shù)據(jù)初始化模塊：體系的格點(diǎn)劃分,內(nèi)存數(shù)組分配等；

2) 組態(tài)生成模塊：利用pseudo-heat-bath 算法生成滿足分布的組態(tài)；

3) 組態(tài)計(jì)算模塊：基于抽樣的組態(tài)計(jì)算膠球的物理性質(zhì)；

4) 輔助分析模塊：組態(tài)驗(yàn)證、統(tǒng)計(jì)計(jì)時(shí)、并行寫組態(tài)等.

2 并行計(jì)算關(guān)鍵問題

QCD 理論是一個(gè)基于SU(3)規(guī)范對(duì)稱性的理論.在格子上,每一條晶格連線上的膠子場(chǎng)的值由一個(gè)3 乘3 的SU(3)復(fù)矩陣描寫：Uμ(x).組態(tài)生成程序從一個(gè)給定的場(chǎng)分布開始,計(jì)算每條連線(也叫l(wèi)ink 變量)的周遭環(huán)境,這里涉及大量矩陣計(jì)算,并且每個(gè)格點(diǎn)需要做獨(dú)立性分析以適合并行化計(jì)算處理.

我們的程序按照pseudo-heat-bath 算法來更新格子上的每個(gè)link 變量的值,按以下幾率分布將一個(gè)link 變量的值從Uμ(x)更新為Uμ(x)’：

其中,A刻畫該link 周圍的環(huán)境,由其他link 變量對(duì)應(yīng)的矩陣值的乘積給出.計(jì)算A的過程中會(huì)涉及到相應(yīng)link (link 構(gòu)成訂書針的形狀,被稱為“staple”)相乘的計(jì)算,其示意圖如圖2.

圖2 link 周圍環(huán)境(staple)計(jì)算示意圖

要計(jì)算圖中的link 變量對(duì)應(yīng)的A,需要對(duì)其近鄰以及次近鄰的格點(diǎn)上相關(guān)的link 變量的乘積(即3×3 復(fù)矩陣的矩陣乘積)進(jìn)行求和,即對(duì)上圖中的8 個(gè)“staple”進(jìn)行求和：

其中,Ai為每個(gè)staple 內(nèi)涉及的link 的乘積,以圖2中的Ⅰ為例,其計(jì)算公式如下：

由此可以看出,這里涉及大量3×3 復(fù)矩陣的矩陣乘積計(jì)算優(yōu)化問題.

由于要更新一個(gè)link 變量,只與其近鄰和次近鄰的相關(guān)link 有關(guān),而與其他link 無關(guān).因此這些相互無關(guān)的link 之間可以并行處理.傳統(tǒng)的格點(diǎn)QCD 程序一般采用blocking和even-odd 算法進(jìn)行并行.該算法將相應(yīng)的link 分配到不同的block 上,根據(jù)block 的x、y、z、t的坐標(biāo)之和的奇偶分為兩組,每組內(nèi)可以進(jìn)行并行計(jì)算.然而我們的膠球模擬中涉及到了近鄰以及次近鄰的相互作用,因此需要首先將24 個(gè)block 組合為一個(gè)hyper-block,之后對(duì)hyper-block 進(jìn)行evenodd 的并行計(jì)算,來生成相應(yīng)的組態(tài).一個(gè)二維結(jié)構(gòu)格子的Hyper-blocking 與even-odd 算法示意如圖3.

圖3 even-odd 算法示意如圖

圖3中同色的虛線框之間可以進(jìn)行并行計(jì)算.

另外格點(diǎn)QCD 需要重復(fù)以上將整個(gè)格子上的所有l(wèi)ink 的值更新若干次,達(dá)到熱平衡,熱平衡之后繼續(xù)不斷更新格子上的link 的值,使之成為新的組態(tài),而組態(tài)的存儲(chǔ)空間占用也比較大,因此將新的組態(tài)導(dǎo)出保存也是并行化時(shí)需要考慮的問題.

3 并行分析與優(yōu)化策略

3.1 瓶頸與熱點(diǎn)分析

為了實(shí)現(xiàn)進(jìn)行大規(guī)模的格點(diǎn)QCD 計(jì)算,我們首先需要分析程序運(yùn)行可能的熱點(diǎn).原型程序采用消息傳遞接口(Message Passing Interface,MPI)來實(shí)現(xiàn)多節(jié)點(diǎn)并行處理,采用基于格點(diǎn)區(qū)域計(jì)算任務(wù)劃分的方式來實(shí)現(xiàn)并行處理,需要在各個(gè)MPI 進(jìn)程之間進(jìn)行格點(diǎn)邊界通信.通信占比與格子的劃分方式有關(guān),具體地說與子格子的面積體積比有關(guān).更小的面積體積比有更小的通信占比.我們知道同樣體積的情況下,正方體比長(zhǎng)方體有更小的面積體積比,因此導(dǎo)致了在通常情況下,x、y、z方向的格子規(guī)模相同,且劃分方式相同,t方向格子規(guī)模與劃分方式可以與x、y、z方向不同.因而在格子規(guī)模V和計(jì)算節(jié)點(diǎn)規(guī)模N確定的情況下,確定劃分方式即為求以下方程的正整數(shù)解：

其中,lt是t方向的子格子大小,lxyz是x、y、z方向的子格子大小.該方程的正整數(shù)解一般是確定的,因而不能靈活的協(xié)調(diào)通信與計(jì)算的比例,來獲取更優(yōu)的性能.為了更加靈活的分配計(jì)算任務(wù),獲取通信與計(jì)算的平衡,我們可以加入多線程(如OpenMP)的并行處理,更加充分的利用計(jì)算資源.

基于原型程序的初步測(cè)試,我們發(fā)現(xiàn)在常用的計(jì)算規(guī)模下,計(jì)算耗時(shí)占比60%以上.此程序的計(jì)算中,涉及大量3×3 復(fù)矩陣的矩陣乘積計(jì)算,是此程序的計(jì)算熱點(diǎn).

另外,運(yùn)行時(shí)輸出組態(tài)文件無疑也是一個(gè)性能瓶頸,尤其是在格子規(guī)模較大的情況下.對(duì)于964 大小的格子,組態(tài)文件大小為40 GB 以上.然而原型程序采用了串行輸出組態(tài)文件的方法,耗時(shí)極大,因此也需要采取優(yōu)化措施.

因此通過分析與測(cè)試,為了實(shí)現(xiàn)程序有效并行處理,我們考慮以下方面：

1)平衡通信與計(jì)算占比

2)復(fù)矩陣矩陣乘等數(shù)值計(jì)算函數(shù)

3)輸出組態(tài)文件

以下分別介紹這幾方面的計(jì)算瓶頸或計(jì)算熱點(diǎn)的優(yōu)化設(shè)計(jì).

3.2 MPI+OpenMP 并行優(yōu)化

針對(duì)膠球模擬的特點(diǎn),在xyzt四個(gè)維度上均可以進(jìn)行并行.程序可以設(shè)置各個(gè)方向分別的劃分方式.假設(shè)格子的體積設(shè)置為(Nt,Nx,Ny,Nz),在t、x、y、z方向分別劃分為mt、mx、my、mz份,則每個(gè)MPI 進(jìn)程上的格子體積為(Nt/mt,Nx/mx,Ny/my,Nz/mz),劃分方式記為(mt,mx,my,mz).以(2,2,2,2)劃分方式為例,4 個(gè)維度的并行分割方式如圖4所示.

圖4 MPI 并行劃分方式(2,2,2,2)示意

基于前述的分析,對(duì)于相同規(guī)模的格子,使用MPI+OpenMP 的并行方式,可以有效增加單個(gè)MPI 進(jìn)程上格子的規(guī)模,從而降低子格子的面積體積比,從而減少通信的占比.因此我們采用了動(dòng)態(tài)負(fù)載平衡的方法,在MPI 進(jìn)程內(nèi)部進(jìn)行OpenMP 的并行.圖5以二維結(jié)構(gòu)示意單個(gè)MPI 進(jìn)程內(nèi)部2 個(gè)OpenMP 線程情況下的OpenMP 并行情況：

圖5 OpenMP 并行示意圖

3.3 矩陣操作的向量化優(yōu)化

在程序中,涉及到的計(jì)算主要是3×3 的復(fù)矩陣的矩陣乘計(jì)算,即兩個(gè)3×3 的復(fù)矩陣相乘,以及復(fù)矩陣與矩陣的厄密共軛矩陣相乘.具體計(jì)算如下所示：

我們針對(duì)復(fù)矩陣矩陣乘操作,設(shè)計(jì)了多種不同的計(jì)算方法,以下將分別介紹：

1) 我們基于結(jié)構(gòu)體的數(shù)組轉(zhuǎn)化為數(shù)組的結(jié)構(gòu)體(AoS to SoA)的思想,設(shè)計(jì)了數(shù)組長(zhǎng)度可變的結(jié)構(gòu)體數(shù)組,用來存儲(chǔ)中間數(shù)據(jù),提高了計(jì)算效率.經(jīng)過測(cè)試,結(jié)構(gòu)體內(nèi)數(shù)組長(zhǎng)度為8 時(shí)的效果最好.

2)我們將復(fù)矩陣矩陣乘的數(shù)學(xué)算式完全展開,利用編譯器的自動(dòng)向量化功能進(jìn)行相關(guān)代碼的向量化.

3) 我們針對(duì)AVX-512 指令集架構(gòu),利用Intel Intrinsics 設(shè)計(jì)了手動(dòng)向量化的算法.對(duì)于這樣的矩陣計(jì)算,我們可以將v矩陣的行向量實(shí)部和虛部[c1,d1,c2,d2,c3,d3]組裝為一個(gè)向量,針對(duì)實(shí)部和虛部采用_mm512_fmaddsub_pd()函數(shù)進(jìn)行計(jì)算.

3.4 輸出組態(tài)文件的優(yōu)化

原始的組態(tài)文件按照一般多維數(shù)組的存儲(chǔ)方式,以行主元的方式進(jìn)行存儲(chǔ),即依次按照z、y、x、t下標(biāo)由小到大進(jìn)行存儲(chǔ).原程序依據(jù)該存儲(chǔ)方法由0 號(hào)進(jìn)程進(jìn)行數(shù)據(jù)收集及輸出到文件的工作.以42大小的格子,(2,2)劃分方式為例,4 個(gè)MPI 進(jìn)程的組態(tài)文件輸出方式如圖6所示.

圖6 MPI 并行寫入優(yōu)化前算法示意圖

對(duì)于N個(gè)進(jìn)程的情況,原始的輸出方式需要進(jìn)行的通訊量為(N-1)×Vsub,其中Vsub為子格子的大小.而且需要進(jìn)行N×mt×mx×my次不連續(xù)的寫操作.

我們修改了組態(tài)文件的數(shù)據(jù)排布方式,將整個(gè)文件劃分為N 部分,并在每個(gè)部分內(nèi)按照數(shù)據(jù)連續(xù)的方式進(jìn)行排布,依次實(shí)現(xiàn)了多個(gè)MPI 進(jìn)程并行輸出組態(tài)文件的算法.采用此種方法后,理論上沒有通信操作,也沒有不連續(xù)的內(nèi)存寫操作,可以極大的提高輸出組態(tài)文件的性能.并行寫采用了MPI-IO,MPI 中可以使用 mpi_file_open 讓所有進(jìn)程都打開同一個(gè)文件,使用mpi_file_write_at 函數(shù)獨(dú)立互不干擾地寫每個(gè)進(jìn)程自己負(fù)責(zé)的數(shù)據(jù).圖7展示了優(yōu)化后的數(shù)據(jù)寫入方式：

圖7 MPI 并行寫入優(yōu)化后算法示意圖

4 實(shí)驗(yàn)結(jié)果與分析

4.1 實(shí)驗(yàn)準(zhǔn)備

我們?cè)凇疤旌? 號(hào)”超級(jí)計(jì)算平臺(tái)以及KNL 平臺(tái)上測(cè)試了我們的代碼,測(cè)試平臺(tái)的信息如表1所示.

表1 軟硬件環(huán)境清單

其中“天河2 號(hào)”超級(jí)計(jì)算機(jī)不支持AVX 指令向量化,而Intel KNL 支持AVX512 指令向量化.

4.2 OpenMP 多線程加速測(cè)試

基于前述3.1 節(jié)的整體分析與3.2 節(jié)提出的OpenMP多線程加速方案,我們?cè)O(shè)定了兩種格子規(guī)模和計(jì)算規(guī)模的情況,考察MPI 單節(jié)點(diǎn)上采用OpenMP 多線程調(diào)整計(jì)算與通信占比以優(yōu)化并行效率的效果.同時(shí)我們也考察了優(yōu)化前后程序強(qiáng)可擴(kuò)展性和弱可擴(kuò)展性的對(duì)比情況.

4.2.1 常規(guī)規(guī)模測(cè)試結(jié)果

我們首先針對(duì)常用的格子與計(jì)算規(guī)模,進(jìn)行了程序優(yōu)化的測(cè)試.我們?cè)O(shè)定格子大小為96×483,使用了5184 個(gè)計(jì)算核心,測(cè)試結(jié)果如表2.

表2 常規(guī)規(guī)模測(cè)試結(jié)果

從測(cè)試結(jié)果可以看出,在規(guī)模不變的情況下,原始的MPI 并行并不是最優(yōu)結(jié)果.隨著OpenMP 線程數(shù)的增加,性能逐漸提高,達(dá)到了理想的計(jì)算效率.較小規(guī)模的情況下最高達(dá)到了3.32 倍的加速.

4.2.2 大規(guī)模測(cè)試結(jié)果

我們又選擇了當(dāng)前環(huán)境下可使用的最大計(jì)算規(guī)模進(jìn)行了測(cè)試.我們?cè)O(shè)定的格子大小為964,使用了41 472個(gè)計(jì)算核心,測(cè)試結(jié)果如表3.

表3 大規(guī)模測(cè)試結(jié)果

從測(cè)試結(jié)果可以看出,在規(guī)模不變的情況下,原始的MPI 并行同樣不是最優(yōu)結(jié)果.加入OpenMP 并行后,可以使得格子的劃分更加靈活,可以找到合理的劃分方式,達(dá)到理想的計(jì)算效率.較大規(guī)模的情況下最高達(dá)到了1.37 倍的加速.

4.2.3 劃分方式對(duì)計(jì)算效率的影響分析

根據(jù)前述3.1 與3.2 節(jié)的分析,我們可以知道減少通信占比的方法是增大子格子的面積體積比.在常規(guī)規(guī)模的測(cè)試中,我們可以看出,通過加入OpenMP 并行,可以有效減小子格子的面積體積比,取得了很好的加速效果.對(duì)于大規(guī)模的測(cè)試結(jié)果,我們發(fā)現(xiàn)單純?cè)龃笞痈褡用娣e體積比并沒有導(dǎo)致計(jì)算性能的持續(xù)提高,針對(duì)這一情況,我們分析原因與計(jì)算機(jī)的網(wǎng)絡(luò)結(jié)構(gòu)有關(guān).計(jì)算任務(wù)在節(jié)點(diǎn)上的分配方式需要進(jìn)一步的優(yōu)化研究.

4.2.4 強(qiáng)可擴(kuò)展性測(cè)試

針對(duì)MPI+OpenMP 并行方案,我們測(cè)試了程序的強(qiáng)擴(kuò)展性.我們對(duì)比了只有MPI 并行的情況和MPI 與OpenMP 混合并行的情況.強(qiáng)擴(kuò)展性測(cè)試用例的格子大小為644,最高并行達(dá)到了4096 核.測(cè)試結(jié)果如圖8所示.

圖8 強(qiáng)擴(kuò)展性測(cè)試結(jié)果

結(jié)果顯示在應(yīng)用了OpenMP 多線程加速后,強(qiáng)擴(kuò)展性有了明顯的提高.

4.2.5 弱可擴(kuò)展性測(cè)試

我們針對(duì)MPI+OpenMP 并行方案,測(cè)試了程序的弱擴(kuò)展性.弱可擴(kuò)展性測(cè)試的規(guī)模為：每個(gè)核上分配有84大小的格子,依次測(cè)試了256 核～6144 核的計(jì)算結(jié)果.測(cè)試結(jié)果如圖9.

圖9顯示了各次計(jì)算中單步生成組態(tài)的計(jì)算時(shí)間,在理想情況下,單步時(shí)間應(yīng)該是保持不變的,然而隨著體系的增大,實(shí)際會(huì)有所增加.我們用最小二乘法線性擬合數(shù)據(jù)結(jié)果,得到了兩條擬合曲線的斜率,以此來評(píng)價(jià)程序的弱擴(kuò)展性.結(jié)果顯示在應(yīng)用了OpenMP 多線程加速后,程序的弱可擴(kuò)展性也得到了提高.

圖9 弱擴(kuò)展性測(cè)試結(jié)果

4.3 向量化優(yōu)化測(cè)試

針對(duì)3.3 節(jié)提出的向量化優(yōu)化方案,我們?cè)贙NL平臺(tái)上測(cè)試了AVX-512 向量化優(yōu)化的計(jì)算結(jié)果.我們通過-xMIC_AVX512 編譯選項(xiàng)控制向量化開關(guān),進(jìn)行了向量化與不向量化的測(cè)試.測(cè)試中使用了324大小的格子,運(yùn)行在MPI+OpenMP 的模式下,MPI 綁定CPU采用I_MPI_PIN_DOMAIN=omp：compact,OpenMP設(shè)置親緣性KMP_AFFINITY=compact.測(cè)試基準(zhǔn)為未向量化的Fortran 內(nèi)置的matmul 函數(shù).后續(xù)測(cè)試了采用-xMIC_AVX512 編譯選項(xiàng)進(jìn)行向量化的matmul函數(shù)；采用AoS to SoA 優(yōu)化與-xMIC_AVX512 編譯選項(xiàng)結(jié)合的算法；采用公式展開和編譯選項(xiàng)結(jié)合的算法；以及采用了Intel Intrinsics 手動(dòng)向量化的算法.幾種算法在16MPI + 1OpenMP 的情況下結(jié)果如圖10.

圖10 16MPI + 1OpenMP 各向量化優(yōu)化測(cè)試結(jié)果及加速比

運(yùn)行在16MPI+16OpenMP 的情況下,測(cè)試結(jié)果如圖11.

圖11 16MPI + 16OpenMP 各向量化優(yōu)化測(cè)試結(jié)果及加速比

結(jié)果顯示,幾種向量化算法均優(yōu)于Fortran 內(nèi)置的matmul 函數(shù)的效果,而且都取得了3 倍左右的加速效果.其中基于Intel Intrinsics 的手動(dòng)向量化算法效果最好,達(dá)到了3.17 倍的加速.

4.4 輸出組態(tài)文件優(yōu)化測(cè)試

針對(duì)3.4 節(jié)提出的組態(tài)輸出優(yōu)化方案,我們測(cè)試了未改變數(shù)據(jù)順序的MPI 并行輸出以及修改數(shù)據(jù)順序后的MPI 并行輸出算法,以單個(gè)MPI 進(jìn)程輸出組態(tài)文件的結(jié)果為基準(zhǔn),具體實(shí)驗(yàn)數(shù)據(jù)如下：

表4 MPI 并行寫組態(tài)文件測(cè)試結(jié)果

表4中“串行寫”的數(shù)據(jù)為測(cè)試的基準(zhǔn),并行寫”為修改了數(shù)據(jù)順序后的并行寫算法測(cè)試結(jié)果,“加速比”為其對(duì)應(yīng)的加速比.通過表4我們可以看出,優(yōu)化了數(shù)據(jù)順序的并行寫算法取得了上百倍的加速,最高達(dá)到了153 倍.

5 總結(jié)與未來工作

格點(diǎn)量子色動(dòng)力學(xué)一直是高性能計(jì)算領(lǐng)域的熱門應(yīng)用方向,多年來一直伴隨著超級(jí)計(jì)算機(jī)的發(fā)展,對(duì)高性能計(jì)算的發(fā)展有強(qiáng)有力的驅(qū)動(dòng)力.本文基于典型的格點(diǎn)QCD 膠球模擬計(jì)算程序,采用了多項(xiàng)優(yōu)化技術(shù),對(duì)程序進(jìn)行了性能分析和并行優(yōu)化,并實(shí)現(xiàn)了大規(guī)模的并行測(cè)試,并行規(guī)模超過4 萬核.性能方面,加入OpenMP 多線程并行部分,最高取得了3.32 倍的加速比,在向量化優(yōu)化矩陣計(jì)算部分,最高取得了3.17 倍的加速比,在并行輸出組態(tài)文件部分取得了153 倍的加速比.

優(yōu)化后的程序可更有效地產(chǎn)生大量組態(tài),有利于得到更精確的膠球譜.優(yōu)化后的程序可較容易推廣到夸克與膠子組成的混雜態(tài)的研究,及普通介子和膠球的混合問題的研究,為高能物理的相關(guān)研究提供更加便利的軟件工具.

致謝

感謝中國(guó)科學(xué)院高能物理研究所理論物理室陳瑩老師參與討論和提供幫助,感謝廣州“天河2 號(hào)”超算工作人員提供的測(cè)試支持.