基于Intel Xeon Phi的激光等離子體粒子模擬研究＊

姚文科，杜云飛，吳 強(qiáng)，楊燦群

（1.國(guó)防科學(xué)技術(shù)大學(xué)并行與分布處理國(guó)家重點(diǎn)實(shí)驗(yàn)室，湖南長(zhǎng)沙410073；2.國(guó)防科學(xué)技術(shù)大學(xué)計(jì)算機(jī)學(xué)院，湖南長(zhǎng)沙410073）

1 引言

激光等離子體粒子模擬通常用于描述強(qiáng)激光與等離子體相互作用的非線性過(guò)程，從而預(yù)測(cè)、設(shè)計(jì)、解釋實(shí)驗(yàn)結(jié)果［1，2］。它對(duì)于解釋激光等離子體多種耦合過(guò)程的物理圖像細(xì)節(jié)，發(fā)現(xiàn)新的物理現(xiàn)象具有重要的意義［1～4］。在實(shí)際應(yīng)用中，模擬迭代的時(shí)間步通?？蛇_(dá)數(shù)千萬(wàn)步，模擬的粒子數(shù)目超過(guò)數(shù)百萬(wàn)個(gè)［5，6］，為了在有效的時(shí)間內(nèi)獲得模擬結(jié)果，模擬對(duì)計(jì)算性能有著極高的要求。因此，使用高性能計(jì)算領(lǐng)域的最新成果加速激光等離子體粒子模擬具有重要的現(xiàn)實(shí)意義。

近年來(lái)，異構(gòu)系統(tǒng)已經(jīng)成為高性能計(jì)算領(lǐng)域的一種主流系統(tǒng)［2，7］。通過(guò)集成通用CPU和加速器單元如GPU等，異構(gòu)系統(tǒng)在功耗和性能方面占有優(yōu)勢(shì)。針對(duì)GPU等加速單元與CPU體系結(jié)構(gòu)差異大、編程模型不兼容的問題［7］，Intel發(fā)布了支持X86指令架構(gòu)的Intel Xeon Phi協(xié)處理器［8］。Intel Xeon Phi一方面可以兼容原CPU上的程序庫(kù)、開發(fā)工具以及應(yīng)用程序，另一方面對(duì)于特定類型的計(jì)算任務(wù)，可以提供極高的計(jì)算性能，如Intel Xeon PhiTMCoprocessor 5110P單卡雙精度浮點(diǎn)性能可達(dá)1.01 Tflops［2］。在2013年5月Top 500排名中，采用了CPU－Intel Xeon Phi異構(gòu)結(jié)構(gòu)的TH－2超級(jí)計(jì)算機(jī)排名榜首?？梢灶A(yù)見，基于CPU－Intel Xeon Phi的異構(gòu)系統(tǒng)將在未來(lái)高性能計(jì)算領(lǐng)域中發(fā)揮重要的作用。

基于這一背景，本文研究了三維等離子體粒子模擬程序LARED－P到Intel Xeon Phi的移植。首先采用Native模式對(duì)LARED－P程序中熱點(diǎn)計(jì)算任務(wù)進(jìn)行優(yōu)化研究，通過(guò)采用SIMD擴(kuò)展指令進(jìn)行手工向量化以及數(shù)據(jù)預(yù)取、循環(huán)展開的優(yōu)化，使該計(jì)算任務(wù)獲得了4.61倍的加速。隨后采用Offload模式將程序移植到CPU－Intel Xeon Phi異構(gòu)系統(tǒng)上，并通過(guò)使用異步數(shù)據(jù)傳輸和雙緩沖技術(shù)提升了程序性能。

本文結(jié)構(gòu)如下：第2節(jié)介紹了Intel Xeon Phi的結(jié)構(gòu)特點(diǎn)和LARED－P算法，回顧了相關(guān)工作；第3節(jié)描述了Native模式下LARED－P在Intel Xeon Phi上的SIMD優(yōu)化；第4節(jié)展現(xiàn)了Offload模式下粒子運(yùn)動(dòng)方程的優(yōu)化過(guò)程；第5節(jié)給出了性能測(cè)評(píng)結(jié)果；第6節(jié)對(duì)全文進(jìn)行了簡(jiǎn)要總結(jié)。

2 背景介紹

