寧 聰,傅志敏,王志剛

(1. 河海大學(xué) 水文水資源學(xué)院，江蘇南京 210098； 2. 黃河水利水電開發(fā)總公司，河南濟源 459017)

水庫大壩作為水利樞紐的重要組成部分，在防洪，發(fā)電，灌溉等功能中發(fā)揮著重要作用。水庫大壩失事將危及下游地區(qū)人民的生命財產(chǎn)安全與社會穩(wěn)定[1- 2]。為了對水庫潰壩進行有效的風(fēng)險防范，需在潰后影響區(qū)域進行洪水模擬計算。目前主流的洪水數(shù)值模擬軟件包括丹麥的MIKE11/21，荷蘭的Delft3D和美國的HEC-RAS。其中HEC-RAS為美國陸軍工程兵軍團水文中心開發(fā)的免費河道水力計算軟件，前人利用HEC-RAS在洪水數(shù)值模擬方面已進行大量工作。周毅[3]利用HEC-RAS和GIS平臺模擬了疏勒河地區(qū)2000年一遇洪水在下游區(qū)域的演進情況；賀娟等[4]利用HEC-RAS對長河壩水電站進行了潰壩洪水模擬；吳博等[5]利用HEC-RAS和GIS平臺對小東川河流域的山洪淹沒范圍做出了較為準確的預(yù)測；孫銳嬌等[6]利用HEC-RAS模擬多種工況下某水庫潰壩洪水演進。但上述研究均建立在HEC-RAS一維水動力學(xué)模型上，而 HEC-RAS近期增加二維水動力學(xué)模型，其模型精度高，模型構(gòu)建要求低，適用于洪水在山區(qū)河道沖刷和平原地區(qū)泛濫等多種場景。本文運用HEC-RAS的二維模型模擬不同工況下小井溝面板壩漫頂潰壩后洪水演進，其結(jié)果對水庫防災(zāi)減災(zāi)以及風(fēng)險防范工作具有重要意義。

1 研究方法

HEC-RAS二維水動力學(xué)模型的原理是Navier-Stokes方程的二維簡化形式——淺水方程。其假定水深尺度遠小于另外兩個平面尺度，計算式如下：

連續(xù)方程：

(1)

動量方程：

(2)

擴散波格式的動量方程，其與連續(xù)方程組合的計算速度較完全的淺水方程的快且累計誤差小，適用于河床坡降大的河流[7]。計算式如下：

(3)

式中：H為水面高程(m)；h為水深(m)；V為流速(m/s)；R為水力半徑(m)；q為旁側(cè)入流(m2/s)；g為重力加速度(m/s2)；υt為水平方向運動黏度(m2/s)；cf為河床底部糙率；f為科里奧利系數(shù)；k為垂直方向單位矢量；n為糙率。

HEC-RAS二維水動力學(xué)模型的數(shù)值計算混合了有限體積法和有限差分法。計算網(wǎng)格采用非結(jié)構(gòu)化網(wǎng)格，非邊界網(wǎng)格為正方形，邊界網(wǎng)格為不規(guī)則多邊形，每個網(wǎng)格的邊類似河道斷面，均能提取所在的地形。因而在較低的網(wǎng)格密度下，仍可以提取足夠的地形細節(jié)，保證模型精度。

HEC-RAS在潰壩計算上采用堰流方程計算潰口流量過程。但其堰流系數(shù)應(yīng)取1.76～1.98[7]。計算式如下：

Q=CLH13/2

(4)

式中：Q為潰口流量(m/s)；L為堰長(m)；H1為堰上水頭(m)；C為堰流系數(shù)。

HEC-RAS潰口預(yù)測模塊集成了一系列以往學(xué)者根據(jù)歷史潰壩數(shù)據(jù)推導(dǎo)出的潰口尺寸回歸方程，能夠根據(jù)大壩相關(guān)參數(shù)預(yù)測潰口。本次模擬采用其中能夠預(yù)測面板壩潰口的徐-張方程估計潰口最大尺寸和成形時間[8]。計算式如下：

(5)

(6)

式中：Bt為潰口頂寬(m)；Vw為潰壩時庫容(m3)；B2為綜合系數(shù)，取0.325；Tf為潰口成形時間(h)；Tr為單位時間(1 h)；B5為綜合系數(shù)，取-1.817；hb為最大潰口深(m)；hd為壩高(m)；hr為參考高度，取15 m。

2 模型建構(gòu)

2.1 工程概況

越溪河屬岷江左岸一級支流, 流域大致呈南北向的狹長形，東接沱江右岸支流釜溪河，西鄰岷江。流域地勢北高南低。由于河流穿行于深丘地區(qū)，河谷狹窄，河道成V型，灘多水急，河道平均比降達7‰左右。禮佳以上至正江有約9 km長的峽谷段，岸陡水急，比降更大。越溪河上游植被良好，河源及山頂多成片幼林及灌木叢。

