張貴金，羅舸旋子，朱博淵，陳會芳，王振華，周芳軍

(1.長沙理工大學水利工程學院，湖南 長沙 410114；2.長沙理工大學水沙科學與水災(zāi)害防治湖南省重點實驗室，湖南 長沙 410114；3.湖南涔天河工程建設(shè)投資有限責任公司，湖南 永州 425500)

興修水庫、攔河筑壩有利于發(fā)揮防洪、發(fā)電和灌溉等多重效益，但同時也會帶來低概率、高風險的大壩潰決事故，給下游地區(qū)造成重大安全隱患。對潰壩洪水進行研究可為下游風險預(yù)警和避洪轉(zhuǎn)移提供關(guān)鍵信息，對防洪減災(zāi)意義重大[1]。

漫頂潰決是堆石壩潰壩的主要方式之一[2]，其潰決機理為：漫頂洪水在壩體形成潰口，繼而沖蝕壩料，并使得潰口邊坡失穩(wěn)，直至形成穩(wěn)定的潰口斷面[3-7]。針對堆石壩的潰決過程，當前主要利用數(shù)學模型開展研究，主要有參數(shù)模型、精細化模型和簡化模型三類。參數(shù)模型多采用經(jīng)驗公式計算潰壩相關(guān)參數(shù)，而不考慮筑壩材料的沖淤特性[6]；精細化模型對物理機制的考慮最為充分，但計算時間較長，且受潰壩過程中復(fù)雜因素的影響，其結(jié)果的準確性與簡化模型相比尚不明確，應(yīng)用也較少[6-8]；簡化模型既考慮了大壩的潰決機理，基于流深、剪應(yīng)力、堤壩材料等特性模擬潰口和流量的變化特征，結(jié)果可靠，又能保證計算效率[8-13]，因此在堆石壩潰壩數(shù)值模擬中應(yīng)用最為廣泛[14-17]。

洪水演進研究基于對洪水波算法的改進而不斷成熟，其基本算法包括有限體積法、有限差分法和有限元法，據(jù)此開發(fā)了多種數(shù)學模型[18-22]。由于潰壩洪水淹沒范圍廣，多利用平面二維數(shù)學模型開展研究，其中MIKE 21應(yīng)用較為普遍，采用有限體積法，計算效率和結(jié)果精度均較高[23-25]。

本文闡述堆石壩潰壩和洪水演進模型構(gòu)建的一般方法，并以涔天河大壩為例，計算潰壩洪水過程線和模擬洪水對下游的淹沒風險，研究成果可為水庫安全運行管理及編制應(yīng)急預(yù)案提供參考，并對湘江上游瀟水流域防洪安全、水資源科學利用具有重要價值。

1 研究方法

1.1 潰壩模型

1.1.1計算原理

采用DB-IWHR簡化模型對水流運動、壩料沖蝕、邊坡穩(wěn)定性等因素進行綜合考慮，假定潰口形態(tài)，采用寬頂堰過流形式和沖蝕模型模擬壩料的沖刷，并分析潰口的縱向下切和橫向擴展，從而建立壩體潰決過程[12-17]：

a. 潰口斷面流量。潰口出口水流按自由出流考慮，出流流量Q采用寬頂堰過流公式計算：

Q=CB(H-z)3/2

(1)

式中：C為綜合流量系數(shù)，取值范圍為1.3～1.7[24]；B為潰口寬度；H和z分別為庫水位和潰口底部高程。

b. 潰口豎向侵蝕。在高速潰壩水流中，潰壩流量計算結(jié)果受參數(shù)影響極為敏感，指數(shù)形式的侵蝕率模型計算不易收斂，此處以雙曲線模型表征侵蝕率和剪應(yīng)力之間的關(guān)系：

(2)

式中：e為侵蝕率，mm/s；a、b為筑壩材料侵蝕性系數(shù)[24-25]；ν為實際剪應(yīng)力與臨界剪應(yīng)力的差值。

