孫晨 王彬 顧文靜 魏敏

（國(guó)家氣象信息中心，北京 100081）

0　引言

GRAPES（Global/Regional Assimilation and Prediction System）是中國(guó)氣象局自主研發(fā)的靜力/非靜力多尺度通用數(shù)值預(yù)報(bào)系統(tǒng)，該系統(tǒng)是氣象與氣候研究的基礎(chǔ)和核心。

GRAPES系統(tǒng)的核心部分包括模式動(dòng)力框架和物理過(guò)程，動(dòng)力框架中計(jì)算時(shí)間最長(zhǎng)的為半拉格朗日計(jì)算和求解Helmholtz大型線(xiàn)性代數(shù)方程組，物理過(guò)程中耗時(shí)較大的為輻射和微物理過(guò)程，其中輻射過(guò)程計(jì)算耗時(shí)占到GRAPES模式的40%，所以提升輻射過(guò)程的性能和計(jì)算效率，對(duì)整個(gè)GRAPES模式的性能提升具有重要科研價(jià)值和實(shí)際意義[1-2]。

近年來(lái)，并行計(jì)算技術(shù)日趨成熟，GPU巨大的計(jì)算能力已經(jīng)吸引了越來(lái)越多的科研人員將其作為一個(gè)高效益低成本的高性能計(jì)算平臺(tái)，然而，要在GPU硬件上實(shí)現(xiàn)加速需要通過(guò)對(duì)底層的API進(jìn)行編程來(lái)實(shí)現(xiàn)，程序編寫(xiě)復(fù)雜、難度大，且容易形成高度依賴(lài)特定設(shè)備的代碼。因此，一種基于指令制導(dǎo)計(jì)算的編程模OpenACC應(yīng)運(yùn)而生[3]，OpenACC的編程機(jī)制是在源程序中添加少量編譯標(biāo)識(shí)，編譯器根據(jù)作者的意圖自動(dòng)產(chǎn)生低級(jí)語(yǔ)言代碼，無(wú)需學(xué)習(xí)新的編程語(yǔ)言和加速器硬件知識(shí)，變更迅速掌握。只添加少量編譯標(biāo)識(shí)不破壞原代碼，開(kāi)發(fā)速度快，即可并行執(zhí)行又可恢復(fù)串行執(zhí)行；硬件更新時(shí)，重新編譯一次代碼即可，不必手工修改代碼。OpenACC標(biāo)準(zhǔn)制定時(shí)就考慮了目前及將來(lái)多種加速器產(chǎn)品，同一份代碼可以在多種加速器設(shè)備上編譯、運(yùn)行、無(wú)成本切換硬件平臺(tái)。

本文以GRAPES模式物理過(guò)程中的長(zhǎng)波輻射為加速處理對(duì)象，依托國(guó)家超級(jí)計(jì)算無(wú)錫中心的“神威·太湖之光”高性能計(jì)算機(jī)系統(tǒng)，設(shè)計(jì)一種基于OpenACC的加速算法，最大程度繼承原始代碼的MPI 級(jí)并行，同時(shí)利用數(shù)量眾多的協(xié)處理器加速熱點(diǎn)函數(shù)。

1　“神威·太湖之光”簡(jiǎn)介

2016年6月20日，新一期全球超級(jí)計(jì)算機(jī)500強(qiáng)榜單在德國(guó)法蘭克福舉辦的國(guó)際超算大會(huì)（ISC）上公布，“神威·太湖之光”超級(jí)計(jì)算機(jī)榮登榜首。