本節(jié)主要對(duì)Intel Xeon Phi加速器的結(jié)構(gòu)和LARED－P程序結(jié)構(gòu)進(jìn)行簡(jiǎn)單介紹；同時(shí)，對(duì)其他學(xué)者的類似研究工作進(jìn)行總結(jié)。

2.1 Intel Xeon Phi的結(jié)構(gòu)

Intel在2012年11月發(fā)布了Intel Xeon PhiTM協(xié)處理器，它基于Intel MIC（Many Integrated Core）架構(gòu)。每片Intel Xeon Phi可以提供62個(gè)x86兼容的計(jì)算核心［8］，每個(gè)核包含一個(gè)支持512位SIMD的向量處理單元（VPU）。通過(guò)使用VPU，每個(gè)核可同時(shí)處理8路雙精度或16路單精度浮點(diǎn)運(yùn)算［8］，目前，Intel通過(guò)提供Intel向量化庫(kù)、Intel Intrinsic函數(shù)等方式使用VPU。Intel Xeon Phi的每個(gè)核擁有32 KB的L1數(shù)據(jù)Cache和32 KB的L1指令Cache。此外，每個(gè)核還可以使用512 KB的L2 Cache。不同核的L2 Cache通過(guò)雙向內(nèi)存控制器相連。

Intel提供了多種模式來(lái)使用Intel Xeon Phi，例如Native模式、Offload模式等。

Native模式將Intel Xeon Phi作為一個(gè)獨(dú)立的處理節(jié)點(diǎn)，其上運(yùn)行有微Linux操作系統(tǒng)，內(nèi)存空間作為全部的存儲(chǔ)空間。

Offload模式則將Intel Xeon Phi作為CPU的加速單元。CPU作為主控端，將計(jì)算密集部分傳送到Intel Xeon Phi上進(jìn)行計(jì)算，然后再將計(jì)算結(jié)果返回給CPU。

2.2 LARED－P簡(jiǎn)介

LARED－P程序是北京應(yīng)用物理與數(shù)學(xué)研究所采用粒子云網(wǎng)格法［3，4］實(shí)現(xiàn)的三維等離子體粒子模擬程序。它適應(yīng)于二維和三維Cartesian坐標(biāo)系，本文主要針對(duì)三維模擬進(jìn)行研究。它采用均勻矩形結(jié)構(gòu)網(wǎng)格離散整個(gè)模擬空間，等離子體采用粒子云近似。

LARED－P程序流程如圖1所示，其中主要的三個(gè)任務(wù)為求解電磁場(chǎng)Maxwell方程、靜電勢(shì)Poisson方程和粒子運(yùn)動(dòng)方程。在計(jì)算過(guò)程中，LARED－P以有限差分方法離散Maxwell方程，然后采用蛙跳格式離散粒子運(yùn)動(dòng)方程，最后采用粒子云網(wǎng)格方法（PIC）耦合粒子和電磁場(chǎng)。這三個(gè)方程中粒子運(yùn)動(dòng)方程占據(jù)了計(jì)算時(shí)間的31.8%，粒子云方程占34.7%，粒子重排過(guò)程占27%。

Figure 1 Flow chart of LARED－P圖1 LARED－P程序流程圖

2.3 相關(guān)工作

莫?jiǎng)t堯等人［4］在所開發(fā)的串行LARED－P程序的基礎(chǔ)上引進(jìn)網(wǎng)格與粒子的關(guān)系鏈表，重排粒子位置，消除了性能遷移，提高了單機(jī)性能。而后又開發(fā)了基于MPI環(huán)境的三維粒子模擬程序LARED－P，并采取三維分割、合理安排通信結(jié)構(gòu)，設(shè)置虛網(wǎng)格等手段進(jìn)行優(yōu)化。

劉來(lái)國(guó)等人［6］完成了基于GPU的LARED－P加速技術(shù)的研究與實(shí)現(xiàn)。加速后的LARED－P程序在NVIDIA Tesla S1070的單個(gè)GPU上獲得了相當(dāng)于主頻2.4 GHz的Intel（R）Core（TM）2 Quad CPU Q6600單核的6倍加速比。

