LQCD Dslash在神威·太湖之光上的研究分析與MPI實現*

張 淼，周 宇，陳建海+，何欽銘，徐 順，宮 明

1.浙江大學 計算機科學與技術學院，杭州310012

2.中國科學院 計算機網絡信息中心，北京100190

3.中國科學院 高能物理研究所，北京100049

+通訊作者E-mail:chenjh919@zju.edu.cn

1 引言

格點量子色動力學（lattice quantum chromodynamics，LQCD）是粒子物理最前沿的研究課題之一。自然界存在四種基本相互作用——強相互作用、弱相互作用、電磁相互作用及引力相互作用。物理學家將強相互作用、弱相互作用及電磁相互作用統一起來稱之為“標準模型”。量子色動力學（quantum chromodynamics，QCD）規范理論是標準模型中描述夸克、膠子組成的粒子強相互作用的基本理論，其重要特征是漸進自由，這使得高能強相互作用過程可利用微擾論進行處理，但此方法針對低能強相互作用過程就無能為力了。1974 年，Wilson 創建了格點規范理論[1]，其將連續時空用離散晶格來替代，夸克放置于格點上，它們之間用規范場鏈建立聯系。此理論不借助任何QCD以外的參數或假設從第一原理出發用非微擾方法來解決上述問題。

在用計算機對LQCD 理論進行數值模擬時需使用蒙特卡羅算法，而這一過程是高度計算密集型的，故現代LQCD 數值模擬過程均使用超級計算機來進行。目前，已經有很多工作在大規模超級計算機上實現了LQCD數值模擬。例如：2006年，Vranas等人[2]基于IBM 的BlueGene/L 超級計算機，在32 000 個CPU核心下，模擬QCD實現了12.2 TFLOPS的性能，并以此獲得了2006年的戈登貝爾獎。2010年，Babich等人[3]在32 個由InfiniBand 互連的GPU 上，采用MPI（message passing interface）并行化QUDA 庫，實現了4 TFLOPS的性能。2011年，Ammendola等人[4]用多核CPU（central process unit）、GPU（graphics process unit）及FPGA（field-programmable gate array）構建了一個由42U racks 組成的3D torus 網絡原型系統——QUonG 用于LQCD 計算；同年，Babich 等人[5]通過使用多維度的并行化策略和域分解的預調節器將由post-Monte Carlo 算法實現的LQCD 擴展到256 個GPU上。2013年，Joó等人[6]將Dslash內核移植到Intel Xeon Phi 加速卡上，單節點取得了280 GFlops 的性能，整個求解器獲得了大約215 GFlops 的性能，當將節點擴展到64 個Intel KNC 時，整個求解器獲得了3.6 TFlops 的性能。2017 年，Bonati 等 人[7]使用OpenACC 實現了可移植的LQCD 蒙特卡羅算法代碼。2018 年，Kanamori 等人[8]將Wilson 內核在Oakforest-PACS超級計算機上實現了并行化。

從“天河一號”初登TOP500榜首，到“天河二號”連續3 年登頂，世界首臺10 億億次超級計算機“神威·太湖之光”連續4次奪冠，我國超算硬件架構已達世界先進水平，軟件水平也有了一定的提升，如在神威平臺上實現的千萬核可擴展全球大氣動力學全隱式模擬[9]和高可擴展非線性地震模擬[10]分獲了2016年與2017 年的戈登貝爾獎，這是歷史性的突破。但目前尚無工作針對LQCD 在這些國內超算平臺上進行移植優化。2018年最新一期TOP500榜單，“神威”已被美國最新的“Summit”超級計算機超越。未來，我國不僅要重奪超算硬件架構水平的制高點，還需加快自身軟件生態的建設。故為提升我國自主超算平臺的軟件水平，針對神威平臺的特有架構，本文工作主要是將LQCD中的熱點代碼Dslash進行移植，通過代碼向量化、循環展開、并行文件讀寫、主-從核DMA（direct memory access）通訊、主-主核MPI 通訊等優化手段，實現了Dslash代碼在申威26010芯片上的異構眾核并行,并實現了一定的優化效果，單核組從核陣列優化取得了165倍的加速比，使用MPI連并16 個核組優化取得了25 倍的加速比，同時發現使用MPI 引發了性能瓶頸，通過實驗定位了瓶頸所在，為后續優化工作打下基礎。

