基于Google Earth的分蓄洪區(qū)水沙模擬與演示系統(tǒng)研究

，

(長江科學院 河流研究所，武漢 430010)

自然潰堤或人工分洪情況下，洪水在洪泛區(qū)內(nèi)的演進將給洪泛區(qū)內(nèi)人民生命財產(chǎn)安全帶來巨大威脅，因此模擬洪水在洪泛區(qū)內(nèi)的演進過程，將給防洪決策和防洪搶險提供技術(shù)支撐。目前，許多學者對分蓄洪區(qū)洪水演進開展了數(shù)值模擬研究[1-3]，然而分蓄洪區(qū)地形獲取的不便和模擬成果顯示不直觀等問題，仍然是制約數(shù)學模型發(fā)展的重要因素。近10多年來，地理信息技術(shù)得到快速發(fā)展，形成了許多的地理信息平臺。因此，如何將洪水演進模擬結(jié)果顯示在地理信息平臺之上，是一件非常有意義的事情。

Google Earth軟件(以下簡稱GE)便是其中應用最為廣泛的通用的地理信息平臺之一，它把衛(wèi)星照片或航空照相和GIS信息布置在一個三維的地球模型上，而且免費對大眾開放。本文基于Google Earth軟件平臺，建立了分蓄洪區(qū)水沙演進數(shù)學模型，能快速完成模型數(shù)據(jù)前處理和模擬計算，并實現(xiàn)了計算結(jié)果在三維虛擬場景中的動態(tài)演示。該模型系統(tǒng)適用于緊急條件下的快速計算和直觀展示，計算成果可為分蓄洪區(qū)防洪決策與防洪搶險提供技術(shù)參考。

1 基于非結(jié)構(gòu)網(wǎng)格的數(shù)值求解方法

分蓄洪區(qū)一般不規(guī)則，為此采用非結(jié)構(gòu)網(wǎng)格[4]進行分蓄洪區(qū)的網(wǎng)格劃分是非常合適的。非結(jié)構(gòu)網(wǎng)格可以是三角形網(wǎng)格、四邊形網(wǎng)格或三角形與四邊形混合網(wǎng)格。

模型基本方程采用守恒形式的水深平均平面二維水流基本方程[5]。采用基于同位網(wǎng)格的有限體積算法進行方程的離散[6]，流速矢量基于直角坐標系，節(jié)點布置采用單元中心法；采用SIMPLEC算法進行離散方程的求解，該算法的優(yōu)點在于即便采用較粗的網(wǎng)格也能保持水量和動量的守恒。

2 基于Google Earth的模擬演示系統(tǒng)

本文基于Google Earth軟件進行二次開發(fā)，將數(shù)值模擬成果在由航拍地貌和三維建筑物構(gòu)成的虛擬地球上進行顯示，有助于更好地分析計算成果的合理性，同時還可利用Google Earth軟件縮放、平移、鳥瞰等功能，更直觀地查看計算成果。

Google Earth軟件的二次開發(fā)有2種方式，一種是基于ComAPI；一種是基于KML(全稱為Keyhole Markup Language)?；?ComAPI的開發(fā)方式主要用來控制GE的視角,實現(xiàn)動畫效果，而基于KML的開發(fā)方式主要用來生成地理要素(如點、線、面、多邊形、多面體以及模型等)，實現(xiàn)數(shù)據(jù)的動態(tài)更新等。在實際應用中，需要結(jié)合這2種開發(fā)方式，才能得到比較理想的效果。

筆者采用VB和Fortran語言混合編程方法，對Google Earth軟件進行二次開發(fā)，研制了“基于虛擬地球的水沙水質(zhì)數(shù)值模擬與演示系統(tǒng)V2.0”，系統(tǒng)由前處理、模型庫、后處理及數(shù)據(jù)庫等部分構(gòu)成，模型可實現(xiàn)干河床的水沙演進計算。該系統(tǒng)不僅實現(xiàn)了界面操作的可視化，而且實現(xiàn)了計算結(jié)果在虛擬地球上的動態(tài)演示，系統(tǒng)主界面見圖1。

圖1 可視化界面圖

2.1 基于ComAPI的開發(fā)

GE 公布的API目前是1.0a版，數(shù)量比較少，所能實現(xiàn)的功能也不多，在實際開發(fā)中主要用來加載KML數(shù)據(jù)、控制視角和實現(xiàn)動畫播放。加載KML數(shù)據(jù)命令有：

(1) IApplicationGE. OpenKmlFile(kmlFilePath)，加載指定的KML文件到GE中；

(2) IApplicationGE.LoadKmlData(kmlData)，加載KML文本流到GE中。

2.2 基于KML文件的開發(fā)

2.2.1 坐標轉(zhuǎn)換