“神威·太湖之光”由國(guó)家并行計(jì)算機(jī)工程技術(shù)研究中心研制，是全球第一臺(tái)運(yùn)行速度超過(guò)10億億次/s的超級(jí)計(jì)算機(jī)，峰值性能高達(dá)12.54億億次/s，持續(xù)性能達(dá)到9.3億億次/s。該系統(tǒng)部署在國(guó)家超級(jí)計(jì)算無(wú)錫中心，系統(tǒng)由40個(gè)運(yùn)算機(jī)柜和8個(gè)網(wǎng)絡(luò)機(jī)柜組成。每個(gè)機(jī)柜有4個(gè)超級(jí)節(jié)點(diǎn)，每個(gè)超級(jí)節(jié)點(diǎn)包括32個(gè)節(jié)點(diǎn)板，每個(gè)節(jié)點(diǎn)板上有4個(gè)節(jié)點(diǎn)卡，每個(gè)節(jié)點(diǎn)卡有兩個(gè)節(jié)點(diǎn)，每個(gè)節(jié)點(diǎn)上裝有1個(gè)“申威26010”眾核處理器和32GB的DDR3內(nèi)存，一臺(tái)機(jī)柜就有1024塊處理器。節(jié)點(diǎn)之間通過(guò)基于PCI-E3.0（外設(shè)部件互聯(lián)標(biāo)準(zhǔn)擴(kuò)展3.0版）的神威高速網(wǎng)絡(luò)進(jìn)行互聯(lián)。

“神威·太湖之光”系統(tǒng)全部采用“申威26010”眾核處理器，架構(gòu)如圖1所示。

該芯片主頻為145 GHz，由4個(gè)核組組成，每個(gè)核組包含一個(gè)主核和一個(gè)從核簇，每個(gè)從核簇由64個(gè)從核組成，共260個(gè)處理器核心，同一個(gè)核組內(nèi)的主核和64個(gè)從核共享主存，每個(gè)從核有獨(dú)立的、容量為64 KB的高速緩存（Local Data Memory，LDM），支持DMA（Direct Memory Access）的方式在主存和LDM之間傳輸數(shù)據(jù)。主核作為處理器核心，可以進(jìn)行通信、IO、計(jì)算等操作，同一核組內(nèi)的64個(gè)從核作為加速計(jì)算部件，用來(lái)加速主核代碼中計(jì)算密集部分。

圖1　“申威26010”處理器架構(gòu)圖Fig. 1 Structure of the SW-26010 processor

2　OpenACC＊總述

OpenACC是由OpenACC組織（www.openaccstandard.org）于2011年推出的眾核加速編程語(yǔ)言，以編譯指示的方式提供眾核編程所需的語(yǔ)言功能，其主要目的是降低眾核編程的難度。用于C/C++和Fortran的OpenACC API和指令負(fù)責(zé)把底層的GPU任務(wù)交給編譯的同時(shí)，又提供跨系統(tǒng)、跨主機(jī)CPU和加速設(shè)備的支持[3]。

“神威·太湖之光”計(jì)算機(jī)系統(tǒng)中的OpenACC*語(yǔ)言是在OpenACC2.0文本的基礎(chǔ)上，針對(duì)申威26010眾核處理器結(jié)構(gòu)特點(diǎn)進(jìn)行適當(dāng)?shù)木?jiǎn)和擴(kuò)充而來(lái)的。

2.1　執(zhí)行模型

OpenACC*程序的執(zhí)行模型是在host（主處理器）的指導(dǎo)下，host和device（加速設(shè)備）協(xié)作的加速執(zhí)行模型，如圖2所示。

圖2　OpenACC＊執(zhí)行模型Fig. 2 Execution model of the OpenACC＊

程序首先在host上啟動(dòng)運(yùn)行，以一個(gè)主線(xiàn)程串行執(zhí)行，或者通過(guò)使用OpenMP或MPI等編程接口以多個(gè)主線(xiàn)程并行執(zhí)行，計(jì)算密集區(qū)域則在主線(xiàn)程的控制下作為kernel或parallel（加速構(gòu)件）被加載到加速設(shè)備上執(zhí)行。加速構(gòu)件的執(zhí)行過(guò)程包括：在設(shè)備內(nèi)存上分配所需數(shù)據(jù)空間；加載構(gòu)件代碼（包含相關(guān)參數(shù)）至device；構(gòu)件將所需數(shù)據(jù)從主存?zhèn)鬏斨猎O(shè)備內(nèi)存；等待數(shù)據(jù)傳輸完成；device進(jìn)行計(jì)算并將結(jié)果傳送回主存；釋放設(shè)備上的數(shù)據(jù)空間等。大多數(shù)情況下，host可以加載一系列構(gòu)件，并在加速設(shè)備上逐個(gè)執(zhí)行[4]。

