基于CUDA的并行雷達拼圖算法研究*

韓 豐 高 嵩 薛 峰 李月安

國家氣象中心,北京 100081

提 要：雷達組網(wǎng)拼圖算法是強對流天氣短時臨近預(yù)報系統(tǒng)(Severe Weather Automatic Nowcasting,SWAN)的重要基礎(chǔ)方法之一。提高拼圖算法的效率,不僅可以提升現(xiàn)有SWAN臨近算法序列的時效性,也能更好地應(yīng)用高分辨率雷達數(shù)據(jù),具有重要的實際意義。采用中央處理器(central processing unit,CPU)和圖形處理器(graphics processing unit,GPU)混合架構(gòu)設(shè)計并行雷達拼圖算法,其中CPU負(fù)責(zé)雷達數(shù)據(jù)的解析和調(diào)度GPU并行模塊,GPU負(fù)責(zé)大規(guī)模數(shù)據(jù)的并行計算。通過分析計算統(tǒng)一設(shè)備架構(gòu)(compute unified device architecture,CUDA)算法的并行開銷和拼圖算法的特點,提出并實現(xiàn)了GPU內(nèi)存管理優(yōu)化和數(shù)據(jù)交換流程優(yōu)化方案,提高了組網(wǎng)拼圖算法的效率。對比試驗結(jié)果表明,基于CUDA的GPU并行拼圖算法和SWAN中30線程并行的CPU算法相比,在全國1 km和500 m分辨率的拼圖任務(wù)上,加速比分別達到3.52和6.82。綜上,基于CUDA的并行拼圖算法不僅可以提高SWAN短時臨近算法序列的時效性,也為更高分辨率雷達資料的拼圖提供了技術(shù)支持。

引 言

隨著觀測技術(shù)和預(yù)報業(yè)務(wù)需求的發(fā)展,氣象雷達的分辨率不斷提高。2019年,中國氣象局啟動我國新一代天氣雷達的雙偏振改造和技術(shù)升級,升級后的S波段雷達強度徑向分辨率由1 km提升至250 m。同時,為了滿足精細(xì)化臨近預(yù)報業(yè)務(wù)的需求,北京市、江蘇省、廣東省等地區(qū)先后建立X波段雙偏振雷達網(wǎng)并投入業(yè)務(wù)使用(張哲等,2021;蘇永彥和劉黎平等,2022;Ryzhkov et al,2005)。雷達徑向分辨率的提升對數(shù)據(jù)的快速處理和客觀算法提出了更高的效率要求。

雷達數(shù)據(jù)的三維組網(wǎng)拼圖技術(shù)是雷達應(yīng)用的重要基礎(chǔ)方法之一(程元慧等,2020;張勇等,2019)。在進行雷達組網(wǎng)時,首先需要將回波從球坐標(biāo)轉(zhuǎn)換到等經(jīng)緯度的直角坐標(biāo)系下,再通過最大值、均值、距離反比等方案,在多部雷達的觀測重疊區(qū)進行融合計算(肖艷姣和劉黎平,2006)。拼圖計算過程中,涉及插值、坐標(biāo)系轉(zhuǎn)換、大規(guī)模數(shù)據(jù)網(wǎng)格等問題,計算密集且數(shù)據(jù)吞吐量大,是對計算資源需求量最大的雷達客觀方法之一。以強對流短時臨近預(yù)報系統(tǒng)(Server Weather Automatic Nowcasting,SWAN)(韓豐和沃偉峰,2018)為例,全國1 km雷達組網(wǎng)拼圖算法是臨近預(yù)報算法中耗時最長的算法,且雷達拼圖后續(xù)還有風(fēng)暴識別和追蹤算法(Thunderstorm Identification,Tracking,Analysis and Nowcasting,TITAN)、外推等其他客觀方法。所以提升雷達拼圖的運算速度可以明顯提升SWAN臨近預(yù)報產(chǎn)品序列的時效性。此外,考慮到雷達徑向分辨率的提升,將拼圖的分辨率提升至500 m時,現(xiàn)有的拼圖算法無法在一個觀測周期(6 min)內(nèi)完成計算,不具備業(yè)務(wù)應(yīng)用能力。所以將更先進的并行計算技術(shù)應(yīng)用于雷達組網(wǎng)拼圖算法中,一方面可以提升現(xiàn)有算法序列的時效性,另一方面也可以更好地應(yīng)用高分辨率雷達基數(shù)據(jù),具有重要的實際意義。

