基于一二維耦合技術的潰壩洪水模擬研究

程 濤，何萬成，趙海鏡

(1.中國電建集團北京勘測設計研究院有限公司，北京 100024；2.國網新源控股有限公司，北京 100761)

隨著社會經濟的發展，世界各國對大壩的安全問題都給予了高度的重視。在研究潰壩洪水可能對河道流域下游帶來危險的同時，也在采取積極主動的防范措施[1]。如：美國啟動的國家大壩安全計劃，歐共體則組織了10余個成員國家開展潰壩問題的合作研究，國際大壩委員會則在此基礎上進行整合和協調，發布了ICOLD 111號潰壩洪水分析導則。這對于潰壩洪水的分析計算具有明確的規范和指導作用[2]。

巴丹托魯水電站位于印度尼西亞蘇門答臘島北部的巴丹托魯河上，蘇門答臘島位于亞歐板塊和印度洋板塊邊緣地區，受地震影響頻繁，根據印尼大壩委員會要求，巴丹托魯水電站需開展潰壩洪水風險研究工作。巴丹托魯水電站壩址位于山區與海灘交接處，地形條件特殊，不同下墊面條件對模型精度要求亦不同。目前，國內一二維水動力耦合模型在山區[3- 4]以及平原城市河道[5-7]應用較多，基于山地與海灘交接區域的模擬較少。本文根據巴丹托魯水電站下游地形特點，建立了一二維耦合水動力學模型，針對擬定的潰壩工況進行數值模擬，為后期電站潰壩洪水風險分析和災害損失評估提供技術支撐。

1 模型構建

潰壩數值模擬主要包括潰口模擬和洪水演進模擬兩部分，根據工程特性和地形特點，本次潰口模型選用DamBreak模型，洪水演進模型選用MIKE FLOOD水動力模型，該模型均是在幾十年來大量的工程應用經驗的基礎上持續發展起來的，曾在世界多個國家和地區得到過成功應用[8]。

1.1 模型區域

巴丹托魯水電站采用混合式開發方式，壩型為重力拱壩，庫容1 790萬m3。電站廠房位于壩址下游約17 km處。壩址至出海口河道長度約80 km，其中27 km位于山區，53 km位于海灘地區，考慮到居民分布點主要位于下游海灘地區。本次建模范圍確定為壩址至河道出海口兩岸周邊區域，覆蓋面積約為1 500 km2。

綜合考慮模型適應性和地形特點，山區河道和下游海灘地區分別采用一維和二維模型進行模擬，通過連接工具建立一二維動態耦合水動力學模型。該模型既有一維和二維模型的優點，又避免了采用單一模型時遇到的計算效率、網格精度和準確性方面的問題。模型區域具體分布見圖1。

圖1 模型范圍分布示意

1.2 基本原理

(1)MIKE 11模型。MIKE 11模型控制方程采用一維非恒定流Saint-Venant 方程組，包括垂向積分的物質方程和動量方程，該模型具有穩定性好、計算精度高的特點。

(2)MIKE 21模型。MIKE 21模型控制方程采用平面二維非恒定流N-S方程組，包括水流連續性方程和動量方程。該模型可以計算三角形非結構網格，保證模型邊界的光滑性，同時模型自身可以設置簡單的水工結構物，通過結構物的組合模擬不同區域的復雜地形。

(3)MIKE FLOOD模型。MIKE FLOOD模型采用了一二維動態耦合技術，它可以將一維模型的末端斷面與二維模型里對應一系列地形網格單元連接，在計算過程中，MIKE 11的計算結果作為邊界條件加入MIKE 21網格進行二維計算，同時MIKE 21參與連接的網格平均水位也傳遞回MIKE 11進行一維計算，模擬效果更貼近于實際情況。

1.3 模型建立

一維模型模擬區域采用實測地形數據，根據河道走勢和地形變化沿河道選取了50個斷面作為模型的輸入地形，斷面范圍覆蓋了庫區和壩址下游山區河道，長度約17 km。潰口模型DamBreak在MIKE11模型中以水工建筑物的形式添加至壩址斷面，可結合斷面地形進行潰壩流量計算。

二維模型模擬區域范圍較大，根據雷達地形和SRTM數據，依照網格剖分原則和方法，選用非結構三角網格對模擬區域進行剖分和地形差值。模型整體網格數量約為41萬個，網格邊長范圍10～30 m，同時根據地形數據調整網格尺寸，確保整體網格的平滑度，提高計算穩定性。

2 模型參數

2.1 模擬方案

根據國際上發生的混凝土拱壩潰壩案例數據統計，各機構推薦的潰壩特征參數較為接近，混凝土拱壩潰決形式基本為瞬間全潰[9]。