OpenACC*支持三級(jí)并行機(jī)制：gang、worker、vector。gang是粗粒度并行，在加速設(shè)備上可以啟動(dòng)一定數(shù)量的gang。worker 是細(xì)粒度并行，每個(gè)gang內(nèi)包含有一定數(shù)量的worker。vector并行是在worker內(nèi)通過(guò)SIMD或向量操作的指令級(jí)并行。

在“申威26010”中，一個(gè)運(yùn)算控制核心（簡(jiǎn)稱(chēng)主核）僅控制一個(gè)運(yùn)算核心陣列(加速設(shè)備，簡(jiǎn)稱(chēng)從核陣列)的運(yùn)行，每個(gè)運(yùn)算核心陣列內(nèi)有64個(gè)運(yùn)算核心（簡(jiǎn)稱(chēng)從核），每個(gè)運(yùn)算核心可以運(yùn)行一個(gè)加速線(xiàn)程。默認(rèn)情況下，64個(gè)加速線(xiàn)程被組織成64個(gè)gang、每個(gè)gang內(nèi)一個(gè)worker、worker內(nèi)可以vector并行的邏輯視圖。通常情況下，盡量用滿(mǎn)64個(gè)加速線(xiàn)程可以獲得較好的運(yùn)行性能。

2.2　存儲(chǔ)模型

在CUDA和OpenCL等低層次加速器編程語(yǔ)言中，主機(jī)和加速器存儲(chǔ)器分離的概念非常明確，內(nèi)存間移動(dòng)數(shù)據(jù)的語(yǔ)句占據(jù)用戶(hù)大部分代碼。

“申威26010”中從核和主核共享內(nèi)存，從核可直接訪(fǎng)問(wèn)主存空間，并在從核內(nèi)提供加速線(xiàn)程私有的LDM，加速計(jì)算需要存放到LDM的數(shù)據(jù)由從核控制傳輸。存儲(chǔ)模型如圖3所示。

圖3　OpenACC＊存儲(chǔ)模型Fig. 3 Memory model of the OpenACC＊

OpenACC標(biāo)準(zhǔn)中的數(shù)據(jù)管理功能將不產(chǎn)生實(shí)際的作用，但直接訪(fǎng)問(wèn)主存往往帶來(lái)性能損失，因此需要充分利用從核中的LDM，提升數(shù)據(jù)訪(fǎng)問(wèn)效率。在“神威·太湖之光”中，OpenACC*對(duì)OpenACC標(biāo)準(zhǔn)所做的主要功能延伸和語(yǔ)法擴(kuò)展就是為了解決共享內(nèi)存架構(gòu)下片內(nèi)高速存儲(chǔ)空間的使用問(wèn)題。

“申威26010”中，從一個(gè)加速線(xiàn)程的角度來(lái)看，其可見(jiàn)的數(shù)據(jù)空間有三種[4]：

1）主線(xiàn)程數(shù)據(jù)空間：位于主存，主線(xiàn)程的內(nèi)存空間對(duì)其創(chuàng)建的加速線(xiàn)程直接可見(jiàn)，加速線(xiàn)程可以直接訪(fǎng)問(wèn)相應(yīng)的數(shù)據(jù)；對(duì)于主線(xiàn)程創(chuàng)建的多個(gè)加速線(xiàn)程而言，這部分空間是共享的；程序中在加速區(qū)外定義的變量，均位于該空間內(nèi)。

2）加速線(xiàn)程私有空間：位于主存，每個(gè)加速線(xiàn)程有獨(dú)立的私有空間，程序中使用private、firstprivate子句修飾的變量將存放于該空間內(nèi)。