目前我國的河道地形圖大地基準面一般采用54坐標系或西安80坐標系，而Google Earth軟件是基于WGS84坐標系統(tǒng)，因此在生成KML文件之前需對計算結(jié)果文件進行坐標轉(zhuǎn)換。上述3類坐標系的轉(zhuǎn)換在相關(guān)教科書中有詳細的說明，這里不再贅述。

2.2.2 流速矢量場

計算成果顯示包括流速矢量場和水位、水深、流速等物理量的標量場，同時還需文字標注、圖例。為進行流場的動態(tài)演示，需進行每幀流場的時間控制。

圖2 流速箭頭示意圖

計算域內(nèi)每個點的流速矢量可用箭頭表示，箭頭長短表示流速大小，箭頭方向代表流速方向。每個箭頭可由6個有序點形成的折線表示，如圖2所示。

在KML文件中，每個點的流速矢量可用LineString標簽來定義，代碼如下：

112.138 084 0, 30.284 331 1, 0

112.140 711 6, 30.283 183 1, 0

......

其中字段表示定義標簽；定義該箭頭名稱；定義流速箭頭顏色，需根據(jù)流速的大小預先定義；繪制折線；下面定義6個點的坐標，注意這里采用的是經(jīng)轉(zhuǎn)換后的經(jīng)緯度坐標。

2.2.3 等值線云圖

等值線的繪制與流速箭頭的繪制類似，只需繪制每條折線就行，云圖的繪制需將封閉折線構(gòu)成的區(qū)域進行顏色填充，可用Polygon標簽來定義，代碼格式如下：

112.141 431 5, 30.283 587 2, 0

112.140 788 8, 30.281 146 2, 0

112.144 288 4, 30.282 297 2, 0

112.141 431 5, 30.283 587 2, 0

112.141 431 5, 30.283 587 2, 0

2.2.4 文字標注

流速值、水深值或其他文本的標注可采用Google Earth中的Point地標標簽來定義，代碼如下：

2.2.5 流場圖例

完整的流場顯示需定義圖例，在Google Earth中圖例由ScreenOverlay標簽控制，在動態(tài)演示過程中保持不變，代碼如下：

2.2.6 時間控制

為實現(xiàn)流場的動態(tài)演示，需定義每幀流場的起止時間，采用TimeSpan標簽定義，代碼如下：

3 在荊江分蓄洪區(qū)的應用

荊江分蓄洪區(qū)位于長江南岸，是長江中游防洪系統(tǒng)的主要組成部分，分洪區(qū)面積921 km2，有效蓄洪容量54億m3。荊江分洪工程包括進洪閘、節(jié)制閘和分洪區(qū)圍堤。進洪閘位于太平口東岸，設(shè)計進洪流量8 000 m3/s。荊江分蓄洪區(qū)的功能在于分泄荊江河段的超額洪水，使洪水能夠安全通過，以保護荊江兩岸廣大農(nóng)田和人民生命財產(chǎn)安全。

1954年7月22日至8月1日荊江分洪區(qū)進行了3次分洪運用，分洪總量為122.6億m3。采用建立的模型系統(tǒng)，對1954年荊江分蓄洪區(qū)分洪過程進行了模擬計算。

3.1 計算條件

計算范圍為整個荊江分蓄洪區(qū)，以北閘為進流邊界，南閘為出流邊界。網(wǎng)格采用非結(jié)構(gòu)網(wǎng)格，見圖3。計算網(wǎng)格地形根據(jù)航拍的DEM地形數(shù)據(jù)插值得到。

圖3 分洪區(qū)計算網(wǎng)格布置

將荊江分蓄洪區(qū)分洪過程概化為潰壩過程，潰壩計算關(guān)鍵在于如何進行干河床與動邊界的處理，筆者提出了“干濕雙最小水深法[6]”，能較好反映潰壩水流在干河床上的演進過程。

由于荊江分蓄洪區(qū)內(nèi)缺乏實測的洪水資料，因此，糙率參考相關(guān)資料確定[6]：樹林0.070，旱地0.065，水田0.050，水面0.025。如果某網(wǎng)格內(nèi)含有多種地形，則按照各種地形糙率的加權(quán)平均值確定該網(wǎng)格的糙率。

模擬計算從1954年7月22日8:00開始，考慮為全潰，直至蓄滿為止。進口含沙量假定為零，床沙中值粒徑取0.2 mm。

3.2 計算成果分析

3.2.1 洪水演進計算成果

計算結(jié)果表明，分洪閘下泄流量先緩慢增大，最大可達8 000 m3/s，后快速減小，直至蓄滿所需時間為11 d左右。