石志才等人［9］則針對(duì)LARED－P中粒子運(yùn)動(dòng)方程在GPU上的實(shí)現(xiàn)提出八色分解方案，消除并行訪問的讀寫依賴，擴(kuò)大了并行度，在NVIDIA Tesla M2090上獲得了相對(duì)于Intel Xeon X5670 CPU單核的20.08倍的加速比。

3 Native模式下基于SIMD的優(yōu)化

本節(jié)介紹在Native模式下針對(duì)LARED－P所做的優(yōu)化工作。該算法的實(shí)現(xiàn)部分見文獻(xiàn)［9］。運(yùn)動(dòng)方程的運(yùn)行時(shí)間占整個(gè)算法運(yùn)行時(shí)間的31.8%，在移植過(guò)程中為了提高其數(shù)據(jù)級(jí)并行，對(duì)它進(jìn)行了向量化。在使用編譯器選項(xiàng)－O3－mmic與－O3－novec－mmic對(duì)比中，發(fā)現(xiàn)兩者運(yùn)行時(shí)間沒明顯區(qū)別；加編譯指導(dǎo)＃pragma vector命令后也沒太大提升。所以，我們選擇了手工向量化。

3.1 基于SIMD擴(kuò)展指令的手工向量化

通過(guò)調(diào)用Intel提供的Intrinsic接口函數(shù)對(duì)粒子運(yùn)動(dòng)方程的求解進(jìn)行手工向量化，該過(guò)程主要包括數(shù)據(jù)打包、向量計(jì)算以及數(shù)據(jù)解包等過(guò)程。以更新粒子位置為例，如圖2所示。其中，x、y、z三個(gè)數(shù)組分別存放粒子的x坐標(biāo)、y坐標(biāo)以及z坐標(biāo)，vxj、vyj、vzj三個(gè)數(shù)組則存放粒子新速度的三個(gè)分量。

這樣經(jīng)過(guò)SIMD向量化，程序局部性增強(qiáng)，性能得到可觀的提升。

在更新粒子新速度時(shí)，涉及許多乘、加操作，為了提高效率可以使用乘加指令：＿mm512＿fmadd＿pd（＿mm512d，＿mm512d，＿mm512d）。乘加指令的使用帶來(lái)了可觀的性能提升。

3.2 預(yù)取優(yōu)化

Figure 2 Sketch map of vectorization圖2 向量化示意圖

512位寬度的向量處理給訪存帶來(lái)了巨大挑戰(zhàn)，為了降低訪問延遲，本文采用了預(yù)取指令mm＿prefetch（），減少Cache失效。

每次迭代訪問八個(gè)網(wǎng)格，這種訪問模式可預(yù)測(cè)性給預(yù)取帶來(lái)了便利。本文測(cè)試了幾組不同預(yù)取距離的運(yùn)行結(jié)果，如圖3所示。從圖3可以看出預(yù)取距離為一個(gè)迭代時(shí)性能最好。

Figure 3 Result of prefetch圖3 預(yù)取效果圖

4 基于Offload模式的程序移植

隨后我們利用offload模式將LARED－P移植到CPU－Intel Xeon Phi異構(gòu)系統(tǒng)上。我們將熱點(diǎn)計(jì)算任務(wù)—運(yùn)動(dòng)方程的求解使用Intel Xeon Phi進(jìn)行加速，程序的其余部分仍然在CPU上完成。在移植的過(guò)程中發(fā)現(xiàn)，CPU和Intel Xeon Phi之間的數(shù)據(jù)傳輸是程序性能瓶頸，占總時(shí)間60%。因此，如何降低傳輸開銷成為了本節(jié)工作的重點(diǎn)。我們主要采用了兩種優(yōu)化方法：異步數(shù)據(jù)傳輸和雙緩沖。

4.1 異步數(shù)據(jù)傳輸

當(dāng)使用Offload模式時(shí)，CPU在將數(shù)據(jù)傳給Intel Xeon Phi之后，由于數(shù)據(jù)依賴的原因，需要等待其計(jì)算完畢并將結(jié)果傳送回CPU以后，CPU才能繼續(xù)執(zhí)行。為了減少CPU的空轉(zhuǎn)時(shí)間，我們采用了異步數(shù)據(jù)傳輸?shù)姆绞綄?duì)其進(jìn)行優(yōu)化。