本文組織結構如下：第2章將介紹神威平臺的系統架構與通信模式；第3章詳細闡述了LQCD在神威平臺上的優化過程；第4章對優化效果進行了數據測試并加以分析；第5章對全文內容進行了總結并對未來工作給予展望。

2 神威架構及其通信模式

2.1 系統架構

神威·太湖之光超級計算機采用的是基于申威26010 處理器的異構眾核架構，其處理器架構如圖1所示。每個處理器包含有4個核組（core group，CG），而每個核組又包括一個主核（management processing element，MPE）、64 個 從 核（computing processing element，CPE）組成的8×8計算處理元件陣列，以及一個協議處理單元（protocol processing unit，PPU）和一個DDR3 存儲控制器（memory controller，MC）。4 個核組通過片上網絡（network on chip，NOC）連接，且每個核組有自己的內存空間，主核和從核均通過存儲控制器來訪問。主核是具有256 位向量單元的全功能64 位精簡指令集（reduced instruction set computing，RISC）核心，而從核是精簡過的具有相同長度向量單元的64位RISC核心。每個從核都配置了一個64 KB大小的本地設備存儲器（local device memory，LDM），其可以配置為用戶可編程的高速緩沖區或用于數據自動緩存的軟件模擬緩存。

Fig.1 Architecture of SW26010圖1 申威26010芯片架構

神威·太湖之光超級計算機總共有40 960 個處理器，它們之間通過中央交換機網絡互連，且每256 個處理器組成一個完全連接的超節點。每個申威處理器具有260 個核心，其雙精度峰值性能為3.06 TFlops，整機核心數目超過一千萬個，總體性能達125 PFlops[11]。

2.2 通信模式

神威平臺目前支持三種主要的編程模型，根據應用的特點，本文工作將研究重點放在了MPI+Athread混合編程模型上，其中Athread是由設備供應商提供的面向細粒度眾核編程的輕量級線程庫。利用Athread在從核陣列上進行眾核高效率計算及從核間通訊，MPI用于不同核組的主核間數據傳遞。

基于此，本文工作所涉及的算法和優化策略將主要充分利用申威處理器的三種高級架構特性來獲得高性能，即用戶控制的LDM、直接存儲器訪問（DMA）和寄存器數據通信機制。每個從核上用戶可控制的LDM 可視為私有的用戶可編程存儲空間。DMA 主要用于主存儲器和LDM 之間的高效數據移動，以便用戶可快速從LDM 中訪問數據而不是從外面的主存儲器中。低延遲的片上寄存器數據通信機制使得在同一核組中的從核之間可以進行快速有效的數據交換，而無需從核陣列和主存儲器之間的額外數據移動。利用MPI的消息傳遞機制進行不同核組的邊界數據的傳輸。

3 數學描述與優化方案設計及實現

3.1 數學描述

LQCD 描繪了一個四維(x,y,z,t)的格點陣列，并把每一維的首尾連接起來，如圖2 所示，其為二維的示意圖，類似一個“甜甜圈”。在每個格點上放置表示夸克的費米子場量，其是含有3個色分量和4個旋量分量的格拉斯曼數，這樣每個格點就成為一個旋子，其在計算機中具體操作時以12 個復數構成的向量表示。在相鄰旋子之間的連接上放置表示膠子的規范場量，其可以表示為色空間的3×3復數矩陣。由于是四維空間，且旋子間的連接又分為正、負兩個方向，故每個格點上將會有8 個3×3 的復矩陣，以連接到8個鄰居格點。

Fig.2 2-dimensional schematic diagram of LQCD圖2 LQCD二維示意圖

通過對每個格點附加費米子常量，每個連接附加規范場量來對格點的狀態進行描述。不同的狀態出現的實際概率不盡相同，此概率正比于e-βSg(U)+ψˉM(U)ψ，其中Sg(U)是全部規范場量Uμ(x)的一個泛函，M(U)是與規范場量相關的一個矩陣，其物理含義為全空間中任意兩個局部夸克場之間的相互作用，這二者定義了LQCD在格點上的物理模型。

對四維空間中的所有狀態進行積分，通過歐式空間的路徑積分得到配分函數Z：

