鮑婷婷，焦圣明，殷笑茹，陳景麗，牛靄琛

(江蘇省氣象信息中心，江蘇 南京 210008)

0 引 言

等值線圖是一種重要圖形，在地質、氣象、水文、設計等領域中得到了廣泛應用。在數值分析中，繪制等值線圖是應用最多的基礎技術手段之一。隨著網絡帶寬的日益完善，傳統(tǒng)的等值線圖繪制應用已明顯不能滿足用戶需求，探索基于Web或云計算下的繪制算法成為目前的研究主流。在軟硬件環(huán)境及算法相同的情況下，傳統(tǒng)應用程序在執(zhí)行效率上要比基于Web的應用高，然而在離散點數量多、網格點密度大的情況下，按傳統(tǒng)的串行方法實現等值面生成的這兩類應用均存在響應速度慢、效率低的問題，這一問題同時也廣泛存在于專業(yè)商業(yè)軟件中，使得業(yè)務系統(tǒng)在實時性方面無法得到保證。隨著多核處理器設計技術的快速發(fā)展，針對Web應用研究基于多核平臺下的等值面快速生成并行算法顯得尤為迫切。

等值面的繪制過程一般分為離散數據網格化、等值點計算、等值線追蹤及光滑、等值面提取、地圖裁剪等幾個步驟。當前，國內外就如何提高等值線圖生成效率展開了大量的研究，其可以從兩方面著手：一方面把GPU等硬件資源應用到普通的數值并行計算中；另一方面通過軟件方法提高計算速度，主要是通過并行算法改進程序提升效率[1-2]。有學者利用GPU提供大量的并行數據，應用其并行性來進行通用計算，如使用GPU進行等值面的提取與繪制，但其繪制過程無法完全在圖形硬件上實現[3-4]；文獻[5]在GPU并行架構和CUDA可編程性基礎上提出了一種基于區(qū)間塊搜索的等值線并行提取算法，但提取方法并未給出；文獻[6]將基于非規(guī)則四面體的等值面提取和繪制過程完全轉移至GPU上執(zhí)行，從而降低了CPU計算負荷，提高了等值面提取和繪制速度；文獻[7]提出了一種等值線二值邊界追蹤的算法，以柵格邊界為基礎追蹤等值線；文獻[8]針對等高線繪制中排序耗時問題提出了一種改進的算法，可應用于2D與3D對象；文獻[9]采用分區(qū)法和OPENMP的并行算法處理形成等值填色場，但文獻中實現的等值面只是根據并行計算插值給圖片像素賦予顏色值，生成的等值面是圖像產品而不是矢量數據，在業(yè)務使用上存在一定的局限性；文獻[10]是基于ArcGIS Server的模型來實現等值面生成，其成本高、算法可擴展性小，在等值面圖形復雜度較高的情況下響應性能會受到一定的影響，致使效率不高。

針對以上問題，文中提出了一種基于Fork/Join框架的等值面快速生成并行算法。首先描述算法的基本原理、并行任務及劃分策略，其次以自動氣象站的累積降水量和不同網格下數值預報溫度場的等值面提取為例實現該并行算法，并通過與傳統(tǒng)串行算法比較，從計算性能、等值面繪制效果驗證該算法的效率及優(yōu)勢。

1 Fork/Join框架

Fork/Join框架是一個分治(divide and conquer)策略的并行編程框架，實現了任務定義與任務處理的分離，大大簡化了并行程序設計的復雜度，開發(fā)人員能夠方便地利用多核平臺的計算能力。盡管還沒有做到對開發(fā)人員完全透明，但已經極大地簡化了編寫并行程序的瑣碎工作。

工作竊取算法(work-stealing algorithm)是該框架用來提高性能、保證負載均衡的一種有效方法，其本質上是一種工作任務的調度方法?；舅枷胧抢靡淹瓿扇蝿盏木€程去查看其他線程是否有未處理完的工作，如果有則竊取一部分工作來執(zhí)行，從而提高資源的利用率，減少程序的處理時間。

該框架主要由Fork和Join兩個操作構成，Fork操作是對任務的劃分，把大任務劃分成若干小問題，能方便對簡化后的問題進行處理；Join操作主要是完成對各個部分的運行結果進行合并的任務。其工作原理如圖1所示。

圖1 Fork/Join框架工作原理

2 等值面并行算法

2.1 串行處理算法

目前，規(guī)則矩形網格法和不規(guī)則三角網格法是進行等值線圖繪制比較常見并且較成熟的算法[11]，由于氣象業(yè)務中大部分數據是已經插值好的網格數據，故多采用規(guī)則矩形網格法進行等值線圖繪制?；谝?guī)則矩形網格法的等值面提取算法描述如下：