c. 潰口橫向擴展。分析土體邊坡的穩(wěn)定性，采用雙曲線經(jīng)驗?zāi)Ｐ痛_定概化梯形[26]側(cè)面的傾角增量β(圖1)，計算潰口寬度B：

(3)

式中：B0為潰口初始寬度；Δz為潰口底高增量；h為水流經(jīng)過潰口跌落后的水深；β為梯形初始傾角。

圖1 DB-IWHR模型潰口擴展過程概化

1.1.2模型構(gòu)建

a. 潰壩參數(shù)。模型參數(shù)由庫壩特征參數(shù)、水力學參數(shù)和土力學參數(shù)組成。庫壩特征參數(shù)包括庫容、初始庫水位、入庫流量及壩體潰口初始寬度、初始底高和啟動流速，可分別查工程報告和根據(jù)模型內(nèi)部公式計算得到。水力學參數(shù)包括寬頂堰綜合流量系數(shù)和跌落系數(shù)，結(jié)合工程特征、根據(jù)計算和查閱規(guī)范得到[27-28]。土力學參數(shù)[29]包括土料容重、黏聚力、內(nèi)摩擦角和侵蝕系數(shù)，取決于壩體的材料特性，可對壩體堆石材料取樣分析獲得，其中侵蝕性系數(shù)的取值范圍隨顆粒級配變化較大[26-28]，工程中需從偏安全角度取值。

b. 潰壩工況。為分析極端條件下的潰壩事故，充分考慮惡劣的入庫洪水條件、不同潰壩可能性(包括泄洪設(shè)施異常等意外情況)和水庫下游遭受不同強度洪水的損失程度，擬定的入庫洪水量級一般為校核洪水[30]，對不同初始庫水位(如汛限水位、校核洪水位)下的入庫洪水(如校核洪水+多年平均入庫流量、校核洪水+設(shè)計洪水或多場校核洪水疊加)進行組合，最終得到極端潰壩工況，從而對潰壩洪水過程線、庫水位變化、潰口發(fā)展過程等進行計算。

1.2 洪水演進模型

1.2.1計算原理

本文利用MIKE 21開展二維非恒定流水動力模擬計算的基本方程如下：

(4)

(5)

(6)

式中：h為水深；t為時間；x、y為笛卡爾坐標系坐標；u、v為x、y方向的速度分量；S為源匯項；f為科氏力系數(shù)，f=2ωsinφ；ω為地球自轉(zhuǎn)角速度；φ為當?shù)鼐暥?；g為重力加速度，取9.8 m/s2；η為水位；ρ0為標準大氣壓下水的密度；ρ為水的相對密度；pa為當?shù)卮髿鈮?；τsx、τsy為表面風應(yīng)力在x、y方向的分量；τbx、τby為底部切應(yīng)力在x、y方向的分量；sxx、sxy、syx、syy為輻射應(yīng)力分量；Txx、Txy、Tyy為水平黏滯應(yīng)力項；us、vs為源匯項在x、y方向的水流流速。

MIKE 21對干濕邊界自動進行處理，利用3個水深指標(干水深、淹沒深度、濕水深)來判斷網(wǎng)格單元的淹沒狀態(tài)：當單元上水深小于濕水深(0.1 m)、但尚未露灘時，將該網(wǎng)格上的動量通量設(shè)置為0，只考慮質(zhì)量通量；當單元上水深減小至露灘狀態(tài)時，則忽略該單元，直至其被重新淹沒。該算法能有效避免水深過小時的計算失穩(wěn)現(xiàn)象。

1.2.2模型構(gòu)建

a. 模型邊界及初始條件。大壩壩址處為模型上邊界，下邊界通常為潰壩下泄流量衰減至下游河道安全泄量時的斷面位置。若大壩下游一定距離有水庫、湖泊等大型水體，且潰壩洪水對其水位不會造成明顯影響，可將下邊界取在該水體處或其回水末端[30]。以壩址處的最高洪水位和下邊界水位為端點，按線性遞減的原則確定河道沿程水位，并平推至河道兩岸得到模型的陸地邊界(圖2)。上邊界處給定潰壩洪水過程線，下邊界處可給定出流斷面的水位、流速、水位-流量關(guān)系等。在下游支流入?yún)R處給定相應(yīng)的洪水過程線。模型計算時河道初始條件為其沿程初始水位或水深，兩岸地表未過水。

