李 健，張大偉，姜曉明，向立云

1.中國(guó)地質(zhì)大學(xué)（武漢）海洋學(xué)院，武漢430074

2.中國(guó)水利水電科學(xué)研究院防洪減災(zāi)所，北京100038

3.水利部防洪抗旱減災(zāi)工程技術(shù)研究中心，北京100038

洪水風(fēng)險(xiǎn)圖編制是防汛減災(zāi)非工程措施的重要手段，對(duì)降低洪水災(zāi)害至關(guān)重要。但在實(shí)際應(yīng)用中，洪水、降雨和地形條件復(fù)雜多變，實(shí)際汛情與計(jì)算方案情景差別很大，制約了洪水風(fēng)險(xiǎn)圖的實(shí)際應(yīng)用效果[1]。實(shí)施可針對(duì)多種情景、任意邊界條件、實(shí)時(shí)演示和人機(jī)交互式計(jì)算的洪水演進(jìn)模擬，為防汛決策等提供動(dòng)態(tài)、實(shí)時(shí)的洪水淹沒(méi)信息，具有重要的現(xiàn)實(shí)意義。但是，基于水動(dòng)力計(jì)算的洪水演進(jìn)模擬計(jì)算量較大時(shí)，模擬耗時(shí)延長(zhǎng)，將妨礙洪水淹沒(méi)的實(shí)時(shí)演示和防洪決策，并行化計(jì)算可解決此類(lèi)問(wèn)題[1]。

洪水演進(jìn)模擬是計(jì)算水力學(xué)領(lǐng)域的一個(gè)重要研究方向，涉及繁雜的數(shù)值算法、研究手段和編程技術(shù)等，目前已幾乎不可能對(duì)其進(jìn)行全面詳盡的綜述。然而，實(shí)際應(yīng)用中多采用基于Godunov 型有限體積算法和網(wǎng)格計(jì)算的洪水演進(jìn)模型，涉及到復(fù)雜地形處理、干濕前鋒追蹤、通量格式和諧性、模型魯棒性等問(wèn)題[2]，對(duì)此國(guó)內(nèi)外學(xué)者從不同角度撰寫(xiě)了大量綜述性論文，包括：洪水模擬的不同物理機(jī)制[3]、洪水演進(jìn)模型的復(fù)雜度分析[4]、洪水淹沒(méi)的不確定性分析[5]、城市洪水淹沒(méi)模擬[6]。近年來(lái)計(jì)算機(jī)硬件和軟件的發(fā)展，尤其異構(gòu)并行計(jì)算技術(shù)的發(fā)展，極大地推動(dòng)了洪水演進(jìn)模型的研發(fā)。但是，目前僅有覃金帛等[7]介紹GPU并行優(yōu)化在水文計(jì)算和水庫(kù)調(diào)度研究中的應(yīng)用，未見(jiàn)有系統(tǒng)性介紹并行化洪水演進(jìn)模擬的綜述性論文。因此，有必要及時(shí)總結(jié)洪水演進(jìn)的并行化研究成果，促進(jìn)今后的洪水演進(jìn)過(guò)程的并行化模擬研究。

考慮到洪水演進(jìn)模擬研究的主流性、前沿性和實(shí)際應(yīng)用價(jià)值，本文將介紹基于網(wǎng)求解淺水方程的并行化洪水演進(jìn)模型，尤其是非結(jié)構(gòu)網(wǎng)格模式。由于粒子法應(yīng)用于實(shí)際河道或流域的洪水模擬時(shí)，存在邊界條件難以定義和計(jì)算量大等問(wèn)題，以及計(jì)算機(jī)圖形學(xué)中的流體模擬，側(cè)重計(jì)算機(jī)可視化和高性能制圖算法，而流體模擬多采用無(wú)條件穩(wěn)定的快速傅里葉變換（FFT）等算法，并不是求解流體方程，因此，本文綜述不包括上述的兩種數(shù)學(xué)模型，粒子法或計(jì)算機(jī)圖形學(xué)領(lǐng)域的并行化洪水演進(jìn)模擬研究綜述可參考文獻(xiàn)[8]。本文將結(jié)合近年洪水演進(jìn)模擬的最新數(shù)值算法和并行化技術(shù)，進(jìn)行文獻(xiàn)綜述、討論和展望，期望為今后快速洪水預(yù)報(bào)研究起到拋磚引玉的作用，為防洪減災(zāi)的實(shí)時(shí)會(huì)商和洪水風(fēng)險(xiǎn)圖的編制提供技術(shù)支撐。