3）加速線(xiàn)程本地空間：位于設(shè)備內(nèi)存，每個(gè)加速線(xiàn)程有獨(dú)立的本地空間LDM，本地空間的訪(fǎng)問(wèn)性能是三種空間中最高的。程序中使用local、copy等數(shù)據(jù)子句修飾的變量將由編譯系統(tǒng)控制全部或局部存放于LDM內(nèi)。主存與本地空間的數(shù)據(jù)交互由加速線(xiàn)程控制。

3　長(zhǎng)波輻射過(guò)程O(píng)penACC優(yōu)化方案

GRAPES模式的輻射模型采用歐洲中期天氣預(yù)報(bào)中心（ECMWF）的長(zhǎng)短波輻射方案，該方案將整個(gè)大氣層劃分成水平方向和垂直方向的三維網(wǎng)格，水平方向?yàn)闄M跨地球表面的二維網(wǎng)格，垂直方向則表示大氣的層數(shù)[5]。數(shù)值天氣預(yù)報(bào)模式是一種非線(xiàn)性化離散計(jì)算模式，其計(jì)算量巨大[6]，最初的都是通過(guò)CPU計(jì)算實(shí)現(xiàn)，GPU出現(xiàn)之后，雖然可以把該算法移植到GPU上去并行執(zhí)行，但編程難度大且易出錯(cuò)。因此，文中設(shè)計(jì)一種基于OPENACC加速方法，以實(shí)現(xiàn)通過(guò)簡(jiǎn)化的編程算法和簡(jiǎn)潔的Fortran語(yǔ)言代碼完成對(duì)長(zhǎng)波輻射過(guò)程并行處理過(guò)程。

3.1　輻射過(guò)程的構(gòu)成

輻射過(guò)程是GRAPES系統(tǒng)中一個(gè)重要的物理過(guò)程，包含云結(jié)構(gòu)的描述和輻射算法（RRTM模塊）兩部分。其中云對(duì)調(diào)解輻射平衡起著至關(guān)重要的作用[7]，本系統(tǒng)在云產(chǎn)生器的基礎(chǔ)上利用McICA進(jìn)行輻射計(jì)算。輻射過(guò)程首先調(diào)用mcica_subcol_lw函數(shù)實(shí)現(xiàn)云結(jié)構(gòu)的描述，之后調(diào)用RRTM模塊進(jìn)行輻射計(jì)算。RRTM模塊初始化函數(shù)rrtmg_lw_init讀取數(shù)據(jù)文件中的輸入數(shù)據(jù)；然后調(diào)用rrtmg_lw_rad函數(shù)計(jì)算輻射傳輸?shù)倪^(guò)程；最后調(diào)用輻射模型子函數(shù)rrtm_lw；rrtm_lw子函數(shù)一次計(jì)算一個(gè)垂直列的輻射傳輸值[6]。它主要包含以下幾個(gè)基本步驟：inatm，準(zhǔn)備大氣剖面；cldprmc，計(jì)算云層光學(xué)深度；setcoef，計(jì)算這種大氣剖面下輻射傳輸算法所需的各種量；taumol，計(jì)算氣體光學(xué)深度和普朗克分?jǐn)?shù)；rtrnmc，計(jì)算晴空和云層的輻射傳輸值。其中，taumol中的16個(gè)子函數(shù)taugb，計(jì)算了全部16個(gè)長(zhǎng)波波域譜帶的氣態(tài)光學(xué)厚度和普朗克（Planck）函數(shù)。輻射過(guò)程的函數(shù)包結(jié)構(gòu)如圖4所示。

圖4　長(zhǎng)波輻射過(guò)程的函數(shù)結(jié)構(gòu)Fig. 4 Structure diagram of the long-wave radiation function

3.2　并行優(yōu)化方案