氣象部門常用中央處理器(central porcessing unit,CPU)并行技術(shù)和基于CPU的分布式計算框架來提高運算能力。王志斌等(2009)通過使用共享存儲并行編程(open multiple processing,OpenMP)并行技術(shù)和消息傳遞接口(message passing interface,MPI)編程規(guī)范,實現(xiàn)雷達拼圖系統(tǒng)的兩級并行。吳石磊等(2012)使用多線程技術(shù)對全國組網(wǎng)雷達估測降水系統(tǒng)進行效率優(yōu)化,通過調(diào)整計算負(fù)載,可以將計算效率提升5倍。SWAN系統(tǒng)設(shè)計實現(xiàn)了基于OpenMP的雷達組網(wǎng)拼圖算法,使用30個邏輯線程,生成全國1 km分辨率三維拼圖產(chǎn)品需約 2 min。

經(jīng)過多年的發(fā)展,CPU性能提升已經(jīng)逐漸達到瓶頸,多核應(yīng)用的成本-效益比不斷下降,僅僅通過挖掘CPU性能來提升算法效率已經(jīng)不能滿足現(xiàn)有的需求。圖形處理器(graphics processing unit,GPU)作為一種新的計算資源,由于其區(qū)別于CPU的硬件架構(gòu),特別適合解決大規(guī)模數(shù)據(jù)的計算問題。同時GPU并行計算接口的快速發(fā)展,也使得GPU在通用計算領(lǐng)域的應(yīng)用越來越廣泛(錢奇峰等,2021;黃驕文等,2021)。

目前國內(nèi)外也有一些使用GPU加速氣象算法的研究成果。Vu et al(2013)、Mielikainen et al(2013)、Govett et al(2017)、Sch?r et al(2020)、王藝儒(2020)、Wang et al(2021)和肖灑等(2019)等學(xué)者將GPU并行技術(shù)應(yīng)用于數(shù)值模式,改寫相關(guān)動力框架和計算模塊,在保持計算精度的同時,顯著提高模式的運行效率。Olander et al(2021)介紹了GPU在臺風(fēng)定強深度學(xué)習(xí)模型中的應(yīng)用。Rukmono et al(2018)、王興等(2015)和胡榮華(2015)將GPU應(yīng)用于雷達數(shù)據(jù)處理和臨近預(yù)報算法,得到了10倍以上的加速比。

由于GPU程序主要使用圖形應(yīng)用程序接口(application programming interface,API)進行編程,如氣象信息綜合分析處理系統(tǒng)(Meteorological Information Combine Analysis and Process System)MICAPS4(高嵩等,2017)使用開放式圖形語言(open graphics language,OpenGL)開發(fā)氣象數(shù)據(jù)渲染引擎。但是圖形API并不適用于通用計算任務(wù)。因此,比較主流的GPU并行算法一般通過一些標(biāo)準(zhǔn)化的并行編程接口實現(xiàn)。常用的并行計算接口有開放運算語言(open computing language,OpenCL)如王興等(2015),開源加速接口(open acceleration,OpenACC)如Govett et al(2017)、Sch?r et al(2020)、肖灑等(2019),計算統(tǒng)一設(shè)備架構(gòu)(compute unified device architecture,CUDA),如Wang et al(2021)等。本文設(shè)計基于CUDA的雷達組網(wǎng)拼圖算法,通過分析CUDA算法的并行開銷和組網(wǎng)拼圖算法的流程,提出高效的CPU+GPU的并行拼圖算法實現(xiàn)方案。通過和現(xiàn)有SWAN中的CPU并行組網(wǎng)拼圖算法對比,獲得較好的性能優(yōu)化和實際效果。