1 洪水演進(jìn)模擬的并行化需求

近年來(lái)洪水演進(jìn)模型被廣泛應(yīng)用于評(píng)估洪水造成的經(jīng)濟(jì)損失，洪水演進(jìn)模型包括1D、2D 和3D 模型，在復(fù)雜度和預(yù)報(bào)能力等方面存在較大差異。1D模型可提供的信息有限，一般僅用于模擬復(fù)雜河網(wǎng)或與2D 模型耦合，實(shí)現(xiàn)區(qū)域洪水模擬，但存在人為處理邊界條件的問(wèn)題[9]。2D 模型采用結(jié)構(gòu)網(wǎng)格或非結(jié)構(gòu)網(wǎng)格求解淺水方程，精度高但計(jì)算量增大[10-12]，可采用自適應(yīng)網(wǎng)格加密（AMR）技術(shù)[13]和3D 模型[14]的方法來(lái)克服，但實(shí)際應(yīng)用中的計(jì)算量任然較大。因此，基于網(wǎng)格計(jì)算的并行化洪水演進(jìn)模型的研發(fā)和應(yīng)用勢(shì)在必行。

2 洪水演進(jìn)模型的并行化模式

常見(jiàn)的并行模式見(jiàn)表1，區(qū)別包括：并行層次（數(shù)據(jù)集、線程和進(jìn)程）、并行粒度（細(xì)粒度和粗粒度）、硬件支持、軟件實(shí)現(xiàn)和編寫(xiě)并行化代碼的編程語(yǔ)言。本文涉及的硬件、軟件和并行計(jì)算的計(jì)算機(jī)術(shù)語(yǔ)，僅列出對(duì)應(yīng)的中文全稱(chēng)和英文縮略詞，具體含義可參考文獻(xiàn)[15-16]，后文中介紹的洪水演進(jìn)模型使用的并行模式與表1 所列相對(duì)應(yīng)。如果計(jì)算完全在中央處理單元（CPU）中完成，稱(chēng)為同構(gòu)并行計(jì)算，如果計(jì)算部分執(zhí)行于CPU，部分執(zhí)行于其他硬件設(shè)備，如圖形處理單元（GPU），稱(chēng)為異構(gòu)并行計(jì)算。從并行層級(jí)的角度，表1列出了從最低層級(jí)的數(shù)據(jù)并行，到線程級(jí)和進(jìn)程級(jí)并行。

表1 常見(jiàn)的并行技術(shù)列表

3 并行化洪水演進(jìn)模擬研究現(xiàn)狀

3.1 同構(gòu)并行化洪水演進(jìn)模型

洪水演進(jìn)模型根據(jù)空間網(wǎng)格模式，可分為結(jié)構(gòu)網(wǎng)格和非結(jié)構(gòu)網(wǎng)格；根據(jù)時(shí)間離散格式，可分為顯格式和隱格式；根據(jù)對(duì)控制方程的不同簡(jiǎn)化，可分為求解儲(chǔ)水單元、擴(kuò)散波、淺水方程和靜水壓力假設(shè)的Navier-Stokes方程等類(lèi)型的模型；根據(jù)數(shù)值方法，可分為有限差分法（FDM）、有限體積法（FVM）和有限單元法（FEM）等模型。針對(duì)不同類(lèi)型的洪水演進(jìn)模型，需要采用不同的并行模式。本文將根據(jù)以上分類(lèi)，從并行計(jì)算的角度，選擇有代表性的模型列于表2，其中不包括未得到廣泛應(yīng)用或測(cè)試或未并行化的洪水演進(jìn)模型[11-12，17-29]。另外，表2所列模型還考慮了代碼的可獲取性，方便讀者選擇使用，檢驗(yàn)洪水淹沒(méi)模擬結(jié)果的可重復(fù)性及模型對(duì)比分析。