OPENACC并行化的唯一目的是利用從核充足的硬件資源來(lái)提高程序運(yùn)行速度，縮短運(yùn)行時(shí)間。整個(gè)長(zhǎng)波輻射過(guò)程中最耗時(shí)間的是RRTMG_LW主過(guò)程，計(jì)算并行化的目標(biāo)是將循環(huán)迭代步分散到多個(gè)不同的線(xiàn)程上執(zhí)行，這些線(xiàn)程運(yùn)行在多個(gè)加速器核心上，從而將計(jì)算任務(wù)由CPU轉(zhuǎn)移到加速器上，減輕CPU的負(fù)擔(dān)。

結(jié)合程序代碼的分析，確定需要眾核加速的代碼段，并根據(jù)源代碼數(shù)據(jù)結(jié)構(gòu)進(jìn)行預(yù)處理，對(duì)相關(guān)數(shù)組和循環(huán)進(jìn)行調(diào)整，以適合OpenACC加速。優(yōu)化的基本思想是確定將程序中的標(biāo)量和關(guān)鍵數(shù)組盡可能多的放到局存LDM中，并設(shè)法提高程序中數(shù)據(jù)的傳輸效率，將需要使用眾核加速的循環(huán)的前后使用OpenACC加速編譯指示進(jìn)行標(biāo)注，其中PARALLEL表明該部分代碼是需要加速執(zhí)行的并行代碼；LOOP表示下面緊跟著的循環(huán)需要在加速線(xiàn)程間進(jìn)行并行劃分；DATA表示主從核之間的數(shù)據(jù)拷貝。

rrtm_lw函數(shù)包含一層循環(huán)，涵蓋inatm、cldprmc、setcoef、taumol和rtrnmc五個(gè)子函數(shù)調(diào)用，每個(gè)子函數(shù)內(nèi)部有多個(gè)do循環(huán)分塊，最外層的iplon循環(huán)閾值為[1，90]，根據(jù)OPENACC以及申威“26010”處理器硬件特點(diǎn)，可以采取的并行化方案有如下幾種：

3.2.1從核并行化方案1：最外層直接并行

結(jié)合輻射過(guò)程的計(jì)算模型為單柱計(jì)算模型，核心迭代針對(duì)每個(gè)單柱串行進(jìn)行多個(gè)物理過(guò)程的計(jì)算，每個(gè)單柱之間不存在數(shù)據(jù)依賴(lài)關(guān)系?？梢栽谠搶友h(huán)外層添加OpenACC加速指示語(yǔ)句，單個(gè)從核將獨(dú)立運(yùn)行循環(huán)內(nèi)部的所有函數(shù)，完成全部計(jì)算后將結(jié)果傳回主核，代碼結(jié)構(gòu)示意如圖5所示，此時(shí)從核使用率為70.3%。

圖5　iplon層直接并行代碼結(jié)構(gòu)示意圖Fig. 5 Structure of direct parallel methods on the iplon layer

本段核心代碼中，5個(gè)子函數(shù)中包含眾多變量，多為二維數(shù)組，少數(shù)為三維數(shù)組（標(biāo)量及一維數(shù)組產(chǎn)生局存壓力較?。Ｆ渲休斎霐?shù)組22個(gè)，包含溫度，壓力，痕量氣體濃度等多個(gè)物理量，總大小約為22×140×60(or 61)×8=1478400 byte；計(jì)算過(guò)程中的中間變量超過(guò)40個(gè)，總大小約為2688000 byte；輸出數(shù)組多為一維數(shù)組或標(biāo)量。即使通過(guò)調(diào)整代碼使計(jì)算過(guò)程中的大多數(shù)中間變量空間能及時(shí)釋放，棧容量需求也遠(yuǎn)超局存大?。?5536 byte），以0.5°分辨率計(jì)算過(guò)程為例，無(wú)論是添加routine函數(shù)將部分從核函數(shù)的棧數(shù)組放到主存，還是使用分塊語(yǔ)句指導(dǎo)數(shù)據(jù)依次拷入拷出，都會(huì)產(chǎn)生極大的通訊開(kāi)銷(xiāo)，影響運(yùn)行速度。