小井溝水庫是一座大(2)型水庫，防洪標準為100年一遇設(shè)計，2 000年一遇校核。排洪工程建筑物級別5級，設(shè)計防洪標準為10年。小井溝水庫攔河大壩位于越溪河上游，禮佳場上游小井溝峽谷河段上，地理位置在東經(jīng)104°10′，北緯29°24′，壩址控制流域面積587 km2，占全流域面積的22%，河長92 km，設(shè)計流域為西北東南向的狹長形。東高西低，東為榮、威高地西緣，分水嶺高程為700～800 m，西邊分水嶺為500～600 m[9]。

2.2 模型構(gòu)建

2.2.1地理建模 地理建模的地形數(shù)據(jù)來源于小井溝流域10 m分辨率DEM。建模中的各種因素概化為3個模型要素，分別為上游水庫、攔河壩和下游控制流域。上游水庫在HEC-RAS中設(shè)為線性水庫，輸入小井溝水庫庫容曲線[9]，如圖1所示。下游控制流域設(shè)為單個二維網(wǎng)格區(qū)域[10-11]，網(wǎng)格密度設(shè)為5 m×5 m。河床主槽糙率取0.027，漫灘糙率取0.032。在設(shè)置網(wǎng)格邊界時，保證計算區(qū)域囊括洪水可能淹沒區(qū)域，考慮越溪河中下游多為平原和丘陵地區(qū)的特點，最終地理建模效果如圖2所示。

圖1 小井溝水庫水位庫容關(guān)系Fig.1Relationship between water level and storage capacity

圖2 地理建模Fig.2Geometric data

HEC-RAS能夠模擬水工建筑物的過流能力，在小井溝攔河壩壩體設(shè)置相應(yīng)的泄洪洞和溢洪道，其上設(shè)置閘門，大壩相關(guān)參數(shù)如下:壩頂高程431.60 m,壩底高程344.00 m,壩頂長263.00 m,壩頂寬8.00 m,壩底寬350.94 m,總庫容為1.66億m3,死庫容為0.35億m3。

2.2.2潰壩參數(shù) 潰壩參數(shù)主要是潰口最大尺寸和成形時間。小井溝攔河壩作為一座面板壩，其潰口發(fā)展并非土石壩潰口的線性增長模式，而是伴隨著洪水的淘蝕造成的上游面板階段性折斷呈現(xiàn)出分段逐級潰決的特點。上游面板從頂部開始，其后部由于洪水沖刷而被逐漸掏空，隨后該部分面板折斷。面板如此逐級折斷且折斷速度加快，直至發(fā)展成線性增長。根據(jù)面板壩分段逐級潰決特點在HEC-RAS中修正潰口發(fā)展曲線，使之契合面板壩的潰壩過程[12-14]，如圖3所示。利用HEC-RAS自有的潰口預(yù)測模塊，預(yù)測潰口最大尺寸如圖4所示，具體參數(shù)如表1所示。

圖3 潰口發(fā)展曲線Fig.3Breach progression curve

圖4 最大潰口尺寸Fig.4Maximum breach dimensions

中心坐標/m底寬/m底部高程/m左邊坡坡度右邊坡坡度堰流系數(shù)成形時間/h模式起潰高程/m120694000.50.52.60.47漫頂431.6

圖5 設(shè)計洪水過程線Fig.5Design flood hydrograph

潰壩事故中往往存在各種突發(fā)狀況，如泄洪建筑物可能無法工作而使?jié)挝：υ龃蟆１敬文M中的閘門分別取完全關(guān)閉和最大開度兩種工況。閘門完全關(guān)閉模擬潰壩洪水到來時，閘門因故障不能開啟導(dǎo)致水庫無法泄洪的場景；閘門最大開度則模擬潰壩洪水到來時，泄水建筑物最大能力泄洪的場景。對兩種工況下潰壩洪水演進結(jié)果差異進行對比，考量水庫泄洪對潰壩危害的影響，計算工況為：工況1，逐漸潰，潰口底寬69 m，底部高程400.00 m，起潰水位431.60 m，閘門完全關(guān)閉；工況2，閘門最大開度，其余同工況1。

2.2.3邊界條件 上游邊界條件為水庫入流，取小井溝水庫2 000年一遇設(shè)計洪水過程線[9]，如圖5所示。下游邊界條件取越溪河下游天然河床比降0.42‰。水庫初始水位設(shè)為校核洪水位430.75 m，模擬入庫洪水到來后水庫水位壅高過壩頂造成漫頂潰壩的場景。此外，鑒于越溪河河床平均比降較大，本次計算選用擴散波格式的淺水方程加快計算速度。