表2 一些有名的同構(gòu)并行化洪水演進(jìn)模型列表

較早的洪水模型多使用結(jié)構(gòu)網(wǎng)格和FDM，如LISFLOOD-FP和FloodMAP模型，結(jié)構(gòu)網(wǎng)格模型數(shù)據(jù)結(jié)構(gòu)簡(jiǎn)單，易于并行化編程，其中LISFLOOD-FP模型采用了OpenMP、MPI 和Clearspeed 三種并行機(jī)制，而Flood-MAP 采用多線程技術(shù)，LISFLOOD-FP 和FloodMAP 兩種模型采用顯格式近似求解黎曼問(wèn)題的數(shù)值通量的方法，允許使用的計(jì)算時(shí)間步長(zhǎng)較小，必須并行化才能應(yīng)用于精細(xì)模擬區(qū)域洪水淹沒(méi)。需要指出的是，不同的并行化洪水演進(jìn)模型，求解的控制方程，使用的數(shù)值方法和并行化機(jī)制均不同，籠統(tǒng)地比較不同模型的并行計(jì)算效率并無(wú)意義，但并行化結(jié)構(gòu)網(wǎng)格模型研究為之后的并行化模型研發(fā)奠定了基礎(chǔ)。

近年來(lái)，洪水演進(jìn)模擬已由結(jié)構(gòu)網(wǎng)格計(jì)算轉(zhuǎn)向非結(jié)構(gòu)網(wǎng)格計(jì)算的模式開(kāi)發(fā)，而FVM在單元界面上近似求解黎曼問(wèn)題，具有較強(qiáng)的洪水激波捕捉能力，如MIKE21[17]、CalTWiMS[18]、ANUGA[19]、TUFLOW[20]和ParBreZo[21]等模型，但缺點(diǎn)是：

（1）對(duì)網(wǎng)格質(zhì)量要求高，網(wǎng)格的正交性嚴(yán)重影響數(shù)值格式的計(jì)算精度，進(jìn)而影響計(jì)算穩(wěn)定性。

（2）采用顯格式，當(dāng)網(wǎng)格尺寸很小時(shí)，受CFL條件限制，允許使用的計(jì)算時(shí)間步長(zhǎng)很小，嚴(yán)重影響計(jì)算效率。

（3）數(shù)值通量格式涉及大量數(shù)組循環(huán)計(jì)算。

（4）采用靜態(tài)網(wǎng)格，可能在干濕地形劇烈變化的區(qū)域網(wǎng)格密度不足，影響計(jì)算精度。

采用FEM離散方法可克服第（1）點(diǎn)問(wèn)題，并行化計(jì)算可克服第（2）、（3）點(diǎn)問(wèn)題，其中基于OpenMP 的細(xì)粒度線程并行可加速計(jì)算數(shù)組循環(huán)計(jì)算，結(jié)合區(qū)域分解后的MPI 進(jìn)程并行，可取得更好的加速計(jì)算效果?；贏MR動(dòng)態(tài)網(wǎng)格技術(shù)可克服第（4）點(diǎn)問(wèn)題，實(shí)現(xiàn)自動(dòng)加密局部網(wǎng)格密度，跟蹤干濕前鋒位置。但改進(jìn)的方法也存在以下局限性：