綜上所述,本文為了提高SWAN臨近客觀算法的時效性,更好地應(yīng)用高分辨率雷達資料,將GPU并行技術(shù)應(yīng)用于雷達組網(wǎng)拼圖算法,設(shè)計了基于CUDA的并行拼圖算法。通過算法性能分析,提出并實現(xiàn)了一種基于GPU內(nèi)存管理和數(shù)據(jù)流程優(yōu)化的方案。對比試驗表明,基于CUDA的并行算法能獲得超過3.5以上的加速比。

1 基于CUDA的GPU并行計算模型

1.1 CUDA編程模型

CUDA是NIVIDA公司發(fā)布的通用并行計算平臺和編程模型,它提供了GPU編程的簡易接口,如初始化GPU空間,主機和GPU設(shè)備間的數(shù)據(jù)傳輸?shù)??；贑UDA編程可以構(gòu)建GPU計算的應(yīng)用程序,而無需關(guān)注底層的圖形API。

CUDA是一個異構(gòu)編程模型,可以協(xié)同調(diào)度CPU和GPU,共同完成一個計算任務(wù)。在CUDA編程模型中,計算資源被分為主機端(host)和設(shè)備端(device),其中主機端表示CPU及其內(nèi)存,設(shè)備端表示GPU及其內(nèi)存。CUDA程序中既包含主機端模塊,又包含設(shè)備端模塊,它們分別在CPU和GPU上運行。

1.2 CUDA并行算法流程

如前文所述,GPU具有獨立的內(nèi)存系統(tǒng),所以CPU在調(diào)度GPU進行并行計算時,首先需要將數(shù)據(jù)從CPU內(nèi)存拷貝到GPU內(nèi)存,最后再將計算結(jié)果拷貝回CPU內(nèi)存。表1是一個典型CUDA算法和CPU算法的執(zhí)行流程。對比兩者可見,CUDA算法多了一次內(nèi)存的分配,一次內(nèi)存的釋放和兩次數(shù)據(jù)的搬運?？蓪⑦@四次操作認(rèn)為是CUDA并行算法所必須的并行開銷。下文將通過組網(wǎng)拼圖試驗,按照表1中的步驟,分析基于CUDA的雷達組網(wǎng)拼圖算法效率瓶頸,并給出優(yōu)化方案。

2 基于CUDA的雷達組網(wǎng)拼圖算法設(shè)計

2.1 基于CPU并行的雷達組網(wǎng)拼圖算法設(shè)計

在雷達組網(wǎng)拼圖算法中,首先需要每個單站分別進行三維離散化處理,最后將所有站的三維離散化結(jié)果拼接到一個大的數(shù)據(jù)網(wǎng)格中。在多雷達重疊部分,按照雷達融合策略進行取值。目前常用的策略為最大值法、距離反比權(quán)重法等。在SWAN中的CPU并行拼圖算法和本文設(shè)計的CUDA并行算法中,均采用肖艷姣和劉黎平(2006)提出的三維離散化方案,并按照距離反比權(quán)重法進行重疊區(qū)融合。

SWAN中目前使用的CPU并行拼圖算法包括兩個并行塊。并行塊1:由于每個站雷達數(shù)據(jù)都是獨立離散化的,則設(shè)計并行方案,在多線程中并行處理不同站的數(shù)據(jù);并行塊2:在離散化后數(shù)據(jù)組網(wǎng)拼圖過程中,由于不同垂直層之間也是相互獨立的,則設(shè)計并行方案,并行處理多層的數(shù)據(jù)拼接。并行塊1和并行塊2為順序執(zhí)行,即先完成全部站點的離散化之后,再進行多站拼圖。

雷達三維離散化插值分為兩步,首先計算出平面直角坐標(biāo)系中點P(x,y,z)在以雷達站為中心的球坐標(biāo)系中的位置(r,θ,φ),隨后將雷達基數(shù)據(jù)用雙線性插值的方法插值到直角坐標(biāo)系中?？梢圆捎妙A(yù)先計算坐標(biāo)轉(zhuǎn)換查找表的方式,加快雷達數(shù)據(jù)空間插值的計算速度。當(dāng)硬盤讀取速度大于查找表計算速度時,優(yōu)化是有效果的。