另外，子函數(shù)cldprmc與setcoef只有簡(jiǎn)單的賦值與計(jì)算過(guò)程，在原本的運(yùn)行過(guò)程中耗時(shí)很少，但與前后子函數(shù)存在數(shù)據(jù)依賴(lài)，所以此方案下，同樣需要拷貝程序與數(shù)據(jù)在從核運(yùn)行，運(yùn)行時(shí)間大量增加。

3.2.2從核并行化方案2：循環(huán)內(nèi)移

對(duì)于氣象模式，運(yùn)算過(guò)程中的中間變量數(shù)組通常個(gè)數(shù)很多并且數(shù)據(jù)量極大，由于局存壓力需要多次多步拷貝，數(shù)據(jù)通信會(huì)造成較大的時(shí)間開(kāi)銷(xiāo)影響加速效果，例如長(zhǎng)波輻射過(guò)程的rrtmg_lw函數(shù)。

此時(shí)可以通過(guò)循環(huán)內(nèi)移，對(duì)大函數(shù)進(jìn)行分塊，內(nèi)移工作從rrtmg_lw函數(shù)調(diào)用的函數(shù)inatm開(kāi)始，依次對(duì)cldprmc、setcoef、taumol和rtrnmc子函數(shù)進(jìn)行循環(huán)內(nèi)移達(dá)到子函數(shù)分塊的目的，之后，對(duì)于適合從核并行的inatm、taumol和rtrnmc子函數(shù)按OPENACC規(guī)則添加導(dǎo)語(yǔ)及子句使其在從核運(yùn)行，不適合從核運(yùn)行的cldprmc、setcoef子函數(shù)仍在主核運(yùn)行，代碼結(jié)構(gòu)示意如圖6所示。此時(shí)，數(shù)據(jù)拷貝通訊量縮短至方案1的20%左右，與此同時(shí)，可根據(jù)每個(gè)子函數(shù)的最外層循環(huán)數(shù)，合并外層循環(huán)，增加從核運(yùn)行指示語(yǔ)句，最大化擴(kuò)大子函數(shù)運(yùn)行的并行度。

圖6　循環(huán)內(nèi)移代碼結(jié)構(gòu)示意圖Fig. 6 Diagram of the loop inside displacement

循環(huán)內(nèi)移使程序原本的計(jì)算順序發(fā)生了改變，由從核單線(xiàn)程串行計(jì)算子函數(shù)1-5變?yōu)橄扔?jì)算所有子函數(shù)1，保存計(jì)算結(jié)果后計(jì)算子函數(shù)2，以此類(lèi)推。所以，代碼中大部分中間變量以及輸出結(jié)果需要升維用于解決數(shù)據(jù)依賴(lài)，保證正確性。

循環(huán)內(nèi)移后，局存壓力顯著降低，特別是rtrnmc函數(shù)，只有6個(gè)二維輸入變量和2個(gè)二維中間變量，經(jīng)過(guò)計(jì)算和試驗(yàn)，只需最低12個(gè)循環(huán)分塊即可正確運(yùn)行，通訊次數(shù)降低，運(yùn)算速度顯著提高。同時(shí)，inatm函數(shù)作為大氣剖面準(zhǔn)備模塊，輸入數(shù)據(jù)較多，中間變量較少，循環(huán)內(nèi)移后，從核運(yùn)行時(shí)所需的循環(huán)分塊個(gè)數(shù)以及數(shù)據(jù)拷貝次數(shù)也明顯降低，運(yùn)行速度有所提高。

3.2.3從核并行化方案3：循環(huán)分列