（1）AMR網(wǎng)格可能會(huì)產(chǎn)生正交性較差的網(wǎng)格，一般結(jié)合FEM模型使用[22]，已廣泛應(yīng)用于洪水或海嘯傳播過(guò)程的模擬，但增大了并行化的難度，特別是采用MPI 并行時(shí)，靜態(tài)非結(jié)構(gòu)網(wǎng)格可使用METIS[23]或SCOTCH[24]等庫(kù)實(shí)現(xiàn)荷載均衡的網(wǎng)格區(qū)域分解，并可降低分區(qū)間的數(shù)據(jù)通信量，提高并行計(jì)算效率，而AMR動(dòng)態(tài)網(wǎng)格密度調(diào)整會(huì)引起荷載不均衡，需要進(jìn)行分區(qū)網(wǎng)格的動(dòng)態(tài)調(diào)整，可使用Zoltan庫(kù)[25]實(shí)現(xiàn)，但編程的難度會(huì)增大。

（2）AMR非結(jié)構(gòu)網(wǎng)格還會(huì)破壞數(shù)據(jù)結(jié)構(gòu)的局部性，這對(duì)GPU異構(gòu)并行計(jì)算效率的影響尤其顯著。

（3）局部高密度網(wǎng)格會(huì)引起計(jì)算不穩(wěn)定問(wèn)題。

（4）在非結(jié)構(gòu)網(wǎng)格上難以構(gòu)建高階的空間離散格式。

（5）非結(jié)構(gòu)網(wǎng)格隱格式模型，最終形成的具有稀疏系數(shù)矩陣的線性方程組求解量巨大，成為計(jì)算瓶頸，且系數(shù)矩陣的數(shù)據(jù)局部性對(duì)求解效率影響顯著。

近年來(lái)，基于FEM 的洪水演進(jìn)模型研究得到廣泛開(kāi)展（見(jiàn)表2），盡管FEM較FVM在編程上難度更大，但FEM 適應(yīng)低質(zhì)量網(wǎng)格的能力更高，因此計(jì)算穩(wěn)定性較FVM 更好。FEM 可分為連續(xù)迦遼金（CG）法和間斷迦遼金（DG）法。CG法結(jié)合MPI并行化，一般采用隱格式離散，可高效地模擬洪水演進(jìn)，但CG法求解對(duì)流占優(yōu)的不連續(xù)問(wèn)題時(shí)，存在質(zhì)量不守恒和格式的耗散性過(guò)強(qiáng)的問(wèn)題，一般采用流線迎風(fēng)Petrov-Galerkin（SUPG）法解決，如RMA2[26]和Telemac2D[27]模型，也有基于求解通用波動(dòng)連續(xù)性方程（GWCE）的CG法來(lái)提高離散格式的質(zhì)量守恒性，如ADCIRC模型[28]。DG法不要求單元界面處滿足連續(xù)性條件，即具有FEM的靈活性又具有FVM的激波捕捉能力，同時(shí)可方便構(gòu)建高階格式，如Fluidity[22]和DGSWE[29]模型。值得注意的是，F(xiàn)EM 模型的數(shù)學(xué)原理抽象，編程難度較大，模型計(jì)算量大，需要并行化才能高精度地模擬洪水演進(jìn)過(guò)程。FEM模型結(jié)合基于網(wǎng)格質(zhì)量指標(biāo)判斷的AMR 動(dòng)態(tài)網(wǎng)格技術(shù)，實(shí)現(xiàn)洪水演進(jìn)過(guò)程中的干濕前鋒的準(zhǔn)確捕捉，可明顯提高城市洪水和局部沖刷[30]等問(wèn)題的模擬精度。求解大型線性方程組時(shí)，可有效利用多種并行模式的線性方程組求解器來(lái)提高計(jì)算效率，如PETSc 庫(kù)[31]，同時(shí)可使用SCOTCH 庫(kù)中的Cuthill-McKee 算法或Gibbs-Poole-Stockmeyer 算法，實(shí)施對(duì)網(wǎng)格單元和節(jié)點(diǎn)編號(hào)重排序，進(jìn)一步改善非結(jié)構(gòu)網(wǎng)格模型的求解效率。需要指出的是，高階DG法中的高斯求積運(yùn)算較大，需要合理設(shè)計(jì)網(wǎng)格分辨率，使計(jì)算結(jié)果在精度和計(jì)算效率上達(dá)到平衡。

