山區河流潰壩洪水演進分析

王 雯， 董嘉銳， 楊 杰， 李 鵬， 李占斌， 劉基興， 朱戰齊

(1.西安理工大學 省部共建西北旱區生態水利國家重點實驗室， 陜西 西安 710048; 2.大唐陜西發電有限公司石泉水力發電廠， 陜西 石泉 725200)

1 研究背景

山區洪水具有流速急，水量集中，突發性強的特點。一直以來，山洪災害對世界上諸多山地國家的公共安全和社經濟發展構成巨大威脅[1]；山區又是建造大壩的天然基地，我國各類水庫數量從新中國成立前的1 200多座，增加到98 000多座，已形成初步的防洪減災工程體系，當壩體存在潰壩風險時，突發的潰壩洪水造成的安全威脅不可忽視。

針對潰壩后的洪水演進過程，周興波等[2]基于潰口擴展模型和一維潛水波模型對白格堰塞湖潰決洪水進行了分析和對比，為堰塞湖應急處置提供了參考。馬利平等[3]通過集成 HLLC 近似黎曼求解器的二維水動力模型，對支溝潰壩洪水進行了模擬，說明支溝潰壩洪水對主河道行洪具有一定影響。楊志等[4]通過耦合模型對黑河金盆水庫潰壩過程進行了模擬，結果表明，耦合模型可提高模型計算的準確性。王曉玲等[5]通過三維水動力模型，對東武仕水庫潰壩洪水演進進行了數值模擬，并將傳統的壩區附近模擬范圍擴大到下游城市建筑群整體范圍。這些研究均側重于潰決下泄流量過程[6]，數值計算方法[7-8]以及洪水波演進[9-10]等方面，由于山區河流蜿蜒曲折，有明顯的區間匯流現象，河道寬窄變化，邊界條件復雜，目前對于長距離山區河流潰壩洪水演進的過程及相關問題研究尚顯不足。鑒于此，本文以石泉大壩至喜河大壩段40 km長度河道為模擬對象，通過擬定大壩潰口形狀及潰壩時刻，獲得潰壩下泄流量過程線，建立區段平面二維計算模型，對典型山區型河道潰壩洪水演進過程進行了模擬，為防洪規劃和調洪決策提供依據。

2 數學模型及數值求解方法

2.1 平面二維洪水演進數學模型

本文基于不可壓縮Reynolds值均布的Navier-Stokes方程，滿足Boussinesq假定和靜水壓力假定，建立洪水演進模型。

二維非恒定淺水方程組[11]為：

(1)

(2)

(3)

式中：t為時間，s；h=η+d為總水深，m，其中d為靜止水深； m、η為水面高度，m；u,v分別為x,y方向速度分量，m/s；f為柯氏力系數，f=2ωsinφ，其中ω為自轉角速度、φ為緯度；g為重力加速度，m/s2；ρ為水的密度，kg/m3；sxx、sxy、syy分別為輻射應力通量；Txx、Txy、Tyx、Tyy為黏滯應力項;S為源項矢量;us,vs為源項流速。

(4)

水平黏滯應力Tij的表達式為：

(5)

在空間上對控制方程采用非結構化網格有限體積法離散，可保證連續性方程和動量方程的守恒；采用顯性歐拉法對時間進行離散。

模型基本方程可表示為如下形式：

(6)

式中：U為守恒型物理向量；F為通量向量；S為源項；I為無黏性通量;V為黏性通量。

2.2 模型參數

計算河段劃分采用非結構三角形網格，采用動邊界處理技術[12]處理干濕邊界。河道糙率反映了計算河段的形態變化、邊界條件等因素的綜合影響，計算所采用的河道糙率主要由實測水流資料率定確定，分段、分高程對河道糙率進行了調整，使得5年一遇及20年一遇洪水模擬條件下的水位與實測水位偏差小于1%，調整后河道糙率范圍為0.020～0.035。水平渦黏系數采用Smagorinsky公式計算，模型中采用值為0.28。垂向渦黏系數采用對數定律公式計算。

3 潰壩洪水演進過程模擬

3.1 潰壩方案及流量確定

石泉大壩主壩壩型為混凝土空腹重力壩，最大壩高65 m，壩頂長度353 m，壩基巖石為石英片巖。壩體從左至右編號：1～3壩段為左岸非溢流壩段，4～6壩段為廠房壩段，7～23壩段為泄洪、排沙、溢流壩段，24～29壩段為右岸非溢流壩段。

本次計算采用瞬時潰壩模式，根據大壩安全監測及運行管理人員建議，潰壩范圍擬定為8～23壩段，潰壩缺口底部高程擬定為387 m，潰口寬度為170.5 m，潰口深度控制為29 m；另一方案為非溢流壩段潰決，潰口寬度為52 m，潰口控制深度為39 m，深度到達壩體建基面。潰決時刻為庫區水位達到擬定潰壩水位時刻，此時庫區來流為入庫洪水(設計洪水位410.29 m，校核洪水位413.67 m)。最大下泄流量采用寧利中等[13]提出的解析解求得，疊加入庫洪水及庫容計算潰壩下泄流量過程線，計算結果見圖1。計算工況見表1。

3.2 研究區域及網格劃分

二維數學模型計算河段總長約40 km，起于石泉水庫大壩，上至城關鎮，途經石泉縣城、池河鎮、后柳鎮、上至城關鎮，喜河鎮，老喜河鎮、下至喜河水電站，見圖2。計算區域地形數據根據區段河道帶狀1∶1000地形圖進行構建，局部地區采用實測斷面數據進行修正，計算網格節點數125 006個，網格單元數235 905個。為了盡可能準確的反映區域流場，對河岸及模型邊界處網格進行加密處理，計算區域局部河段網格分布見圖3。