根據巴丹托魯電站壩型特點和泄流能力分析，發生1 000年一遇洪水情況下，大壩泄流能力依然滿足泄洪要求，重力型拱壩遭遇洪水潰壩風險較小，但巴丹托魯水電站的所在地蘇門答臘島位于亞歐板塊和印度洋板塊邊緣地區，受地震影響頻繁，遭遇地震潰壩的風險相對較大。因此，本次潰壩洪水分析計算工況為地震潰壩工況，計算方案特征參數見表1。

表1 模擬方案特征參數取值范圍

2.2 邊界條件

一維模型上游邊界采用多年平均流量106 m3/s，一維模型下游邊界給定出口斷面的實測水位流量關系。二維模型上下游地區設置為開邊界，兩側水流未到達邊界設置為陸地邊界，二維模型上游流量邊界通過耦合工具MIKE FLOOD給定，下游水位邊界采用恒定海面高程0 m，模型計算時間步長為1 s。

河道斷面糙率根據實測資料推算，取值為0.05～0.07。下游灘地范圍較大，參考不同類型下墊面的糙率取值范圍，綜合實際土地利用情況確定灘地模型糙率，綜合糙率取值為0.07。模型經天然河道水位、流量資料率定，模擬精度較高，可用于壩址下游洪水演進計算。

3 模擬成果

3.1 潰口流量過程

圖2與圖3分別展示了潰口處洪峰流量過程線以及水庫水量變化過程線。由于巴丹托魯水電站庫容較小，整個水庫水體下泄過程僅持續了約28 min，受潰壩引起的洪水波影響，潰口流量和庫水位在潰壩發生后出現短期震蕩現象，其中潰口最大洪峰流量為28 300 m3/s，隨著庫水位逐漸下降，潰口洪峰流量也逐漸減小。

圖2 潰口流量過程

圖3 水庫水量變化過程

3.2 洪水演進過程

由于村莊和農田主要分布在下游灘地地區，因此圖4統計了潰壩洪水進行過程中下游地區最大淹沒水深，成果顯示洪水在上游山區段水深變化較為顯著，瞬時局部最大水深接近20 m，隨著洪水向灘地行進，洪水大量漫出主河槽，對河道兩岸的部分村莊造成了影響，但由于灘地面積較大，淹沒水深遠小于山區河道，灘地淹沒水深在0～3 m范圍內。

圖4 最大淹沒水深分布

為了更直觀地展示潰壩洪水沿程的變化過程，從壩址到下游灘地選取了5個典型斷面進行了水力要素統計，典型斷面水力要素統計成果見表2。

表2 典型斷面水力要素統計

由圖4及表2可以看出，潰壩發生初期，壩址斷面洪峰流量和行進流速最大，隨著向下游演進，洪峰流量逐漸坦化，在山區河道行進約48 min即到達位于下游的巴丹托魯村斷面。巴丹托魯村下游地區地形較為平坦，洪水經過巴丹托魯村后大量漫灘，由于水庫庫容較小，下泄水量有限，潰壩水體釋放后下游河槽流量急劇減小，同時流速降低，下游海灘地區河道受潰壩洪水影響較小。

3.3 敏感區域淹沒分析

巴丹托魯電站下游地區人煙稀少，根據潰壩洪水淹沒水深計算結果，主要影響區域為位于山區和海灘交匯處的巴丹托魯村，該村沿河分布有大量房屋和公路、橋梁，下游海灘地區沿岸村莊受洪水影響程度相對較輕。巴丹托魯村斷面最大洪峰流量5 880 m3/s，最高洪水位73.2 m，洪峰持續時間20 min，參考實測村子分布高程數據，確定出各設施遭遇潰壩洪水時的影響范圍，受影響設施統計結果見表3。

表3 敏感區域淹沒設施統計

4 結 語

本文以印度尼西亞巴丹托魯水電站為例，對壩址下游地形條件復雜的水電站潰壩洪水數值模擬技術進行了探討，總結出如下幾點體會，供同類工程參考。

(1)大區域、復雜地形條件下潰壩洪水數值模擬需同時考慮計算精度與計算效率問題。本文建立的適用于不規則計算域和地形的潰壩洪水一二維耦合模型，引入了最新的動態耦合技術，避免了采用單一模型時遇到的計算效率、網格精度和準確性方面的問題。

(2)潰壩洪水在不同區域衰減速度差異較大，潰壩發生初期，山區河道洪峰流量和行進流速最大，山區河段洪峰坦化幅度相對較小，隨著向下游演進，潰壩水體進入海灘區域后，洪水發生漫灘，下游河槽流量急劇減小，同時流速降低，下游海灘地區受潰壩洪水影響較小。

(3)潰壩時的水流運動過程常伴隨有泥沙輸移等其他物理過程，泥沙沖淤與水流運動會產生相互作用，后期可以水動力模型為基礎，建立適用于潰壩模擬的水沙數學模型，實現潰壩洪水多過程耦合模擬。