算法：規(guī)則矩形網格法等值面提取算法。

輸入：某種氣象要素的離散數據集U，等值線間隔數組K；

輸出：等值面矢量數據集。

Step1:將離散點氣象數據集利用一定的插值方法進行網格化處理；

Step2:計算網格上的等值點，并將等值點連接成線段存貯到線段組隊列；

Step3:利用隊列中的線段組根據等值間隔追蹤出所有等值線，并進行光滑處理；

Step4:通過邊界線把非閉合等值線處理成閉合等值線，對等值面進行標記識別；

Step5:按目標區(qū)域范圍進行等值面裁剪；

Step6:算法結束。

2.2 可并行過程

并行算法是一些可同時執(zhí)行的諸進程的集合，這些進程互相作用和協(xié)調動作從而達到給定問題的求解[12]。并行基本思想就是在一段時間內同時處理兩件或兩件以上的作業(yè)，這兩件或兩件以上的作業(yè)可以是同一性質也可以是不同性質的，如果要并行計算就必須保證過程之間的獨立性，即不包括數據依賴、控制依賴、數據交換等關系。從整體上分析程序的并行化工作，綜合考慮整個計算任務在所有計算資源上的分配，盡量做到負載平衡，可以取得較大的性能提升[13]。就等值面提取算法而言，其處理流程的各個步驟是有一定依賴關系的，它們不能同時并行處理串行算法中的各個步驟，而是要按一定的先后順序執(zhí)行。因此，并行算法只能作用在可以實施并行計算的單獨步驟上，從而減少這些步驟的作業(yè)時間，達到提高程序執(zhí)行效率的目的。串行算法里有四個步驟可以實施并行計算，分別為：

(1)離散點數據網格化處理。離散數據量的多少、網格設定的大小以及不同的插值方法都會影響到算法的處理速度。在傳統(tǒng)方法下，網格密度越大，插值精度越高，算法花費在格點值的計算時間就越長，反之網格密度越小，計算速度越快，但輸出的等值面失真度就越高。精度和效率都是不可缺少的，為了同時保持高精度和高效率，一個有效的解決辦法就是使用并行算法來加快網格點的計算速度。

(2)等值點計算。根據步驟1處理好的格點資料計算出有效的等值線間隔的最大值和最小值，確定出等值點判別循壞的起始值，在并行計算等值點時可能會減少當前等值線段組屬于哪條等值線的標識時間，它在最差狀態(tài)與上面不進行預處理時的循環(huán)次數相同。

(3)等值線追蹤與光滑。從步驟2線段組中提取出所有不同間隔的等值線值，分別建立新的隊列，將原先的線段組擺放到各自相對應的新隊列中去，接著對這些新隊列進行并行等值線追蹤，最后按文獻[14]中的方法對每條等值線進行穿過原點的光滑處理。

(4)等值面標記識別。僅根據等值線的標記值是無法判斷由其構成的等值面等級，需要利用格點上的值對等值面進行判別，標記識別的目的是確定使用什么顏色對該等值面進行顏色填充。文中對標記識別采用一種簡單實用的辦法，就是先對所有等值面進行面積大小排序，排序完成后從等值面邊界線任取一點計算出相對應的格點，利用該格點及周邊四個格點進行判別，如果該格點在等值面內，就以此格點的值換算出等值面的標記值。通過邊界線把非閉合等值線處理成閉合等值線，把所有排好序的等值面按順序存放到數組隊列，然后把這一大數組隊列劃分成數各小數組隊列并行標記識別。

當用戶對實時性要求比較高時，若操作響應時間超過3秒，則其交互體驗效果將會變差。因此，對上述四個步驟進行并行化處理后，計算速度明顯提升，有效地保證了客戶端用戶交互響應的時效性。

2.3 任務劃分策略

通常算法中循環(huán)處理操作是程序執(zhí)行中蘊含并行性最豐富的一種結構，上述四個步驟是串行程序的并行化，就是將可以并行化的循環(huán)移植到并行環(huán)境中進行計算的過程。一個大任務操作在劃分成許多小任務時，也不是并行小任務越多處理速度就越快，需要從硬件資源和計算任務性質統(tǒng)籌考慮，并行任務線程池的線程數主要依據硬件資源設定，計算性質對并行閾值(threshold value)取值也十分重要，在任務規(guī)模不大的情況下采用串行的解決方法可能會更好。上述步驟1、2計算量級相似，其劃分小任務的閾值基本一致但不能太小，可以設定在3 000～4 000；步驟3隊列數組不會太大但是一個較耗時的操作，分解的并行小任務數量和數組長度一樣；步驟4需要幾何空間，計算閾值不宜太大，可以設定在10～30。這四個并行計算都使用了工作竊取算法。