首先采用一種CPU提前計(jì)算的方式減少CPU的空轉(zhuǎn)時(shí)間。當(dāng)需要調(diào)用Intel Xeon Phi進(jìn)行粒子運(yùn)動(dòng)方程加速時(shí)，CPU以異步方式將數(shù)據(jù)傳輸?shù)絀ntel Xeon Phi上，并提前進(jìn)行與粒子運(yùn)動(dòng)方程計(jì)算結(jié)果無(wú)關(guān)的任務(wù)的處理，如圖4所示。然而，由于大部分的計(jì)算任務(wù)都與運(yùn)動(dòng)方程相關(guān)，因此這種方式取得的效果并不明顯。

隨后采用了一種提前傳輸?shù)姆绞絹?lái)減少數(shù)據(jù)拷貝時(shí)CPU的等待時(shí)間。當(dāng)粒子運(yùn)動(dòng)方程相關(guān)的數(shù)據(jù)已經(jīng)準(zhǔn)備完畢時(shí)，CPU立刻以異步的方式將數(shù)據(jù)傳給Intel Xeon Phi后，繼續(xù)以正常流程進(jìn)行任務(wù)處理，這樣增加了CPU和Intel Xeon Phi并行執(zhí)行的時(shí)間，從而隱藏了部分?jǐn)?shù)據(jù)傳輸時(shí)間。

Figure 4 Sketch map of asynchronous mode圖4 異步傳輸示意圖

4.2 雙緩沖

為了進(jìn)一步減少數(shù)據(jù)傳輸開銷，我們采用了雙緩沖技術(shù)，在CPU端開辟兩個(gè)發(fā)送緩沖區(qū)和兩個(gè)接收緩沖區(qū)。所需數(shù)據(jù)分組傳輸，相比之前MIC端需要接收全部數(shù)據(jù)之后才開始計(jì)算，在使用了雙緩沖之后，帶來(lái)的明顯好處是：以流水線方式執(zhí)行傳輸和計(jì)算，可以將絕大部分?jǐn)?shù)據(jù)傳輸時(shí)間隱藏，如圖5所示。

Figure 5 Sketch map of double buffer圖5 雙緩沖示意圖

5 性能評(píng)測(cè)

本節(jié)對(duì)LARED－P算法進(jìn)行性能測(cè)試，測(cè)試平臺(tái)的具體配置如表1所示。我們采用的編譯器為Intel C＋＋Compiler 13.0，MPI版本為3.0。

Table 1 Configuration of the test platform表1 測(cè)試平臺(tái)的配置

首先測(cè)試了手動(dòng)向量化和數(shù)據(jù)預(yù)取取得的性能提升。測(cè)試網(wǎng)格數(shù)為64＊128＊128，平均每個(gè)網(wǎng)格包含10個(gè)粒子。圖6為沒有向量化、只有向量化和向量化加預(yù)取三種情況下程序的運(yùn)行時(shí)間，從圖6中可以看出，向量化加速效果不錯(cuò)，加速比可達(dá)到4.61，而預(yù)取的加速效果沒有向量化明顯，但是仍有19.7%的性能提升。因?yàn)镸IC卡內(nèi)存有限，無(wú)法再擴(kuò)大規(guī)模，而當(dāng)線程數(shù)接近64時(shí)，運(yùn)行時(shí)間占整個(gè)程序時(shí)間的比例已經(jīng)非常小了，所以只測(cè)試到了64線程。CPU上運(yùn)行16線程時(shí)總時(shí)間為1.058 s，而MIC上64線程運(yùn)行時(shí)間為0.717 s，所以在MIC上每個(gè)核運(yùn)行一個(gè)線程的情況下，相對(duì)單塊CPU有1.48倍的加速。

Figure 6 Result in native mode圖6 Native模式下實(shí)驗(yàn)結(jié)果

接著測(cè)試了Offload模式下的異步輸入和雙緩沖的效果，如圖7所示。從圖7中可以看出，異步輸入和雙緩沖都對(duì)運(yùn)動(dòng)方程產(chǎn)生了一定的性能提升，相比于傳統(tǒng)Offload分別有9.8%和21.8%的性能提升。

6 結(jié)束語(yǔ)