3.2 異構(gòu)并行化洪水演進(jìn)模型

GPU 異構(gòu)并行計(jì)算是利用顯卡中的眾多核心實(shí)現(xiàn)大量浮點(diǎn)運(yùn)算的線程并行，非常適合用于加速FVM 洪水演進(jìn)模型中數(shù)值通量涉及的數(shù)組循環(huán)計(jì)算，GPU異構(gòu)并行的研究文獻(xiàn)在近年來(lái)快速增長(zhǎng)。下面選擇有代表性的GPU異構(gòu)并行洪水演進(jìn)模型列于表3[17，20，32-41]。

表3 一些有名的異構(gòu)并行化洪水演進(jìn)模型列表

早期的GPU異構(gòu)并行洪水演進(jìn)模型，均采用CUDA C語(yǔ)言編程，且GPU的硬件架構(gòu)可有效加速結(jié)構(gòu)網(wǎng)格的FDM 模型的計(jì)算效率，如JFlow[32]和Flood2D[33]模型。由于FVM 模型在求解不連續(xù)問(wèn)題方面的優(yōu)勢(shì)，異構(gòu)并行的結(jié)構(gòu)網(wǎng)格FVM 洪水演進(jìn)模型也大量開(kāi)展，如TUFLOW[20]和SWsolver[34]模型。非結(jié)構(gòu)網(wǎng)格模型的數(shù)據(jù)結(jié)構(gòu)較為復(fù)雜，異構(gòu)并行編程的難度更大，且數(shù)據(jù)結(jié)構(gòu)也會(huì)影響異構(gòu)并行計(jì)算效率。因此，異構(gòu)并行非結(jié)構(gòu)洪水演進(jìn)模型較異構(gòu)并行結(jié)構(gòu)網(wǎng)格模型的研究要晚，但已成為前沿研究課題，此類(lèi)模型包括商業(yè)軟件MIKE21、免費(fèi)軟件但不公開(kāi)代碼的BASEMENT[35]以及開(kāi)源代碼的VOLNA[36]模型。

異構(gòu)并行非結(jié)構(gòu)網(wǎng)格洪水演進(jìn)模型存在的局限性有：（1）CUDA 編程復(fù)雜，會(huì)破壞串行程序的代碼結(jié)構(gòu)，代碼維護(hù)成本高；（2）CUDA SDK版本進(jìn)化和GPU硬件設(shè)備更新速度快，代碼更新難度大；（3）CUDA程序只能在Nvidia 公司的顯卡上執(zhí)行，兼容性差；（4）非結(jié)構(gòu)網(wǎng)格模式下，網(wǎng)格單元、節(jié)點(diǎn)和邊的連接關(guān)系必須顯式給出，這將導(dǎo)致大量的間接尋址訪問(wèn)，對(duì)數(shù)據(jù)傳輸帶寬要求高，這在同構(gòu)并行計(jì)算中問(wèn)題不大，因?yàn)镃PU內(nèi)存的帶寬為Gb/s 量級(jí)，而GPU 顯存的帶寬一般為Mb/s 量級(jí)，因此必須使用有利于提高間接尋址訪問(wèn)效率的數(shù)組形式。