2.4 并行算法實現

Java從JDK1.7版本開始引入Fork/Join框架編程模式，采用這個并行框架最大的優(yōu)勢是能夠高效快捷地實現程序跨平臺功能。Fork/Join框架里的java.util.concurrent包提供ForkJoinTask和ForkJoinPool兩個類來實現算法的并行，通常情況下不需要直接繼承ForkJoinTask類，而只需要繼承它的子類RecursiveAction和RecursiveTask。RecursiveAction類用于沒有返回結果的任務，RecursiveTask用于有返回結果的任務。上述步驟1～步驟3的并行算法都是繼承有返回任務的RecursiveTask類實現。ForkJoinTask需要通過ForkJoinPool來執(zhí)行，任務分割出的子任務會添加到當前工作線程所維護的雙端隊列中，進入隊列的頭部,當一個工作線程的隊列里暫時沒有任務時，它會隨機從其他工作線程的隊列的尾部獲取一個任務，原理實現如圖2所示。

通過對上述等值面提取算法及其并行性的分析，得出算法核心過程的任務可分割性，提取了各個任務中的可并行部分，得到基于Fork/Join框架編程模式下實現的并行算法流程示意圖(見圖3)。文中基于Fork/Join框架編程模式實現的并行算法最核心的偽代碼如下：

public class TaskName extends RecursiveTask {//V是返回結果的變量

//變量定義

……

public TaskName(//傳入參數){

//類初始化

}

@Override

protected V compute( ){

//定義返回的變量

V sum；

if(問題規(guī)模<閾值){

//使用串行模式解決或選擇其他更優(yōu)算法解決

……

}

else{

//將任務Task進行分解,分解成Task1、Task2…

TaskName Task1=new TaskName();

圖2 Work-stealing(工作竊取)實現

圖3 等值面提取并行算法流程示意圖

TaskName Task2=new TaskName();

……

//將任務Task1、Task2…提交給線程池執(zhí)行

Task1.fork();

Task2.fork();

……

//如果任務有返回結果，收集結果合并

V v1=Task1.join();

V v2=Task2.join();

sum=v1+v2;

}

return sum;

}

}

等值面提取算法是通過編寫com.jandy.isolines類包來實現的，其中IsoLines是主類，定義了一系列屬性和方法供程序調用，最重要的方法是parallelMakeIsolines(),它包含了四個并行過程。在com.jandy.isolines.parallel子包下，這四個并行過程的類就是利用上述偽代碼思想編碼完成，只不過各自類的輸入任務、計算方法和輸出結果不一樣而已。這個包單獨編譯后形成的jar文件還可以提供其他Java開發(fā)人員當作組件使用，提高了代碼的復用性。

3 實 驗

3.1 算法處理時間

等值面提取算法的計算時間是可以進行估算的。

設算法處理時間為T，串行算法處理時間為：

T=Ti+Tg+Tc+Ts+Tp+Tw+Tm

(1)

其中，T為處理總的運行時間，Ti為讀取數據處理時間，Tg為離散數據網格化處理時間，Tc為等值點計算時間，Ts為等值線追蹤和光滑時間，Tp為非閉合等值線處理成閉合等值線和面積排序時間，Tw為等值面標記識別時間，Tm為地圖裁剪時間。

并行算法處理時間為：

(2)

其中，N是處理器核個數，算法通過把各個任務分配到N個核處理器上以相同方法同時進行計算，從而減少程序的計算時間，提高執(zhí)行效率。實驗數據分析主要針對可以并行的過程，不能并行的過程執(zhí)行時間相似且極小，基本可以忽略不計。

根據上述串行算法處理時間(T1)與并行算法處理時間(T2)得到的兩者所耗時間比(加速比)來評價該并行程序的性能，即，并行執(zhí)行與串行執(zhí)行時的計算性能的比值[15]:

(3)

3.2 實驗平臺

等值面提取算法的性能測試是在軟硬件環(huán)境相同的實驗平臺上進行的，通過串、并行算法實例運行，比較兩者的計算時間。采用的硬件平臺是一臺配備Intel(R) Xeon(R) E7-4820 v2 主頻2.00 GHz的8核CPU處理器和128 G內存的X3860服務器,操作系統(tǒng)為64位Redhat Linux 6.0，Java虛擬機為64位JDK1.8。為了避免其他軟件的干擾，對操作系統(tǒng)及測試必要的軟件都進行了全新安裝。