本文基于Intel Xeon Phi協(xié)處理器實(shí)現(xiàn)了LARED－P程序的移植。在移植的過(guò)程中，本文綜合運(yùn)用了Native和Offload兩種編程模式：首先運(yùn)用Native模式對(duì)LARED－P程序中熱點(diǎn)計(jì)算任務(wù)進(jìn)行優(yōu)化研究，通過(guò)采用SIMD擴(kuò)展指令使該計(jì)算任務(wù)獲得了4.61倍的加速；然后運(yùn)用Offload模式將程序移植到CPU－Intel Xeon Phi異構(gòu)系統(tǒng)上，并通過(guò)使用異步數(shù)據(jù)傳輸和雙緩沖技術(shù)分別獲得9.8%和21.8%的性能提升。

Figure 7 Result in Offload mode圖7 Offload模式下實(shí)驗(yàn)結(jié)果

［1］ Zhu Shao－ping，Zhang Wei－yan.Overview of computer simulation on laser fusion in China［J］.Journal of the Korean Physical Society，2006，49：33－38.

［2］ Chang Tie－qiang.Laser－plasma interaction and laser fusion［M］.Changsha：Hunan Science and Technology Press，1998.（in Chinese）

［3］ Ma Yan－yun，Chang Wen－wei，Yin Yan，et al.An object－oriented 3D parallel simulation program PLASIM3D［J］.Chinese Journal of Computation Physics，2004，21（3）：305－311.（in Chinese）

［4］ Mo Ze－yao，Xu Lin－bao，Zhang Bao－lin，et al.Parallel computing and performance analysis for 2D plasma simulation with particle clouds in cells method［J］.Chinese Journal of Computation Physics，1999，16（5）496－504.（in Chinese）

［5］ Cao Xiao－lin，Zheng Chun－yang，Zhang Ai－qing，et al.Program development of 3D plasma simulation oriented thousands of processors［J］.Chinese Journal of Progress in Natural Science，2009，1（5）：544－550.（in Chinese）

［6］ Liu Lai－guo，Xu Wei－xia，Yang Can－qun，et al.Accelerating LARED－P algorithm based on GPU［J］.Computer Engineering ＆Science，2009，31（A01）：59－63.（in Chinese）

［7］ Stantchev G，Dorland W，Gumerov N.Fast parallel particleto－grid interpolation for plasma PIC simulations on the GPU［J］.Journal of Parallel and Distributed Computing，2008，68（10）：1339－1349.

［8］ Intel?Xeon PhiTMCoprocessor System Software Developers Guide［R］.SKU＃328207－001EN，2012.

［9］ Yang Can－qun，Wu Qiang，Hu Hui－li，et al.Fast weighting method for plasma PIC simulation on the GPU－accelerated heterogeneous systems［J］.Journal of Central South University，2013，20（6）：1527－1535.

［10］ MPI－2：Extensions to the message－passing interface［EB/OL］.［2012－05－16］.http：//micro.ustc.edu.cn/Linux/MPI/mpi－20.pdf.

［11］ Schulz K W，Ulerich R，Malaya N，et al.Early experiences porting scientific applications to the many integrated core（MIC）platform［C］∥Proc of the 2012 Highly Parallel Computing Symposium，2012：1.

附中文參考文獻(xiàn)：

［2］ 常鐵強(qiáng).激光等離子體相互作用與激光聚變［M］.長(zhǎng)沙：湖南科學(xué)技術(shù)出版社，1998.

［3］ 馬燕云，常文蔚，銀燕，等.三維面向?qū)ο蟮牟⑿辛Ｗ幽M程序PLASIM3D［J］.計(jì)算物理，2004，21（3）：305－311.

［4］ 莫?jiǎng)t堯，許林寶，張寶琳，等.二維等離子體模擬粒子云網(wǎng)格方法的并行計(jì)算與性能分析［J］.計(jì)算物理，1999，16（5）：496－504.

［5］ 曹小林，鄭春陽(yáng)，張愛清，等.面向數(shù)千處理器的三維等離子體粒子模擬程序研制［J］.自然科學(xué)進(jìn)展，2009，19（5）：544－550.

［6］ 劉來(lái)國(guó)，徐煒遐，楊燦群，等.基于GPU的LARED－P算法加速［J］.計(jì)算機(jī)工程與科學(xué)，2009，31（A01）：59－63.