圖2 潰壩洪水演進模型計算邊界確定方法

b. 地形處理。從地理空間數(shù)據(jù)云平臺(URL: http://www.gscloud.cn/)獲取研究區(qū)域的DEM數(shù)據(jù)(30 m×30 m)，考慮到潰壩對下游的淹沒范圍大，該精度能滿足計算要求和保證結(jié)果的可靠性。利用ArcGIS10.2將地形數(shù)據(jù)轉(zhuǎn)換到WGS-84坐標系統(tǒng)，并換算至1985國家高程基準。為較好刻畫研究區(qū)域邊界和地形條件，采用非結(jié)構(gòu)化三角形網(wǎng)格，并使得網(wǎng)格尺寸小于地形測點的空間距離，以保證地形插值的精度。由于河流走向蜿蜒曲折，河道內(nèi)淺灘、深槽并存，且潰壩洪水沿河道行進，需對河道網(wǎng)格進行加密和細化處理。地形插值方法采用克里金方法，該方法對復(fù)雜地形的處理效果較好，近年來得到普遍應(yīng)用[31-34]。最后，將插值得到的連續(xù)地形與地圖影像進行比對，修正失真區(qū)域。

c. 河道及地表糙率。對于實際發(fā)生潰壩的工程，可利用歷史洪水和潰壩洪水對河道及地表糙率進行率定[31]。若工程未發(fā)生潰壩事件，下游河道糙率可采用歷史洪水資料進行率定，而兩岸地表的糙率則需根據(jù)植被分布、房屋阻水等情況查閱相關(guān)手冊[35]、并參考同類工程選取。

2 算 例

涔天河水庫是湘江支流瀟水流域上游以防洪為主的Ⅰ等大(1)型水利樞紐，對下游城鎮(zhèn)具有重要防洪保護任務(wù)，其擴建工程的堆石壩主體工程于2017年竣工，亟須進行潰壩風險評估，故將其作為本研究的算例。

基于潰壩沖蝕機制及寬頂堰過流特性計算潰口初始寬度、底高及啟動流速；寬頂堰綜合流量系數(shù)根據(jù)經(jīng)驗公式[35]計算得到，并與已有研究成果[26]中的取值范圍比對，跌落系數(shù)在已有研究成果[26]的基礎(chǔ)上、結(jié)合本工程特性選??；土料容重、黏聚力和內(nèi)摩擦角經(jīng)現(xiàn)場調(diào)研分析和查閱《涔天河水庫擴建工程正常蓄水位蓄水安全工作報告》[36]得到。最終確定DB-IWHR模型的計算參數(shù)如表1所示。

表1 DB-IWHR模型計算參數(shù)

涔天河水庫防洪庫容加安全超高庫容共計6億m3，校核洪水24 h洪量為5.7億m3[34]，在汛期，按汛限水位運行且考慮泄洪洞正常工作時，水庫可容納4場以上校核洪水疊加，不考慮泄洪設(shè)施和水電站引水洞出流，也能完全蓄滯一場校核洪水。因此，為研究極端條件下的涔天河潰壩事件，假設(shè)泄洪設(shè)施遭遇地質(zhì)災(zāi)害停止工作(如受庫區(qū)霧江滑坡體垮塌堵塞)，在初始水位為汛限水位和校核洪水位的條件下，入庫洪水過程分別設(shè)定為校核洪水疊加多年平均流量和兩場校核洪水疊加的情況。同時，根據(jù)涔天河筑壩材料的侵蝕性[37]，侵蝕性系數(shù)a、b在常用范圍1.2～1.5和0.001～0.01[26-27]內(nèi)分別選取中值(代表侵蝕性小)和最小值(代表侵蝕性大)。最終擬定計算工況如下(兩種工況中校核洪水均采用24 h洪量計算)：

工況1。水庫初始水位為汛限水位(310.50 m)，遭遇一場校核洪水(洪峰流量為9 300 m3/s)后再以多年平均流量(84 m3/s)入庫，庫水位逐漸上升直至漫頂潰壩，取a=1.35，b=0.005。

工況2。水庫初始水位為校核洪水位(320.27 m)，遭遇兩場校核洪水(洪峰流量為9 300 m3/s)入庫，庫水位逐漸上升直至漫頂潰壩，取a=1.20，b=0.001。

圖3顯示，兩種工況均在庫水位達到最大值325 m時開始發(fā)生潰壩，下泄流量逐漸增大，直至分別達到峰值8 275 m3/s和41 528 m3/s，隨后流量逐漸減小至84 m3/s和4 610 m3/s，潰壩總歷時分別為100 h和24.5 h，庫水位最終降低至251 m和248 m，潰口深度達到72 m和75 m。

圖3 涔天河水庫潰壩流量、庫水位及潰口底部高程變化

模型上邊界取涔天河水庫壩址，下邊界取下游雙牌水庫回水末端道縣縣城防洪控制斷面，并根據(jù)河道沿程水位推算陸地邊界，最終得到計算范圍(圖4(a))。計算范圍內(nèi)共有網(wǎng)格單元46 409個，河道區(qū)域進行網(wǎng)格加密，網(wǎng)格面積控制在100～2 000 m2，滿足河道洪水計算的網(wǎng)格要求[29]。壩下游為山地、平原交錯的復(fù)雜地形，采用克里金方法插值并校正后的連續(xù)地形真實度較高(圖4(a))。上邊界給定兩種工況下的潰壩洪水過程線(圖3)，下邊界給定道縣防洪控制斷面的水位流量關(guān)系。河道沿程左岸有兩條較大支流入?yún)R，為充分考慮洪水風險，采用支流歷史最大洪水過程線，并通過試算使得支流入?yún)R洪峰流量與壩址下泄洪峰流量在支流河口處遭遇。計算初始時刻壩下游河道水深為初始水面以下水深，兩岸地表未過水。將干流大水年汛期實測日均流量系列和支流實測同期流量系列疊加，通過調(diào)整河道糙率，使得洪水演進至下邊界附近道縣水文站時流量系列的模擬結(jié)果與實測成果吻合(圖4(b))，最終確定河道糙率為0.035。河道兩岸植被主要為雜草和灌木，房屋在縣城區(qū)域較為密集，經(jīng)查閱《水力計算手冊(第二版)》[35]和對比同類工程[38]最終取地表綜合糙率為0.05。