3.3 實驗結果與分析

實驗中分別使用了離散站點的累積降水量和不同網格大小的數值預報溫度場這兩類資料進行測試和分析，之所以選擇這兩種類型的資料是因為它們在實際的氣象業(yè)務應用中極具代表性。為了評估并行計算對等值面處理效率的提高，使用加速比作為衡量標準，執(zhí)行時間單位為毫秒。

(1)離散數據。

離散數據使用的是江蘇省氣象局全省地面自動站2015年6月某日的累積降水量數據，對71個基本站和含加密站的1 441個所有站的兩組數據分別進行了實驗。網格大小設定為201×201，經度116°～122°，緯度30°～36°，空間分辨率遠小于加密站分布密度能夠滿足插值精度，生成了7級降水分布圖，測試不同站點數量對插值速度的影響，見表1。

表1 71個基本站和含加密站的1 441個所有站的等值面提取串并行性能對比

從表1可以看出，上述離散數據的等值面提取算法四個并行過程Tg、Tc、Ts和Tw，無論是串行還是并行最耗時間的過程是離散數據網格化處理過程Tg，在71個基本站這一組數據中Tg的串、并行執(zhí)行時間分別占總處理時間的96.34%和90.88%；在含加密站的這一組數據中分別為99.46%和94.98%；隨著參與插值計算站點個數的增加，Tg過程的計算時間呈增加趨勢但加速比獲得顯著提高，但其他三個過程的加速比就不那么明顯甚至會降低。因此對于小網格、等值線稀疏的場景下更適合使用串行算法，并行算法可能會適得其反，增加了處理時間。針對此類業(yè)務應用是采取串行還是并行需要具體對待，但離散數據網格化處理過程Tg采用并行算法處理是毫無疑問的。圖4(a)是以71個基本站數據進行插值計算的累積降水量分布圖，圖4(b)是以含加密站的1 441個所有站數據進行插值計算的累積降水量分布圖。圖4(b)更能反映出較小尺度下的降水差異性，提高了空間精度，與圖4(a)相比在處理時間上要長些，但采用并行算法處理后生成時間控制在2秒以下能夠滿足業(yè)務要求。

圖4 降水分布對比

(2)網格數據。

網格數據是使用美國、歐洲、日本和中國的數值預報中心輸出的2米高溫度場格點預報產品作為實驗數據，產品數據都是采用氣象信息綜合分析處理系統(tǒng)(meteorological information comprehensive analysis and process system，MICAPS)標準的第四類文件格式進行存儲，十分方便業(yè)務人員進行讀寫操作。提取2016年3月06日16時四家數值預報中心全球模式不同分辨率下的溫度場數據，等值線間隔設為5 ℃，溫度范圍設定在-60～40 ℃，對比串行與并行算法的耗時，生成20級等值面，測試不同網格大小對加速比的影響，見表2。

表2 不同網格大小等值面提取串并行性能對比

圖5 不同網格下的串、并行處理時間及加速比

根據表2結果繪制如圖5所示的折線圖。在網格規(guī)模較小的情況下，并行算法的性能突出不明顯，但網格規(guī)模加大時，并行算法效率提升的效果十分顯著。從表2中可以看出，在1 440×721網格下達到最大加速比5.308 1，網格大小從720×361到1 440×721，串行處理時間增加了5.6秒，而并行處理時間只增加了1.0秒，串行算法的處理時間增長幅度遠高于并行算法。

圖6為利用文中并行算法繪制的網格大小為1 440×721的溫度要素等值面效果圖，從讀取數據到繪制完成用時不到2秒，而用專業(yè)氣象業(yè)務軟件MICAPS繪制耗時近15秒。

圖6 1 440×721網格溫度要素個例

4 結束語

通過對常用的兩類氣象資料等值面生成的特點分析，研究了基于Fork/Join框架下多核并行計算數據處理的方法。實驗結果表明，該方法在實現復雜計算量下的等值面提取時較大地提高了數據處理速度，等值面的生成時間控制在2秒內，能夠滿足氣象業(yè)務需求，在生成效率上具有優(yōu)越性，為其他相關業(yè)務的等值面分析提供了一種借鑒思路。特別是針對實時性要求極高、基于Web應用的等值面分析時，在服務端使用并行計算是一種較好的解決方案。后續(xù)工作將利用此算法構建云服務平臺，為注冊用戶開放API云接口，便于用戶在自己的Web應用中接入使用。