2.2 基于CUDA的雷達組網(wǎng)拼圖算法設(shè)計

2.2.1 CUDA并行粒度分析

上述CPU并行拼圖算法是合理而且有效的,主要是因為CPU一般只有幾十個邏輯核心,且每個核心的計算能力很強,可以處理復(fù)雜的任務(wù)。將不同的站點和層次分配到單獨的線程里面,既可以充分利用所有的核心,還能夠充分發(fā)揮CPU核心的性能,減少頻繁創(chuàng)建、銷毀線程的開銷。

但是GPU一般有幾十個流式多處理器,每個流式多處理器可以同時啟動上千個CUDA計算單元,且每一個計算單元的計算能力較弱,不能完成一個站點的三維離散化任務(wù)。所以要充分發(fā)揮GPU的并行性能,需要設(shè)計格點級的并行方案。

2.2.2 CUDA并行拼圖算法初步設(shè)計

圖1給出了CUDA并行拼圖算法初步設(shè)計流程圖。算法采用CPU+GPU并行混合的設(shè)計方案。其中,CPU負(fù)責(zé)雷達基數(shù)據(jù)的解析和最終拼圖產(chǎn)品的輸出,GPU負(fù)責(zé)每一步插值算法的并行計算。整個計算流程分為兩個部分,單站的三維離散化插值和多站拼圖。在單站三維離散化插值時,由CPU作為主要的調(diào)度程序,負(fù)責(zé)解壓、解析雷達基數(shù)據(jù),并調(diào)用GPU中的三維插值模塊,實現(xiàn)并行計算。不同雷達站之間的離散化部分采用串行的方式依次進行。在多站拼接階段,依然使用CPU作為調(diào)度程序,調(diào)用GPU中的空間插值模塊,實現(xiàn)多站拼圖的并行計算。

圖1 CUDA并行組網(wǎng)拼圖算法流程設(shè)計Fig.1 Process design of radar mosaic algorithm based on CUDA

2.2.3 效率分析

圖2給出了GPU中單站三維離散化并行模塊的執(zhí)行流程。對照表1,GPU并行計算中包括內(nèi)存管理、數(shù)據(jù)傳輸和計算模塊三個部分。表2給出了一個單站三維離散化實例(1 km分辨率)和多站拼圖(213個站)實例的運行時間分布。表中可以看出單站的三維離散化總計耗時0.483 s,考慮到多站雷達數(shù)據(jù)三維離散化之間采用串行方式,以表2的執(zhí)行速度,處理全國大約200個站的數(shù)據(jù)大概需要80～90 s,再算上CPU解析處理雷達基數(shù)據(jù)的時間,GPU并行算法和CPU并行算法相比,沒有明顯的效率優(yōu)勢。由表2可以看出單站離散化過程中有超過99%的運行時間都用在了內(nèi)存管理和數(shù)據(jù)傳輸方面。所以從優(yōu)化算法和提高效率的角度考慮,應(yīng)提高內(nèi)存管理和數(shù)據(jù)傳輸部分的效率。多站拼圖并行模塊的時間分布與此類似,不做贅述。

圖2 基于CUDA的單站三維離散化并行模塊Fig.2 Parallel module of radar archives discretization based on CUDA

表2 基于 CUDA的組網(wǎng)拼圖算法中的耗時分配(單位:s)Table 2 Time consuming in each step of radar mosaic algorithm based on CUDA (unit: s)

2.2.4 改進方案

通過效率分析,本文從內(nèi)存管理和數(shù)據(jù)傳輸兩個方面提出CUDA組網(wǎng)拼圖算法改進方案。