針對(duì)以上問(wèn)題，研究者們給出了解決方案，例如Zhang等[37]使用OpenACC指令集的異構(gòu)并行方法，該方法類(lèi)似于OpenMP指令的同構(gòu)線程并行方法，對(duì)代碼結(jié)構(gòu)幾乎沒(méi)有破壞，但是只能用PGI編譯器編譯。還可采用OpenCL 語(yǔ)言編程，提高模型對(duì)硬件的兼容性，如HiPIMS 模型[38]。還有研究者使用特定域語(yǔ)言（DSL）的程序庫(kù)，如OP2庫(kù)[39]，可將串行版本的程序，快速轉(zhuǎn)換為不同的同構(gòu)和異構(gòu)并行化的程序，如BASEMENT 和VOLNA 模型，DSL 轉(zhuǎn)換可讓科學(xué)家和工程人員專(zhuān)注本專(zhuān)業(yè)領(lǐng)域的研究，而不必關(guān)注并行化編程本身，同時(shí)可應(yīng)對(duì)軟件和硬件架構(gòu)快速更替的問(wèn)題。數(shù)組的排列形式分為Structure of Array（SoA）排列和Array of Structure（AoS）排列兩種，SoA排列和AoS排列對(duì)同構(gòu)并行計(jì)算效率的影響不大，但在異構(gòu)并行計(jì)算時(shí)，SoA 排列具有顯著優(yōu)勢(shì)，可提高間接尋址訪問(wèn)效率，同時(shí)實(shí)施非結(jié)構(gòu)網(wǎng)格編號(hào)的重排序。

隱格式的洪水演進(jìn)模型求解的瓶頸就是求解線性方程組，可調(diào)用異構(gòu)并行的求解器，這樣可利用具有異構(gòu)并行計(jì)算資源的集群系統(tǒng)，如Fluidity[40]和CCHE2D[41]模型。值得關(guān)注的是，采用異構(gòu)并行化的系數(shù)矩陣預(yù)處理的Krylov 空間求解器（如PETSc 庫(kù)）和多重網(wǎng)格求解器（如AMGX 庫(kù)[42]），可進(jìn)一步提高求解效率，已應(yīng)用于一些計(jì)算流體力學(xué)商業(yè)軟件，但尚未見(jiàn)到用于隱格式洪水演進(jìn)模擬的研究。

3.3 基于其他并行機(jī)制的洪水演進(jìn)模型

除了基于MPI、OpenMP的同構(gòu)并行和基于CUDA的異構(gòu)并行洪水演進(jìn)模型外，還有基于其他加速硬件設(shè)備和并行模式的洪水演進(jìn)模型，此類(lèi)并行化模式并未成為主流，但涉及最底層的數(shù)據(jù)并行或是對(duì)主流硬件加速研究的補(bǔ)充，也列于表4[11，34，43-44]?；谔厥馓幚砥鳎ㄈ鏑ELL和Clearspeed處理器）的并行化洪水模型，如SWSolver[34]和LISMIN-CS[11]模型。還有采用單指令多數(shù)據(jù)（SIMD）指令集（如SSE 和AVX）的數(shù)據(jù)并行，如CITEEC[43]和NUFSAW2D[44]模型。值得關(guān)注的是，NUFSAW2D 模型采用AVX 數(shù)據(jù)矢量化、MPI 和OpenMP 的混合并行，達(dá)到近1 000 倍的加速比，混合并行充分利用了并行計(jì)算硬件資源，在今后的洪水演進(jìn)模型研發(fā)中應(yīng)重視混合并行技術(shù)。由表4還可以看出：數(shù)據(jù)并行的洪水演進(jìn)模型都采用結(jié)構(gòu)網(wǎng)格，在非結(jié)構(gòu)網(wǎng)格模型的數(shù)據(jù)并行和混合并行研發(fā)方面的加速計(jì)算效率還有待進(jìn)一步研究。

表4 基于其他并行機(jī)制的洪水演進(jìn)模型列表

4 討論