任意物理量O(U,Gq(U))的期望值可寫為：

其中，Gq(U)=M-1(U)代表夸克傳播子。

使用蒙特卡羅算法進行式（2）積分的計算僅需兩個步驟：第一步，隨機地生成規范場組態Uμ(x)，并使其出現的概率為P(U)∝e-βSg(U)detM(U)；第二步，計算可觀測量值。其中，由于難以直接計算detM(U)與M-1(U)，因此通過采用隨機偽費米子場的方法來估算detM(U)，從而將其轉化為求解M-1(U)，即稀疏矩陣的求逆操作。M-1(U)的物理含義是為了表示任意兩個局部夸克場之間的傳播振幅，對于求解M-1(U)，每次只需計算它的若干列而非整個逆矩陣，問題便轉化為求解線性方程組：

其中，b表示四維格點空間，通過x的求解，方便對M-1(U)的計算。

對于求解上述方程組，只需在{b,Mb,M2b,M3b,…}張開的Krylov子空間中迭代逼近即可達到滿足要求的精度。在如共軛梯度法等Krylov 子空間算法中，關鍵的操作是矩陣乘以向量和向量內積的計算，其中矩陣向量乘是主要的計算熱點，而向量內積涉及全局規約操作。

對于矩陣向量乘主要是計算以下的公式，即Wilson費米子矩陣：

其中，μ表示4 個時空方向；x和y是兩個4-矢量時空坐標；Uμ(x)表示從x點向μ方向的規范場連接，其為3×3的復矩陣，稱為“組態”，U?是它的轉置復共軛；γμ是4×4的稀疏復矩陣，形式為：

此費米子矩陣對角元部分的計算相對簡單，非對角元上的計算才是關鍵部分，定義：

Wilson費米子矩陣每經過一次迭代計算，空間中的所有旋子就會更新自身的狀態，上式定義的Dslash即為重點優化的熱點部分。

3.2 代碼結構

本文重點優化的包括Dslash 熱點部分在內的Wilson費米子矩陣主要由3個函數所構成，它們的依賴關系如圖3 所示（箭頭指向的是被調用者）。這3個函數分 別 為mr_solver(…)、m_wilson(…)、dslash(…)，它們的代碼結構如下：

Fig.3 Dependences of 3 important functions圖3 3個主要函數的依賴關系

LQCD的整個計算流程可簡化為：在求解線性方程組的迭代條件下，將表示費米子場量和規范場量的復矩陣通過Wilson費米子矩陣調用Dslash運算以遍歷更新所有格點的數據，然后進行規約運算，最后計算殘差，輸出結果解。如圖4展示了上述基本流程。

Fig.4 Flow of LQCD圖4 LQCD計算流程

3.3 單核組優化方案

從式（6）可以看出，Dslash 運算可以歸為一種四維的9-points stencil運算，如圖5所示為二維示意，即每計算一個格點的數據均需要訪問4個維度上相鄰8個格點的數據及8個方向上的連接的數據。

Fig.5 2-dimensional schematic diagram of 4-dimensional 9-points stencil圖5 四維9-points stencil的二維示意圖

關于stencil 計算，目前也已有很多的研究成果。如：Datta等人[12]針對stencil計算開發了一套行之有效的優化策略，并以此構建了一個自動調優的機制，其可通過搜索優化策略組及其參數來選出最優的策略，從而最大限度地減少運行時間，同時提高性能。此機制在Intel Clovertown、IBM QS22 PowerXCell 8等多核架構上進行了驗證，并取得了很好的效果。Ho等人[13]通過利用本地存儲器和異步軟件預取技術來優化三維Lattice Boltzmann（LBM）stencil 模型，并在Kalray MPPA-256 多核處理器上獲得了33%的性能提升。Bondhugula 等人[14]提出了一種全新的平鋪技術——Diamod Tiling，其解決了現有自動平鋪框架流水啟動和負載不均衡的問題，使得可以并發啟動且保證負載均衡，從而最大化stencil 計算的并行化。Ao 等人[15]在算法層面提出了一種計算-通信重疊技術，以減少進程間的通信開銷；一種局部感知的阻塞方法，以增強數據局存性，從而充分利用片上并行性；以及一種針對不同線程中共享數據的聯合訪問模式。以此，他們將三維非靜力大氣模擬模型在神威·太湖之光上進行了實現，并獲得了26 PFlops 的stencil計算性能。