2.3 結(jié)果分析

洪水水深、流速和滯留時間是衡量潰壩洪水影響的3個主要指標。其中依據(jù)洪水水深差異可以劃定不同的風(fēng)險等級[15]，洪水風(fēng)險分級見表2[15]。此外，根據(jù)洪水流速大小同樣能夠給出相應(yīng)的預(yù)期損害[16]，如表3所示。

表2 洪水風(fēng)險等級Tab.2 Flood hazard classification等級水深/m風(fēng)險等級H1≤0.5很低H2>0.5～1.0底H3>1.0～2.0中等H4>2.0～5.0高H5>5.0極高表3 洪水流速等級Tab.3 Flood velocity classification分類流速/(m·s-1)損害等級1≤0.2極低2>0.2～0.5低3>0.5～1.0中4>1.0高

HEC-RAS能夠輸出潰壩洪水下游演進各時刻的水深和流速分布圖、模擬時段內(nèi)的洪水滯留時間、抵達時間和退水時間分布圖等。計算完成后導(dǎo)出最大水深和流速分布圖(將各淹沒處的最大水深或流速疊加)，以及閾值1.5 m的洪水滯留時間分布圖，在GIS軟件中加以處理[17]，獲得潰壩洪水風(fēng)險分布、最大流速分布和潰壩洪水滯留時間分布(見圖6～8)。

由圖6可知，越溪河上游為深V河谷，潰壩洪水沒有向兩岸明顯的泛濫。中下游為平原丘陵地區(qū)，河道蜿蜒曲折，兩岸地勢平緩，洪水壅出河道淹沒周邊地區(qū)。而人員聚居區(qū)多集中在中下游河道兩岸地勢平坦地帶，在這一區(qū)域，工況1洪水淹沒范圍較工況2大，洪水風(fēng)險等級約在H5，也比工況2的風(fēng)險等級(約H4)高一級。說明水庫全力泄洪能夠有效削弱潰壩洪水在潰后影響區(qū)域造成的洪水風(fēng)險。

圖6 潰壩洪水風(fēng)險分布

由圖7可知，上游越溪河河道臨近水庫，洪水演進快，主槽中流速在10 m/s以上。中下游越溪河河道主槽流速衰減到5～10 m/s，河岸漫灘在2～5 m/s。在人員聚居的兩岸淹沒區(qū)域，兩種工況的最大流速均小于1 m/s，影響程度為中等?？梢妰煞N工況下潰壩洪水流速不是對下游地區(qū)造成危害的主要因素。

圖7 潰壩洪水最大流速分布

洪水水深超過1.5 m預(yù)期會造成嚴重的生命財產(chǎn)損害[16]。由圖8可知，兩種工況下河道兩岸淹沒區(qū)域洪水滯留時間均在24 h以下，且自河岸向外逐步衰減。緊臨河岸的淹沒區(qū)域洪水滯留時間在12～24 h，在稍遠區(qū)域滯留時間下降到了2～8 h，局部地區(qū)更短。工況2兩岸淹沒范圍小且緊臨河岸，在這一區(qū)域，兩種工況下洪水滯留時間并無顯著差異。

圖8 潰壩洪水滯留時間(閾值1.5 m)

綜上所述，水庫全力泄洪主要能夠削減洪水在中下游河道兩岸淹沒區(qū)域的水深和淹沒范圍，從而降低潰后影響區(qū)域洪水危害，同時為下游地區(qū)人員轉(zhuǎn)移疏散爭取時間。此外，結(jié)果也表明了水庫日常檢查的重要性，尤其是帶閘門的泄洪建筑物，必須充分保證工作狀態(tài)和泄洪能力，為防范事故發(fā)生提供安全保障。

3 結(jié) 語

本文采用HEC-RAS模擬不同工況下小井溝攔河面板壩漫頂潰壩后洪水下游演進。結(jié)果表明，HEC-RAS二維模型適用于模擬潰壩洪水向河道兩岸平原丘陵地區(qū)泛濫的，并能夠?qū)⒛M結(jié)果輸出為可視化圖層。若將結(jié)果圖層與社會行政，經(jīng)濟區(qū)劃圖疊加制作信息更豐富的洪水風(fēng)險圖，可為后續(xù)防洪應(yīng)急預(yù)案的編制和財產(chǎn)損失的估算提供有力支持；其次，通過模擬了閘門不同開度下潰壩洪水的演進過程，表明了潰壩事故發(fā)生時水庫全力泄洪對潰壩洪水危害的削弱作用。