圖4 洪水演進模型構(gòu)建及糙率率定

洪水向下游演進，兩種工況下的最大淹沒面積分別為232 km2和253 km2(圖5)。洪水演進受地形影響較大，在河道兩岸低洼地帶發(fā)生側(cè)向漫溢，支流河口處，洪水向支流上游發(fā)生倒灌(圖5)。最大流速與最大水深分布規(guī)律基本一致，主河道附近以及受高地環(huán)繞的低洼地區(qū)、特別是計算范圍內(nèi)河道地形較低的中間區(qū)段，淹沒水深和流速均較大，其余地區(qū)水流分散、淹沒水深和流速較小(圖5)。近壩段水深和流速較大，可能沖毀壩下廠房(圖5)。相較于工況1，工況2對下游區(qū)域的淹沒風險顯著增大，表現(xiàn)為淹沒水深及流速整體增大，尤其是計算范圍內(nèi)上游區(qū)段的淹沒水深大幅增加，河道流速及兩個側(cè)向漫溢區(qū)域的流速增大明顯，其余區(qū)域亦有增大(圖5)。潰壩洪水向下游迅速傳播，到達各城鎮(zhèn)斷面(斷面位置見圖4(a))的時間沿程變晚，發(fā)生側(cè)向漫溢后，干流河道的洪峰流量明顯減小，而支流入?yún)R后，干流河道的洪峰流量顯著增大(圖6)。相較于工況1，工況2洪水前鋒到達各斷面的時間大幅提前，各斷面的洪峰流量顯著增大(圖6)。以上計算結(jié)果符合潰壩洪水演進的一般規(guī)律。