根據Dslash運算的特點，目前實現的Dslash代碼是根據從核LDM 的大小專門設計的，希望盡可能地把每個從核需要處理的數據存放在它上面的64 KB大小的LDM 上，從而最小化主存儲DMA 開銷。在此設計下，每個核組負責84大小的子格子空間，每個從核負責64 個格點，且這64 個點采用指定z和t的一個xy平面，如圖6所示。

3.3.1 代碼向量化及循環展開

Fig.6 Schematic diagram of 4-demensional lattice data spreading on CPE cluster圖6 四維格點數據在從核陣列上的分布示意圖

在實際開發過程中，發現申威處理器針對簡單循環以外的其他結構并沒有提供自動向量化（如SLP（successive linear programming）算法）功能，因此在無法利用編譯器進行向量化的情況下只能手動向量化代碼。申威處理器每次能夠處理4個浮點數值，且旋量的4個分量剛好可以一并處理，于是修改了數據的存儲格式，將旋量指標放到最內層，這樣矩陣元素在計算時被逐個處理，轉置也不會帶來額外的開銷。Dslash中有8項不同的操作，不同的項需要考慮轉置或不轉置、1±γμ的符號等情況進行分別處理。另外，SU(3)矩陣有18個浮點數，而向量的存儲是4個浮點數對齊的，這就意味著如果某個方向上的矩陣是對齊的，那么隨后的一個就剛好不是對齊的，因此對齊與不對齊的情況也要分別處理。通過使用宏嵌套定義來分別對8個項生成對應的函數定義，從而對數據對齊或不對齊等各種不同的情況分別進行處理。

除了Dslash以外，其他涉及旋量場操作的函數也進行了向量化，主要是把循環內的操作改寫成向量的形式。

指令流水線是否一直通暢并不斷發射指令是影響從核計算效率的重要因素之一。循環是計算中最常見的結構，每次循環都需要進行分支判斷，而條件判斷會采用分支預測技術，但申威編譯器的分支預測功能尚不健全，若直接采用其默認的分支預測功能就會造成流水線的不通暢，因此直接將Dslash中的64個時空點相關的循環全部展開。代碼徹底展開后增加了代碼的靈活性，其中有關和主存進行DMA傳輸的代碼部分是異步的，這樣就可以通過巧妙地安排代碼次序來實現隱藏通訊時間。

3.3.2 從核寄存器通訊機制

從核寄存器通訊機制是申威處理器的一個非常特殊的機制。在實踐過程中，感覺此機制使得從核上的程序并不像線程，而更接近于進程，這就直接改變了很多編程的基本思路。寄存器通訊只能發生在同行或同列的從核上，因此安排每個從核上的計算任務時，這種通訊關系必須滿足。也就是說，對于任何一個時空點而言，相鄰的8個時空點要么在同一個從核上，要么在同行或同列的其他從核上。若不滿足這個條件，數據傳輸就需要中轉，代碼因此會變得很復雜，效率也會大大降低。如前所述，在目前的版本中，采用簡單的分配方案，將從核的行和列映射到四維模型的z和t維度，然后給每個從核分配一個xy平面的8×8的薄片。此方案的好處是代碼清晰直觀，但效率不是最好的，一方面由于薄片的面積較大，另一方面是與主存通信的DMA 數量不平衡。在后續版本中，將采用參數化代碼自動生成方案，為每個從核生成其特定的任務代碼，使得從核間協調工作，并使用遺傳算法、模擬退火算法等隨機算法來優化從核的代碼設計，使得從核上的任務操作順序最優。基于薄片方案，在開發過程中了解到寄存器通信涉及到3 套隊列，分別是“發送”“路上”“接收”，且它們有各自不同的大小和共用規則，非常復雜。如果接收從核的接收隊列滿了，路上隊列也滿了，那么有些數據就會停留在發送從核的發送隊列中，在同步過后，不同從核的隊列轉移進路上隊列的次序并不能保證是最初的發送次序，從而導致次序錯亂。這個原理弄清楚后，讓后一個發送從核與當前接收從核同步，保證發送的時候接收隊列是空的，從而確保發送次序的穩定。