改進1:以計算代替存儲?？紤]到GPU計算的耗時明顯少于內(nèi)存管理的耗時,盡可能減少數(shù)據(jù)存儲和設(shè)備間的數(shù)據(jù)拷貝,進而減少內(nèi)存管理開銷。通過表2可以看出,單站離散化中的空間插值表計算和三維離散化插值計算總計用時不到0.001 s。由于存儲查找表需要在GPU中申請單獨的大數(shù)組,其內(nèi)存管理時間遠遠大于計算時間。由此,將空間查找表計算和三維離散化插值合并成一個并行計算模塊,格點的坐標(biāo)系轉(zhuǎn)換關(guān)系直接存儲在GPU中。通過測試,應(yīng)用改進1后,分配空間的耗時由0.269 s下降至0.202 s。

改進2:增加數(shù)據(jù)結(jié)構(gòu)復(fù)用,避免重復(fù)申請和釋放GPU內(nèi)存空間?？紤]到CPU以串行的方式調(diào)度GPU進行不同站的三維離散化計算,在不同站處理時,部分?jǐn)?shù)據(jù)結(jié)構(gòu)可以直接進行復(fù)用,減少GPU內(nèi)存管理開銷。在算法開始時,可以一次性申請全部內(nèi)存空間,并在最后統(tǒng)一進行空間釋放,減少內(nèi)存管理的操作次數(shù)。通過測試,應(yīng)用改進2后,200個站的內(nèi)存申請耗時0.428 s,內(nèi)存釋放耗時0.192 s,效率提升非常明顯。

改進3:以存儲代替?zhèn)鬏敗＝M網(wǎng)拼圖算法的第二步是將單站離散化之后的數(shù)據(jù)融合拼接到拼圖網(wǎng)格上,需要將單站離散化結(jié)果從單站模塊搬運到多站拼圖模塊。此處考慮到由于兩個并行模塊都在GPU中完成,可以將單站離散化結(jié)果臨時存儲在GPU內(nèi)存中,拼圖模塊直接從GPU內(nèi)存讀取相關(guān)數(shù)據(jù)進行拼圖。應(yīng)用改進3后,可以省去單站三維離散化模塊的數(shù)據(jù)返回耗時和多站拼圖模塊的數(shù)據(jù)傳輸耗時。

表3給出了改進前后的GPU并行開銷對比。在四個對比項目中,除了CPU—GPU傳輸數(shù)據(jù)部分是必要開銷外,其余三部分都得到了明顯優(yōu)化,并行開銷總計時間也得到了顯著下降。

表3 改進前后GPU組網(wǎng)算法并行開銷對比(單位:s)Table 3 Comparison of parallel overhead of GPU networking algorithms before and after optimization based on CUDA (unit: s)

2.2.5 基于CUDA的并行組網(wǎng)拼圖算法設(shè)計

綜上,圖3給出了CUDA并行拼圖算法的最終設(shè)計,圖中n表示雷達站點數(shù)。和2.2.2節(jié)相比,最終設(shè)計增加了對于GPU內(nèi)存的管理和CPU—GPU數(shù)據(jù)傳輸?shù)膬?yōu)化,即上文提到的三個優(yōu)化方案。GPU完成單站三維離散化計算后,僅發(fā)送回調(diào)消息,通知CPU開始下一個雷達站數(shù)據(jù)的處理。拼圖模塊直接從GPU內(nèi)存中讀取數(shù)據(jù),完成計算,最后將拼圖產(chǎn)品返回給CPU端。通過數(shù)據(jù)流和內(nèi)存管理的優(yōu)化,組網(wǎng)算法一共只有n次CPU—GPU數(shù)據(jù)拷貝,1次GPU—CPU的數(shù)據(jù)拷貝,1次申請空間和1次釋放空間。

圖3 基于CUDA的并行拼圖算法設(shè)計Fig.3 Flow chart of radar mosaic algorithm based on CUDA

3 組網(wǎng)拼圖算法對比試驗

3.1 對比試驗介紹

SWAN實時業(yè)務(wù)逐6分鐘生成全國(18°～55°N、70°～135°E)三維(垂直20層)拼圖產(chǎn)品,每個時次大約處理210個雷達(約2 GB)的數(shù)據(jù)。以此業(yè)務(wù)場景為基礎(chǔ),本文設(shè)計兩組對照試驗,驗證CPU并行和GPU平行拼圖算法的效率。試驗內(nèi)容如表4所示,其中CPU并行算法采用SWAN系統(tǒng)中的拼圖算法(30個線程并行)。試驗在同一臺服務(wù)器上進行,其中CPU為E5-2690,內(nèi)存64 G,顯卡為GTX1080,顯存8 G。耗時采用5次運算取均值的方法。