首先，不同并行模式的洪水演進(jìn)模型的研發(fā)技術(shù)難度、開(kāi)發(fā)周期和維護(hù)成本存在較大差異，見(jiàn)表5，同構(gòu)線程并行（OpenMP）和異構(gòu)線程并行（OpenACC）僅需在串行代碼中添加一定的指令即可，編程難度最低，代碼的修改量最小，且不需要第三方程序庫(kù)，采用具備一定編譯功能的編譯器即可編譯和運(yùn)行代碼，其開(kāi)發(fā)周期較短。集群并行模式需要使用第三方的MPI 庫(kù)實(shí)施計(jì)算節(jié)點(diǎn)間數(shù)據(jù)通信，并實(shí)施網(wǎng)格分區(qū)間的數(shù)據(jù)收發(fā)、進(jìn)程管理等，非結(jié)構(gòu)網(wǎng)格區(qū)域分解還需要實(shí)施荷載均衡的第三方程序庫(kù)，編程難度中等，且需要在集群系統(tǒng)上運(yùn)行。GPU異構(gòu)并行和數(shù)據(jù)并行，需要特定的SDK，對(duì)串行代碼的修改量大，甚至要重新采用CUDA語(yǔ)言編程移植，編程難度最大，代碼維護(hù)成本最高。

表5 不同并行化洪水演進(jìn)模型的總結(jié)

其次，采用MPI+OpenMP 和MPI+CUDA 的混合并行技術(shù)，充分挖掘洪水演進(jìn)模型的細(xì)粒度和粗粒度并行潛力[45]，可更進(jìn)一步提高并行計(jì)算效率，見(jiàn)表5。MPI+CUDA 的異構(gòu)并行稱(chēng)為GPU 集群并行，涉及多個(gè)加速器之間的數(shù)據(jù)交換，如何實(shí)現(xiàn)荷載均衡條件下多個(gè)設(shè)備間的數(shù)據(jù)交換速度，是近年來(lái)的研究熱點(diǎn)，實(shí)現(xiàn)方法有兩種：（1）通過(guò)CPU 主機(jī)實(shí)現(xiàn)多GPU 設(shè)備間的數(shù)據(jù)通信，稱(chēng)為三路數(shù)據(jù)通信；（2）多GPU 設(shè)備間直接數(shù)據(jù)通信，稱(chēng)為雙路數(shù)據(jù)通信。后一種方式的數(shù)據(jù)交換速度最快，可使用Nvidia公司開(kāi)發(fā)的GPUDirect技術(shù)來(lái)實(shí)現(xiàn)。

最后，使用類(lèi)似FORTRAN 和C/C++這樣的靜態(tài)編程語(yǔ)言開(kāi)發(fā)洪水演進(jìn)模型，盡管具有很高的計(jì)算效率，但需要大量的代碼編寫(xiě)，且容易出錯(cuò)，開(kāi)發(fā)效率低；盡管FORTRAN或C++語(yǔ)言的面向?qū)ο缶幊炭纱蟠蠼档途幊塘縖46]，但仍無(wú)法實(shí)現(xiàn)像FEM這樣復(fù)雜模型的快速建模。使用動(dòng)態(tài)編程語(yǔ)言（如Python），可快速實(shí)現(xiàn)FDM、FVM和FEM模型的建立[47-48]，結(jié)合Python語(yǔ)言的MPI、CUDA和OpenCL 等庫(kù)的接口程序[49-50]，實(shí)現(xiàn)并行化計(jì)算。另外，由于FEM 算法易于形成抽象化范式，例如FEniCS項(xiàng)目建立了使用FEM，求解弱形式偏微分方程的統(tǒng)一形式語(yǔ)言（UFL）[51]，見(jiàn)表5，應(yīng)用UFL建立了Python語(yǔ)言編程的Firedrake模型[52]，F(xiàn)iredrake模型實(shí)現(xiàn)了具有Python 語(yǔ)言接口的MPI、CUDA 和OP2 程序庫(kù)[53]、具有SSE和AVX指令集優(yōu)化的COFFEE庫(kù)[54]、有限元自動(dòng)組建的FIAT 庫(kù)[55]、用于線性和非線性問(wèn)題求解的PETSc庫(kù)和非結(jié)構(gòu)網(wǎng)格數(shù)據(jù)管理和編號(hào)的DMPlex庫(kù)[56]之間的無(wú)縫耦合，使用FEM 建立了模擬洪水演進(jìn)的Fluids 模型[57]和模擬海嘯波傳播的Thetis模型[58]，包括前處理、計(jì)算主程序和后處理，僅需要300 行Python 代碼。使用DSL 和UFL 可快速實(shí)現(xiàn)對(duì)串行代碼的并行化以及快速實(shí)施新的算法，上述不同并行洪水演進(jìn)模型的優(yōu)缺點(diǎn)及適用場(chǎng)景的具體總結(jié)見(jiàn)表5。