3.3.3 主-從核DMA異步通信

從核LDM與主存之間的數據傳輸采用的是異步DMA。神威架構定制的Athread 庫包含了athread_get()等API（application programming interface）接口，這些接口一般都是一些宏，包裝了編譯器實現的“內置函數”，如dma()等。這些內置函數一般會被發射為若干行匯編代碼，大致相當于嵌入匯編。這三層接口的效率由慢到快，在程序的不同部分會使用不同層級的接口，其中對計算要求較高的地方采用直接嵌入匯編代碼或內置宏接口的方式調用所需功能。在DMA 操作方面，經測試Athread 的效率明顯較低，因此直接采用dma_set_*()系列宏接口來實現DMA 傳輸。在主核和從核相互傳輸數據時，雙方需要保證對“進度”的統一認知，這可以通過在主存上設置信號量來實現。具體編寫代碼時，無論是普通變量還是指針變量，都需要用volatile 關鍵字來修飾（其中對于指針變量的修飾，放置位置不同含義也不同），以保證緩存一致性。對于信號量的更新操作，申威處理器提供了3個指令，它們各自有相應的包裝宏來調用，在實際使用時，重寫了這些宏，以保證代碼正確性。

3.4 基于狀態機模型的多核組MPI通訊

在整個LQCD程序優化過程中，設計的整體思路是從核負責計算，主核負責通信。故在從核優化的基礎上，使用MPI消息傳遞機制來進行核組間通信，以擴展并行規模。

核與核之間的通信主要出現在熱點計算dslash函數中，每個核組負責存儲84數量的格點，并負責本區域內格點的計算。邊界從核在對四維空間邊界上的格點進行計算時，需要從存儲相鄰格點的相鄰核組中的從核中獲取數據，同時也需要將自身的邊界數據發送給對方。由于是四維空間，每個維度上都需要交換數據，且分前、后兩個方向，那么就需要8次發送數據與8次接收數據，每次發送和接收的數據量為四維空間某個邊界上的數據，即83數量的格點。整個計算需要多次執行dslash函數，因此相應的主核之間也需要進行多次MPI通信。

3.4.1 建立笛卡爾通信域

由于每個核組對應著四維空間的數據，那么核組之間就有了空間關系，為了更好地表示這種空間位置關系，建立了對應維度和維度長度的MPI 笛卡爾四維torus通信域，即此通信域的維度為4(x、y、z、t)，且每個維度有其對應的進程數。這樣在通信中就可以很方便地通過進程自身在笛卡爾通信域中的坐標（MPI_Cart_coords）找到其對應邊界相鄰進程的坐標（MPI_Cart_shift）。

3.4.2 非阻塞輪詢監聽及持久通信

在一次dslash 執行過程中，需要進行8 次數據發送和接收操作，如果使用阻塞函數（MPI_Send/MPI_Recv）來進行發送和接收，那么就需要很繁瑣的設計來避免進程間的死鎖。同時存放在從核陣列中的格點數據的狀態在同一時刻不盡相同，因此主核并不能隨時對邊界數據進行發送，亦不能隨時提供從核所需要的邊界數據?；诖耍瑢Πl送和接收采用非阻塞函數（MPI_Isend以及MPI_Irecv）[16]，并考慮到需要多次進行dslash計算，對非阻塞發送和接收采用持久通信方式（MPI_Send_init 和MPI_Recv_init）以減小通信開銷。但是采用非阻塞方式，就需要等待或者查看非阻塞函數是否完成，還涉及到主從核如何彼此相互通知數據已交予對方，因此采用模擬網絡socket 通信中的I/O 復用機制select 方式來輪詢監聽發送和接收數據事件。算法如下：

將8 個方向的數據的發送和接收看成是16 個事件（8 個發送事件、8 個接收事件）并進行監聽。主核層面設立8種狀態來進行事件描述，分別為發送過程的4種狀態NOT_READY（從核并沒有準備好主核需要發送的數據）、NOT_SEND（從核數據已準備好，但主核尚未發送）、SENDING（數據發送中）、SEND_SUCCESS（數據發送成功）和接收過程的4 種狀態NOT_RECV（還未接收數據）、RECVING（接收中）、RECV_READY（主核接收完成，等待從核接收）、RECV_SUCCESS（從核接收完成）。