圖1 各工況潰壩下泄流量過程線

表1 各工況計算參數

圖2 計算范圍水系圖及重點防護對象分布圖

圖3 研究區域局部網格結構圖

3.3 結果分析

3.3.1 潰壩洪水演進過程 潰壩洪水演進的關鍵研究內容之一是洪峰的傳播到達時間，精確掌握潰壩洪水到達時間，可有效進行防洪預警，減少洪水損失。工況2和6的特征斷面平均流速隨時間變化圖見圖4和5。由圖4、5可看出，溢流壩段潰壩時，由于潰壩增加流量相對正常泄洪流量較小，下游特征斷面流速增大不顯著；當非溢流壩段潰壩，斷面越靠近大壩，流速峰值越明顯，潰壩發生于洪水過程第45 h，石泉大壩下游G210線石泉漢江大橋處斷面平均流速可達5.28 m/s，但由于石泉縣城段河道順直且河寬較寬，壩體下游5 km處河道自然束窄且形成彎道，為一天然卡口河道，壩體至彎道間形成滯洪區，彎道后十天高速以及安陽鐵路漢江6號橋河段流速增大較為平緩，斷面平均流速峰值相比潰壩下游1 km范圍內削減明顯。由各斷面峰值出現的時間可獲得潰壩傳播過程，潰壩水流在研究區域整體傳播時長為43 min，由于河道支流匯入較多，且大洪水條件下主流河水倒灌入支流，遲滯并削減了洪峰向下游的傳播。

3.3.2 水位特征分析 選取具有代表性的石泉縣城城區段河道進行說明，圖6、7分別為工況1、工況2潰壩發生40 min時洪水水位分布圖。工況1最高水位為378.0 m，工況2最高水位為382.0 m，高水位均位于饒峰河與漢江交匯口左岸區域。主要原因為：石泉大壩泄洪水流頂沖點位于該部位，部分水流動能轉化為勢能，并逐漸向兩側擴散；下泄洪水在饒峰河河口形成逆時針環流，環流頂托饒峰河支流水位，并在支流形成倒灌；石泉縣城段下游河道急劇束窄，并形成彎道，彎道凹岸頂沖點水位明顯抬升，并在彎段下游逐漸降低，石泉縣城段潰壩產生的洪水波波峰在經過下游彎道卡口段后被削減。

圖4工況2特征斷面潰壩洪水流速隨時間變化過程 圖5工況6特征斷面潰壩洪水流速隨時間變化過程

圖6工況1潰壩后40min石泉縣城洪水水位分布圖 圖7工況2潰壩后40min石泉縣城洪水水位分布圖

圖8、9分別為工況2潰壩發生40、43 min計算區域水位分布圖。圖8與9表明，由于石泉大壩下游11 km處池河匯入，池河匯入方向與漢江主河道正交，且交匯口位于漢江彎段左岸凹岸，一方面使得潰壩洪水流量倒灌池河，另一方面頂沖點處山體嚴重頂托下泄洪水，壅高上游段漢江河道水位及池河水位。計算模型整體分析，由于山區河道的彎曲和不規則形態，洪水演進沿程阻力大，受邊界條件影響顯著，使得潰壩洪峰在傳播過程中迅速坦化，對下游喜河電站影響作用逐漸減小，表明河道的蜿蜒能夠很好地降低遠處的受災。

3.3.3 流速分布 以工況2為例進行整體分析，工況2潰壩發生40 min計算區域流速分布見圖10。

由圖10可看出，潰壩發生后大壩下游河道流速迅速增加，但由于河流較長、彎道眾多且主河道兩側存在有支流及支毛溝，潰壩洪峰流量及能量在傳播過程中由于支流倒灌以及交匯口處的環流消耗水量及能量，使洪峰流量過程逐漸坦化，水流流速減小。除潰口段水流流速較大外，其余河段并無出現6 m/s以上斷面平均流速，各工況河段最大流速介于4～6 m/s之間。

以老喜河鎮段河道為例，工況5潰壩發生40 min老喜河鎮段河道洪水流速矢量分布如圖11。此河段洪水由于支流與主流間產生眾多的擴散回流區，支溝水流與主流進行動量交換，削減主流能量，但水流橫向切割沖刷岸坡，可能導致邊灘、河岸崩塌，河道拓寬，應注意重點位置處的防護。局部河段兩岸束窄處岸坡存在有突出的陡峭巖體，造成岸坡段形成較高的流速值，由圖11可見，局部流速可達5 m/s。

圖8工況2潰壩發生40min圖9工況2潰壩發生43min圖10工況2潰壩發生40min

計算區域水位分布圖 計算區域水位分布圖 計算區域流速分布圖

圖11 工況5潰壩發生40 min老喜河鎮段河道洪水流速矢量分布圖

4 結 論

本文通過建立山區河道平面二維水動力學模型，對模型參數進行了率定，模擬了6種工況下的潰壩洪水演進過程，對復雜的地形進行了精細化建模并考慮區間匯流，得出以下結論：

(1)山區河流形態對潰壩水流具有顯著的影響作用，卡口段上游河段易形成局部的滯洪區，通過調蓄洪水，影響下泄流量過程。河流彎段對調整水流結構起著至關重要的作用，彎道的存在，極大削減了潰壩洪水波波峰；

(2)山區河流的支流以及存在的眾多支毛溝，在潰壩洪水發生的條件下，能有效吸納洪峰流量，在交匯處形成環流，削減主流能量，但應注意回流區可能導致的岸坡沖刷失穩；

(3)各潰壩工況斷面平均流速最大值未超過6m/s，但局部范圍內由于河道邊界形狀突變，可能產生局部高流速區，工程防護時應加以注意。