5 結(jié)論及展望

本文對(duì)具有代表性和實(shí)用性的并行化洪水演進(jìn)模型進(jìn)行了綜述，獲得以下幾點(diǎn)結(jié)論：

（1）基于線程并行（OpenMP）、基于進(jìn)程并行（MPI）和基于GPU 異構(gòu)并行（CUDA）的洪水演進(jìn)模型具有很好的實(shí)際應(yīng)用價(jià)值，不同的混合并行模式（如MPI+OpenMP以及MPI+CUDA）將更能充分利用不同層級(jí)的計(jì)算資源，達(dá)到更高的并行加速計(jì)算效率。

（2）并行化非結(jié)構(gòu)網(wǎng)格模式的洪水演進(jìn)模型開(kāi)發(fā)是發(fā)展趨勢(shì)，但存在數(shù)據(jù)結(jié)構(gòu)復(fù)雜、荷載均衡的區(qū)域分解以及數(shù)據(jù)結(jié)構(gòu)對(duì)異構(gòu)并行計(jì)算效率的影響較大等問(wèn)題，可采用相關(guān)算法針對(duì)不同問(wèn)題給出解決方案。

（3）使用特定域語(yǔ)言（DSL）和動(dòng)態(tài)編程語(yǔ)言（如Python）開(kāi)發(fā)并行化的洪水演進(jìn)模型，可快速實(shí)施不同并行模式的洪水模擬加速效率的科學(xué)比較和新的洪水演進(jìn)模擬算法，應(yīng)該是今后的發(fā)展方向。

基于對(duì)并行化洪水演進(jìn)模型的綜述，對(duì)今后并行化洪水演進(jìn)模型研發(fā)提出幾點(diǎn)展望：

（1）非結(jié)構(gòu)網(wǎng)格模式下的并行化洪水演進(jìn)模型研發(fā)應(yīng)成為著力點(diǎn)，應(yīng)重點(diǎn)研究混合并行機(jī)制（如MPI+OpenMP+CUDA）、多GPU 間的數(shù)據(jù)直接交換技術(shù)（如GPUDirect技術(shù)）、解決并行化中自適應(yīng)非結(jié)構(gòu)網(wǎng)格加密的荷載不均衡問(wèn)題（如使用Zoltan庫(kù)）。

（2）針對(duì)GPU架構(gòu)特點(diǎn)，研發(fā)適合于GPU異構(gòu)并行的網(wǎng)格形式，如Block Uniform Quadtree（BUG）網(wǎng)格[59]，BUG網(wǎng)格的數(shù)據(jù)結(jié)構(gòu)具有連續(xù)存儲(chǔ)、結(jié)構(gòu)緊致的特點(diǎn)，特別是結(jié)合BUG 網(wǎng)格和變時(shí)間步長(zhǎng)技術(shù)[60]，可高速求解具有內(nèi)邊界[61]和復(fù)雜河網(wǎng)的洪水演進(jìn)過(guò)程。探討適于GPU加速的網(wǎng)格形式和數(shù)據(jù)結(jié)構(gòu)仍有很大空間。

（3）基于計(jì)算機(jī)圖形庫(kù)OpenGL的實(shí)時(shí)洪水演進(jìn)可視化和交互式計(jì)算是未來(lái)的發(fā)展趨勢(shì)[62]，但如何實(shí)現(xiàn)并行計(jì)算（包括CPU集群并行和GPU集群并行）環(huán)境下的實(shí)時(shí)可視化和交互式計(jì)算，仍是有待解決的問(wèn)題。