每次輪詢前對發送和接收事件的狀態進行初始化。輪詢時，對每個事件查看其狀態，并執行相應的操作。若某個事件狀態對應的動作未執行完，那么就立即查看下一個事件的狀態；若某個事件的狀態對應的函數或者變量發生變化，那么就改變此事件的狀態并查看下一個事件的狀態。發送和接收事件的狀態都顯示為成功時即標識著一次輪詢的結束。

另外主從核之間需要特定標志來判斷從核已準備好主核需要發送的數據，主核已接收到從核所需的數據。在設計方案中為這些事件各設立了一個變量，用以表示主核需發送的數據有多少從核已經放置到主核上，以及主核已經接收到數據并通知從核進行獲取。

主核間除上述的通過MPI 通訊來傳輸格點數據，還需在求解殘差過程中對每個核組中的殘差數據進行全局規約。

3.5 并行文件讀寫

MPI-2 中提供了一套并行I/O 接口，實現了進程間并行讀寫文件的方法。LQCD 中每個進程對應的數據在數據文件中并不是連續存放的。進程數據分布方式為循環塊分布，通過建立新的數據類型構造符（MPI_Type_create_darray）來定義分布式數組類型。打開文件之后（MPI_File_open），設置進程的文件視圖（MPI_File_set_view），使得每個進程通過各自的文件視圖訪問數據時，文件中只包含視圖對應的數據，且數據之間是沒有空隙的。

4 性能測試與分析

為測試上述兩部分的并行優化效果，現分別設置幾組實驗進行驗證。首先，測試從核優化版本的效果；其次，測試MPI版本的擴展性效果。測試標準為熱點計算部分（即涉及Dslash 運算的上述3 個函數）迭代150次的運行時間。

4.1 單主核版本與單核組從核優化版本

基于84 大小的格點數量，首先將程序在單個主核上編譯成功并運行，這相當于純CPU 版本；然后，將使用從核優化策略改進過的異構并行版本在單個核組，即64 個從核上運行。它們的運行時間如表1所示。

Table 1 Time comparison of single MPE version and single CG version表1 單主核版本與單核組從核版本運行時間對比

從以上實驗數據可以看出，目前的從核優化策略取得了不錯的效果，從核優化版本相較單主核版本的加速比達到了165倍。

4.2 單主核版本與從核優化MPE版本

在從核優化取得不錯加速效果的基礎上，為測試使用MPI 擴展核組規模的性能情況，將程序規模擴展到16 個核組上，每個核組負責的格點數據量仍為84，則總數據量為16×84。與之對應的單主核版本的數據量也擴大到16×84。它們的運行時間如表2所示。

Table 2 Time comparison of single MPE version and MPI version表2 單主核版本與MPI版本運行時間對比

從以上實驗數據可以看出，MPI版本相較單主核版本依然有不錯的運行速度提升，但相比上一節得出的加速比卻有了較大的下降。由于單主核版本的數據量增大了16 倍，理論上運行時間也應增大16倍，實際的數據也符合這一情況（57.73/3.31=17.44,接近16）。但MPI版本相較單核組從核版本而言，整體數據量雖增大了16 倍，但每個核組負責的數據量沒有變，且處理器的數量也相應擴大了16倍，那么理想情況下這兩個版本的運行時間應相差不大，然而實際的數據卻相差較大，MPI 版本慢了99.1%。這說明，MPI 版本在實際運行過程中出現了性能瓶頸，具體原因需要進一步探究。

進一步地，測試了隨著核組規模的擴大，運行時間的變化趨勢，如表3和圖7所示。

Table 3 Running time of MPI version of different CG numbers表3 不同核組數MPI版本的運行時間

Fig.7 Trend of CGs scale and running time圖7 核組規模與運行時間的變化趨勢

從圖7 可分析得出，使用MPI 擴展核組規模后，程序運行時間與單核組版本相比都有較大的差距，且隨著核組規模的擴展，運行時間有增大的趨勢。這進一步說明了使用MPI 誘發了性能瓶頸，且隨著核組規模的擴大，程序的整體性能有下降的趨勢。具體原因將在下一節中詳細討論。

4.3 從核代碼拍數與MPI代碼拍數對比