為最大程度利用異構(gòu)眾核處理器的結(jié)構(gòu)優(yōu)勢(shì)，同時(shí)進(jìn)一步降低局存壓力，減少通訊時(shí)間，可以對(duì)循環(huán)內(nèi)移后的子函數(shù)進(jìn)行循環(huán)改寫(xiě)，擴(kuò)大并行度使每個(gè)從核承擔(dān)更少的計(jì)算任務(wù)，從而降低計(jì)算所需通訊量。通過(guò)分析子函數(shù)循環(huán)結(jié)構(gòu)，對(duì)于迭代層數(shù)較多，耗時(shí)較長(zhǎng)的計(jì)算熱點(diǎn)進(jìn)行進(jìn)一步從核并行化改寫(xiě)。同時(shí)，通過(guò)collapse循環(huán)合并語(yǔ)句，可以使代碼運(yùn)行時(shí)的并行度顯著增加。以inatm函數(shù)為例，改寫(xiě)循環(huán)后k由子函數(shù)內(nèi)層轉(zhuǎn)至函數(shù)外，為保證數(shù)據(jù)正確性，需改寫(xiě)原函數(shù)，部分常量數(shù)據(jù)需要升維以符合最外層k循環(huán)，部分變量數(shù)據(jù)需要轉(zhuǎn)置以適合最高效的訪(fǎng)存方式，新函數(shù)名稱(chēng)以inatm_trans表示，循環(huán)分列前與分列后的代碼結(jié)構(gòu)示意如圖7所示，此時(shí)從核使用率為99.3%。

圖7　inatm函數(shù)循環(huán)分列前后結(jié)構(gòu)示意圖Fig. 7 Diagrams before and behind the loop splitting

循環(huán)分列的目的是通過(guò)減少單個(gè)從核承擔(dān)計(jì)算任務(wù)量的方式降低數(shù)據(jù)通信量，同時(shí)增大并行度且簡(jiǎn)化程序結(jié)構(gòu)，需要在一定程度上改寫(xiě)代碼結(jié)構(gòu)使其重組，合并，擴(kuò)大最外層循環(huán)大小。整個(gè)rrtmg_lw函數(shù)改寫(xiě)合并循環(huán)后，同時(shí)并行的線(xiàn)程增大至5400個(gè)，每個(gè)從核的輸入數(shù)據(jù)量為22×140×8=24640 byte，中間變量大小為7×140×8=7840 byte，總數(shù)據(jù)量已經(jīng)低于局存大小，不再需要數(shù)據(jù)分塊拷貝，重新調(diào)整其通訊策略后，運(yùn)行速度得到進(jìn)一步提升。

對(duì)于某些原本不適合進(jìn)行從核并行的計(jì)算核心段需要循環(huán)改寫(xiě)后再進(jìn)行從核并行，否則只能通過(guò)循環(huán)內(nèi)移分塊跳過(guò)，使其在主核完成計(jì)算過(guò)程，影響模式運(yùn)行效率。例如計(jì)算氣態(tài)光學(xué)厚度和普朗克量的taumol函數(shù)包含16個(gè)子函數(shù)對(duì)應(yīng)lay變量的16種分割方式，不具備從核并行條件，通過(guò)合并成為一個(gè)最外層以lay為變量的大循環(huán)后，順利在從核運(yùn)行且使用率超過(guò)90%。

在實(shí)際情況下，方案1，2，3實(shí)現(xiàn)數(shù)據(jù)并行計(jì)算的不同思路通常需要結(jié)合使用。針對(duì)不同核心函數(shù)不同的計(jì)算結(jié)構(gòu)與數(shù)據(jù)規(guī)模，分析設(shè)計(jì)出各自最高效的通訊策略。本次工作中主要選取了長(zhǎng)波輻射部分的2個(gè)主函數(shù)和7個(gè)子函數(shù)作為核心段，最終實(shí)現(xiàn) “粗粒度MPI并行+細(xì)粒度眾核并行”多級(jí)異構(gòu)并行方案。

3.2.4其他并行優(yōu)化措施