圖4給出不同時刻洪水在分蓄洪區(qū)的演進情況(圖中箭頭表示流速矢量大小和方向，藍色云圖表示水深，藍色越渾則水深越大)。從圖4中可見，系統(tǒng)較直觀地反映了洪水在航拍地貌上的演進過程，在潰口和分洪區(qū)縮窄段，流速明顯增大，最大流速為1.2 m/s，圖中彩色云圖能直觀反映不同時刻洪水淹沒范圍及淹沒水深。

圖4 不同時刻洪水演進流場圖

3.2.2 泥沙輸移計算成果

圖5(a)給出了分洪67 h時分蓄洪區(qū)含沙量分布情況(圖中顏色越紅表示含沙量越大)。從圖4中可見，由于進口未給定含沙量，泥沙主要來自于床面的沖刷，因此流速越大部位(如圖5(a)中深藍色區(qū)域)，含沙量也會越大，最大含沙量可達2 kg/m3左右，位于分洪口門下游附近；而在分洪區(qū)下游，由于流速不斷減小，含沙量也隨之減小。

圖5 分蓄洪區(qū)含沙量和沖淤厚度分布(T=67 h)

圖5(b)給出了分洪67 h后分蓄洪區(qū)沖淤厚度分布情況(圖中紅色表示淤積，藍色表示沖刷，顏色越深則沖淤幅度越大)。從圖中可見，一般大流速部位，床面易發(fā)生沖刷，最大沖深為1.7 m，其余大部分區(qū)域發(fā)生淤積，最大淤厚約1.5 m。

4 結(jié) 論

(1) 建立了基于地理信息技術(shù)的分蓄洪區(qū)平面二維水沙演進數(shù)學模型。模型能直接在Google Earth軟件環(huán)境中確定計算域邊界，采用無結(jié)構(gòu)的三角形、四邊形網(wǎng)格進行計算域的網(wǎng)格劃分，根據(jù)航拍的DEM地形插值得到計算網(wǎng)格地形高程，提高了前處理的效率。

(2) 對Google Earth軟件進行二次開發(fā)，研制了計算結(jié)果后處理模塊，實現(xiàn)了分蓄洪區(qū)水沙演進過程在三維航拍地貌中的動態(tài)演示，并能與Google Earth同步實現(xiàn)縮放、鳥瞰、漫游等多種查看功能。

(3) 該模型系統(tǒng)在長江荊江分洪區(qū)得到成功運用，復演了1954年分洪后的水沙演進過程，并在Google Earth地貌中進行動態(tài)演示，可直觀查看水深、流速、洪水淹沒范圍、沖深等要素，可為洪水風險評估和應急預案的制定提供參考。

參考文獻：

[1] MACCHIONE F, MORELLI M. Practical Aspects in Comparing Shock-Capturing Schemes for Dam-break Problems [J]. Journal of Hydraulic Engineering, 2002, 129(3)：187-195.

[2] 余明輝，張小峰. 平面二維潰堤水流泥沙數(shù)值模擬[J]. 水科學進展，2001， 12(3)：286-290.(YU Ming-hui, ZHANG Xiao-feng. Horizontal 2-D Uneven Sediment Mathematical Model in Dike Burst[J]. Advances in Water Science, 2001, 12(3): 286-290.(in Chinese))

[3] 王志剛，汪德爟，賴錫軍, 等. 下游為干河床的潰壩水流數(shù)值模擬[J]. 水利水運工程學報，2003，(2):18-23.(WANG Zhi-gang, WANG De-huan, LAI Xi-jun,etal. Numerical Simulation of Flow on Dry Bed for Dam Break[J]. Hydro-science and Engineering, 2003, (2):18-23.(in Chinese))

[4] 張細兵.河道有限元網(wǎng)格自動剖分方法研究[J]. 長江科學院院報，2002，19(3):19-21.(ZHANG Xi-bing. Auto-generating Method of Triangulation Grids for River Reach[J]. Journal of Yangtze River Scientific Research Institute, 2002,19(3):19-21.(in Chinese))

[5] 張細兵，范北林. 潰壩洪水演進平面二維數(shù)學模型初步研究[J]. 長江科學院院報, 2006，23(6):19-21.(ZHANG Xi-bing, FAN Bei-lin. Numerical Simulation for Two-Dimensional Dam-break Flood Flow[J]. Journal of Yangtze River Scientific Research Institute, 2006, 23(6): 19-21.(in Chinese))

[6] 張細兵，歐治華，崔占峰，等. 基于非結(jié)構(gòu)網(wǎng)格的分蓄洪區(qū)水沙演進數(shù)學模型研究[J].長江科學院院報，2011，28(4):75-79.(ZHANG Xi-bing, OU Zhi-hua, CUI Zhan-feng,etal. Unstructured Grid Model of Flow and Sediment Evolution in Flood Diversion and Storage Area[J]. Journal of Yangtze River Scientific Research Institute, 2011, 28(4):75-79.(in Chinese))