為了詳細說明使用MPI 在主核間傳輸數據引發性能下降的原因，需要更準確地分析各部分程序運行的時間。使用神威平臺提供的程序運行拍數計數接口rpcc進行測量，將拍數乘以處理器的主頻即可準確得到該段程序運行的墻鐘時間（wall clock time）。

利用神威平臺提供的計數器接口，分別在從核關鍵代碼和MPI代碼中進行插樁，分別得到這些代碼部分的拍數統計，測試數據如表4和表5所示。

Table 4 rpcc results of CPE codes表4 從核代碼拍數測試結果

Table 5 rpcc results of MPI codes表5 MPI代碼拍數測試結果

表4中，m_wilson_dslash 表示從核Dslash熱點部分計算過程的拍數；slave_doit 表示整個從核代碼運算過程的拍數。

表5中，在3 408次輪詢下，共進行了27 264次的MPI 發送（send）和接收（recv）測試，它們之間的關系為TPC=TSC(orTRC)/8。在此基礎上，測得單次輪詢過程需17 981 513拍數，其中一次調用MPI發送指令需134 716 拍數，測試數據發送完成的過程需8 131拍數，且需大約9 次測試完成此過程；接收過程同理。從核發送到主核上的數據進行轉化需1 791 068拍數。它們之間的關系為：

從以上兩個表格數據中可以看出，MPI代碼運行的平均拍數比從核代碼運行的平均拍數大了3 個數量級，這就充分說明了MPI 代碼占用了大量的運行時間，且從核代碼的運算不能將其隱藏，這是目前LQCD 神威版本最大的性能瓶頸。再具體地分析，MPI 運行拍數中的Data transfering 部分又占據了大量的運行時間，其占比達到了895 534×16/17 981 513×100%=79.7%，這部分的功能主要是為了主核間傳輸數據，將從核上準備好的數據通過DMA傳輸給對應的主核，然后主核根據接收的次序對數據進行重排。

經過實驗數據測試與分析，單核組從核優化取得了很好的效果，使用MPI 進行規模擴展后亦取得了一定的效果，但同時也造成了性能下降，通過利用拍數測量工具對代碼的精確分析，將性能瓶頸定位到了主核對從核數據進行重排轉化的過程，這促使下一步的工作需要改進這部分的算法設計。另一方面，減去這部分的性能障礙，MPI 傳輸過程的拍數也要比從核代碼的拍數高出兩個數量級，因此優化MPI傳輸也是下一步工作的重點。

5 總結與展望

LQCD 是當今高能物理領域最前沿的研究方向之一，在大規模超級計算機上實現其數值模擬，對于一個國家占領科技制高點具有重要意義。在本文中，總結了目前在神威·太湖之光超級計算機上對LQCD 中的Dslash 熱點部分所做的移植和并行優化工作。主要貢獻包含如下幾點：通過分析LQCD的應用特點及數值特征，首次將其在神威平臺上實現了成功移植及運行；通過使用向量化、指令流水線、寄存器通訊機制等手段在申威26010 處理器上實現了異構眾核并行，并實現了不錯的加速比；在實現從核陣列并行化的基礎上，進一步使用MPI 實現了多核組連并運行，以此實現一定的并行規模，并著重分析了在目前版本中其成為性能瓶頸的原因。將LQCD在全自主國產超算平臺上實現并行化進行了有力嘗試。測試數據表明，神威平臺對實現大規模的LQCD并行化具有巨大的潛力，對我國未來在高能物理強相互作用領域的數值模擬能力建設與提升具有重大意義。

接下來將會在以下幾方面繼續進行研究和開發：

（1）在神威平臺上對LQCD 進行全面的眾核優化。目前版本所使用的格點數量尚且較小，后面將優化所采用的stencil 模型，增大數據量，進一步發掘LDM 與寄存器通訊機制的功能，以更加充分地利用從核陣列的并行計算能力，提高運行效率。

（2）通過進一步的分析，極大地消除MPI 通訊瓶頸阻礙，或將使用MPI-2 等其他版本進行開發，以期進一步擴大并行規模，充分挖掘神威平臺的整體計算能力。

（3）在并行規模和效率達到滿意的程度后，將程序與Chroma框架進行整合，便于科研人員使用，為高能物理領域做出些許貢獻。