使用ldm數(shù)學(xué)函數(shù)庫(kù)：在OpenACC加速指導(dǎo)語(yǔ)句的基礎(chǔ)上，使用ldm數(shù)學(xué)函數(shù)庫(kù)，經(jīng)測(cè)試，ldm數(shù)學(xué)函數(shù)精度與默認(rèn)庫(kù)保持一致，計(jì)算性能提升30%～80%不等，進(jìn)一步提升了GRAPES系統(tǒng)整體運(yùn)算效率。

3.2.5結(jié)果與分析

試驗(yàn)采用分辨率為0.5°×0.5°的GRAPES全球模式長(zhǎng)波輻射方案，積分步長(zhǎng)為36，優(yōu)化選項(xiàng)使用-O2。由于cldprmc、setcoef耗時(shí)較少，沒(méi)有進(jìn)行優(yōu)化。加速后，對(duì)比純主核運(yùn)行輸出的全球長(zhǎng)波輻射通量與模式輸出值，經(jīng)驗(yàn)證誤差，相對(duì)誤差的量級(jí)在10－5至l0－5之間，在允許范圍內(nèi)。優(yōu)化結(jié)果對(duì)比如表1所示。

表1　眾核加速結(jié)果Table 1 The result of many-core acceleration

4　展望與不足

隨著編譯器的進(jìn)一步優(yōu)化和硬件技術(shù)的發(fā)展，OpenACC加速與CUDA等在底層技術(shù)實(shí)現(xiàn)上的差距將會(huì)越來(lái)越小，而且支持OpenACC加速指令轉(zhuǎn)換將不僅僅針對(duì)CUDA設(shè)備，還包括其他更多廠(chǎng)商的硬件加速設(shè)備，從而能極大地提高OpenACC加速指令的普適性，增大程序可移植性。

文中以GRAPES全球模式RRTMG輻射模塊長(zhǎng)波過(guò)程為加速對(duì)象，目前只應(yīng)用到為數(shù)不多的熱點(diǎn)函數(shù)，就能得到較好的加速效果，優(yōu)化還有很大的提升空間。

下一步將繼續(xù)深入研究源代碼數(shù)據(jù)結(jié)構(gòu)，加強(qiáng)數(shù)據(jù)預(yù)處理，減少?gòu)暮伺c主存之間的通訊頻率，充分利用從核中高速緩存LDM，提升數(shù)據(jù)訪(fǎng)問(wèn)效率；同時(shí)嘗試多種OpenACC加速指令組和使用，進(jìn)一步提升加速效果。

[1]伍湘君, 金之雁, 黃麗萍, 等．GRAPES模式軟件框架與實(shí)現(xiàn). 應(yīng)用氣象學(xué)報(bào), 2005, 16(4): 539-546．

[2]陳德輝, 沈?qū)W順．新一代數(shù)值預(yù)報(bào)系統(tǒng)GRAPES研究進(jìn)展. 應(yīng)用氣象學(xué)報(bào), 2006, 17(6): 773-777．

[3]曾文權(quán), 胡玉貴, 何擁軍, 等．一種基于OPENACC的GPU加速實(shí)現(xiàn)高斯模糊算法. 計(jì)算機(jī)技術(shù)與發(fā)展, 2013, 23(7): 147-150．

[4]何滄平. OpenACC并行編程實(shí)戰(zhàn). 北京: 機(jī)械工業(yè)出版社, 2016．

[5]郭妙, 金之雁, 周斌. 基于通用圖形處理器的GRAPES長(zhǎng)波輻射并行方案. 應(yīng)用氣象學(xué)報(bào), 2012, 23(3): 348-354．

[6]鄭芳, 許先斌, 向冬冬, 等. 基于GPU的GRAPES數(shù)值預(yù)報(bào)系統(tǒng)中RRTM模塊的并行化研究. 計(jì)算機(jī)科學(xué), 2012, 39(6A): 370-374.

[7]荊現(xiàn)文, 張華. McICA 云—輻射方案在國(guó)家氣候中心全球氣候模式中的應(yīng)用與評(píng)估. 大氣科學(xué), 2012, 36(5): 945-958.