表4 對比試驗Table 4 Comparison experiments of algorithm efficiency

3.2 運行效率分析

本文按照公式1進行算法運行效率分析:

ttotal=tdecode+tstation+tmosaic

(1)

式中:ttotal表示算法的總運行時間,tdecode表示解析全部雷達基數(shù)據(jù)的耗時,tstation表示全部單站完成三維離散化的耗時,tmosaic表示多站拼圖的耗時。

使用加速比SG定量評估CUDA并行算法加速效果:

SG=tC/tG

(2)

式中:tC表示CPU并行算法的運行時間,tG表示CUDA并行算法的運行時間。SG越大,說明CUDA算法加速效果越好。

表5給出了算法效率試驗結(jié)果。對比ttotal可以看出,CUDA并行算法的效率優(yōu)勢明顯,SG分別達到3.52和6.82。對比試驗一、試驗二可以看出,隨著拼圖分辨率的提高,CUDA算法的效率優(yōu)勢也更加明顯。這主要是由于在當(dāng)前拼圖的計算量下,CUDA算法的耗時主要是內(nèi)存管理和數(shù)據(jù)傳輸。當(dāng)拼圖分辨率提高時,GPU內(nèi)存管理和拼圖產(chǎn)品的傳輸開銷變大,并行計算部分的耗時還是僅占不到10%。這說明CUDA并行算法具備應(yīng)對更大規(guī)模計算的潛力。表中tdecode主要包括CPU解壓和解析雷達基數(shù)據(jù)的時間,此步驟在四個試驗中是基本一致的。從試驗結(jié)果也可以看出,四個試驗方案的耗時接近,符合預(yù)期。如果排除tdecode,僅考慮計算部分的加速比,SG分別達到11.8和12.4。表中的試驗二表明,基于CUDA的并行拼圖算法可以在40 s內(nèi)完成全國500 m分辨率拼圖的運算,說明應(yīng)用GPU并行技術(shù)之后,更高分辨率的雷達拼圖具備一定的業(yè)務(wù)可行性。

表5 算法效率試驗結(jié)果(單位:s)Table 5 The result of algorithm efficiency experiments (unit: s)

4 結(jié) 論

本文基于CUDA設(shè)計一種CPU+GPU并行的全國雷達組網(wǎng)拼圖算法,其中CPU負(fù)責(zé)文件解析和流程調(diào)度,GPU負(fù)責(zé)大規(guī)模的插值并行計算。在算法設(shè)計時,考慮到GPU特性和拼圖算法特點,通過優(yōu)化GPU內(nèi)存管理和數(shù)據(jù)流程,減少CUDA算法的并行開銷,和CPU并行算法(30個線程)相比,獲得了超過3.5的加速比。

本文的試驗表明,CUDA并行拼圖算法可以在40 s內(nèi)完成全國500 m分辨率拼圖的實時計算,為高分辨率雷達資料在業(yè)務(wù)中的應(yīng)用提供了一種新方式。本文設(shè)計的雷達組網(wǎng)拼圖算法,可以代表計算密集、數(shù)據(jù)吞吐量大的這一類氣象算法的CUDA并行方案。此類算法都可以通過本文類似的方法實現(xiàn)CUAD的并行化改造。

在今后的工作中,考慮將拼圖算法和其前序的單站質(zhì)控、后序的外推算法進行進一步整合,直接實現(xiàn)基于GPU的SWAN算法序列。此外,配合中國氣象局雷達流傳輸業(yè)務(wù)改革,拼圖算法可以直接對接雷達基數(shù)據(jù)數(shù)據(jù)流,通過在體掃觀測結(jié)束前,預(yù)先加載部分?jǐn)?shù)據(jù)的方式,進一步提高產(chǎn)品時效。