3 應(yīng)急預(yù)案對策

水庫預(yù)留庫容結(jié)合配套泄洪設(shè)施一般足以容納一定防洪標準的入庫洪水(包括校核洪水疊加的情況)，因此，漫頂潰壩工況的發(fā)生概率很小，即使發(fā)生，其潰壩過程和洪水演進所帶來的風險要比本文所考慮的極端情況(沿程糙率取大值、潰壩洪峰與支流歷史最大洪水過程線的洪峰遭遇)小得多，對下游的影響將大幅減弱。基于以上成果，提出應(yīng)急預(yù)案如下：

圖6 重點斷面流量過程線

a. 做好建筑物安全管理，確保庫壩安全運行。要確保汛期泄洪洞啟閉設(shè)備的可靠性和電站正常引水發(fā)電，并根據(jù)降雨預(yù)報不斷優(yōu)化防洪調(diào)度方案；加強庫區(qū)水位觀測和大壩安全監(jiān)測，特別要加強對庫區(qū)邊坡的安全監(jiān)測，做好穩(wěn)定和防護措施，防止出現(xiàn)嚴重地質(zhì)災(zāi)害，影響泄洪洞正常運行。

b. 提高下游保護區(qū)防護標準，制定避洪轉(zhuǎn)移預(yù)案。適當提高防護區(qū)重要建筑物和區(qū)段的防洪標準，確保下游地區(qū)信息傳遞渠道暢通，一旦遭遇不利工況，及時采用網(wǎng)絡(luò)平臺或廣播對下游地區(qū)下達提醒、預(yù)警或撤退時間、線路等信息通知。在編制水庫應(yīng)急預(yù)案時，重點關(guān)注最不利工況下影響嚴重的區(qū)域，并根據(jù)洪水演進分析成果，按照潰壩洪水的到達時間進行生命財產(chǎn)轉(zhuǎn)移。根據(jù)SL 483—2017《洪水風險圖編制導則》[30]，轉(zhuǎn)移批次按洪水前鋒到達時間劃分，<12 h、12～24 h和>24 h的區(qū)域分別為第1批、第2批和第3批轉(zhuǎn)移區(qū)，對照算例中兩種極端工況(圖6)的避洪轉(zhuǎn)移預(yù)案如表2所示。

4 結(jié) 語

漫頂潰壩是堆石壩潰決的主要方式之一，基于寬頂堰過流機制，將潰口形態(tài)概化為梯形斷面，考慮不同極端洪水組合下水流對堆石料的侵蝕能力和潰口發(fā)展過程，通過合理選取侵蝕性系數(shù)等參數(shù)，采用DB-IWHR模型對潰壩洪水過程線進行計算。進一步采用MIKE 21模擬潰壩洪水向下游演進的時空規(guī)律，為最大限度考慮洪水淹沒風險，支流選取歷史最大洪水過程線，并使其洪峰與干流洪峰在入?yún)R口遭遇。采用克里金方法對地形進行插值、并修正失真區(qū)域，以保證下游地形導入精度。涔天河水庫算例的結(jié)果表明洪水對下游地區(qū)的淹沒面積、淹沒水深和流速分布及沿程各斷面的洪水傳播特征均符合潰壩洪水演進的一般規(guī)律，且潰壩流量越大、洪水淹沒風險越大。為規(guī)避洪水風險，需做好庫壩建筑物安全管理，優(yōu)化防洪調(diào)度，提高下游保護區(qū)防護標準，并按照潰壩洪水到達時間制定避洪轉(zhuǎn)移預(yù)案。

表2 涔天河潰壩避洪轉(zhuǎn)移預(yù)案