缺水文資料地區(qū)小流域暴雨洪水數(shù)值模擬研究

付 剛 劉 威 何思明 羅星文 王任國

（1.中鐵二院工程集團有限責任公司，成都 610031；2.中國科學(xué)院水利部成都山地災(zāi)害與環(huán)境研究所，成都 610041）

老撾地處東南亞，中老鐵路磨丁至萬象段北端連接與我國玉溪至磨憨鐵路，南端與規(guī)劃中的泰國曼谷至廊開鐵路相接。該鐵路起于中老邊境口岸磨丁，經(jīng)老撾北部的南塔省、烏多姆塞省、瑯勃拉邦省、萬象省，到達首都萬象。鐵路跨越水系主要為湄公河及其支流，鐵路所經(jīng)地區(qū)河流眾多。在全長424 km 的鐵路中，除去196 km 隧道外，剩余228 km 共設(shè)橋梁165 座，擁有排洪功能的涵洞更是多達510 座。受老撾獨特的地理和氣候條件影響，海洋和大陸對老撾氣候產(chǎn)生強烈影響，使其成為亞洲季風區(qū)之一。年內(nèi)降雨變化大，有明顯的旱季和雨季之分，雨季持續(xù)時間長。鐵路所經(jīng)的北部山區(qū)（磨丁至瑯勃拉邦）年平均降雨量1 000～2 000 mm，而萬榮年平均降雨量可達3 875 mm。加之河道位于中南半島山地，河道縱坡陡峻，導(dǎo)致河流洪峰流量大，陡漲陡落歷時短促，給當?shù)鼐用裆詈徒煌ㄟ\輸帶了巨大的挑戰(zhàn)［1］。

然而，受限于老撾當?shù)鼗A(chǔ)設(shè)施建設(shè)，缺乏足夠的水文設(shè)施支撐相關(guān)水文信息的采集，水文實測資料匱乏和水文資料序列較短成為影響鐵路排洪工程設(shè)計以及運營防排洪工作的重難點。因此，本文通過實測暴雨洪水數(shù)據(jù)，提出缺乏實測資料地區(qū)小流域暴雨洪水數(shù)值模擬系統(tǒng)建設(shè)和研究方法，構(gòu)建山區(qū)小流域三維高精度數(shù)字場景，研究不同降雨模式下植被截流、坡面飽和-非飽和入滲到坡面產(chǎn)匯流等多物理過程時空演化機制，構(gòu)建無資料山區(qū)小流域暴雨洪水數(shù)值預(yù)報物理模型。研究成果對缺實測水文資料地區(qū)的小流域暴雨洪水計算具有重要的參考價值。

1 缺實測水文資料地區(qū)暴雨洪水數(shù)值模擬數(shù)學(xué)模型

小流域暴雨洪水過程涉及降雨入滲、溝道匯流、洪水運動等多個物理階段，各階段內(nèi)在物理過程存在較大差異。因此，需針對暴雨洪水各物理階段分別建立數(shù)學(xué)模型以反映相應(yīng)的過程演化特征。本數(shù)值模型將降雨入滲、坡面匯流、洪水運動劃分為3 個階段，綜合考慮小流域地表土壤覆蓋層、小流域植被、小流域地形等因素對暴雨洪水形成過程的影響，建立小流域暴雨洪水過程的系統(tǒng)框架如圖1所示。

圖1 小流域降雨誘發(fā)暴雨洪水過程研究系統(tǒng)框架圖

小流域暴雨洪水成災(zāi)過程各物理階段的數(shù)學(xué)模型為：

（1）降雨入滲數(shù)學(xué)模型

小流域土壤蓄水是影響洪水水源條件的重要因素，降雨條件下的土壤入滲過程采用以含水率和基質(zhì)吸力耦合變量表達的Richard 方程描述［2］：

式中：t——時間；

θ——土壤含水率；

ψ——土壤內(nèi)部壓力水頭；

z——土壤厚度；

K——土壤滲透系數(shù)；

S——植被根系引起的水分吸收。

考慮到土壤在雨水下滲過程中，底部邊界對入滲過程的擾動，采用兩種不同的邊界條件（第一類Dirichlet 邊界和第二類Neumann 邊界）進行描述：

式中：θs——飽和土壤含水率；

θm——初始土壤含水率；

I——入滲強度。

（2）坡面匯流數(shù)學(xué)模型

小流域土壤經(jīng)歷入滲飽和后產(chǎn)生坡面積水，在地形條件的影響下產(chǎn)生匯流，進而為洪水形成提供水源條件。降雨條件下的坡面匯流過程采用Shallow-water方程描述［3］：

式中：h——坡面流水深；

x、y——坡面流運動方向；

R——降雨強度；

I——土壤飽和入滲強度；

g——重力加速度；

u、v——坡面流運動速度；

zb——地表高程；

Sfx、Sfy——坡面匯流所受到的摩擦阻力。

（3）洪水運動數(shù)學(xué)模型

坡面流匯入溝道的過程中，會帶動坡面松散物源啟動形成含沙洪水，這一過程會對洪水運動特征產(chǎn)生影響，洪水自身運動特征與坡面流存在不同。現(xiàn)場資料分析表明，洪水成災(zāi)過程中存在明顯的由侵蝕引起的體積放大效應(yīng)，同時也造成地形條件的改變。因此，在傳統(tǒng)Shallow-water 方程的基礎(chǔ)上考慮了物源啟動對洪水運動的影響，可描述為［4］：

式中：h——坡面流水深；

x、y——坡面流運動方向；

R——降雨強度；

I——土壤飽和入滲強度；

g——重力加速度；

u、v——坡面流運動速度；

zb——地表高程；

Sfx、Sfy——坡面匯流所受到的摩擦阻力；

c——洪水含沙量；

E——物源侵蝕速率；

D——物源沉積速率；

p——被侵蝕物源飽和度；

ρf——雨水密度；

ρs——被侵蝕物源密度。

（4）子物理模型耦合

小流域暴雨洪水成災(zāi)過程的各物理階段相互銜接，相互影響，存在顯著的階段互饋機制。因此，在建立的子物理模型基礎(chǔ)上，通過暴雨洪水形成過程中各階段之間的互饋原理，確定影響各階段銜接及過程演化的關(guān)鍵因子，對子物理模型開展耦合研究，如圖2所示。對于第一階段，選取降雨強度、土壤狀態(tài)、植被為關(guān)鍵因子，通過第一階段數(shù)學(xué)模型計算得出的土壤飽和度及坡面產(chǎn)流量等變量代入第二階段進行計算。選取地形、植被為第二階段的關(guān)鍵因子，通過第二階段數(shù)學(xué)模型計算得出的坡面匯流量等變量代入第三階段進行計算。選取地形、土壤狀態(tài)為關(guān)鍵因子，通過第三階段數(shù)學(xué)模型計算得出溝道洪水流量、洪水分布及含沙量等變量。

圖2 暴雨洪水形成過程各階段互饋機制及子物理模型耦合圖

（5）暴雨洪水數(shù)值模擬模型平臺構(gòu)建

針對小流域暴雨洪水過程中各階段物理模型方程的結(jié)構(gòu)特征，采用結(jié)合HLLC 黎曼間斷近似解的有限體積算法進行求解。基于MATLAB 語言開展代碼編寫，利用并行計算及網(wǎng)格重劃分技術(shù)提高計算效率，結(jié)合MATLAB 可視化功能，實現(xiàn)了計算數(shù)據(jù)的讀取與展示，以此完成暴雨洪水成災(zāi)全過程的計算模擬及定量化評價。小流域暴雨洪水過程流程如圖3所示。

圖3 小流域暴雨洪水災(zāi)害數(shù)值模擬流程圖

2 缺實測水文資料地區(qū)暴雨洪水監(jiān)測系統(tǒng)

2.1 監(jiān)測系統(tǒng)及監(jiān)測內(nèi)容

選取老撾鐵路2 個小流域作為示范區(qū)，對其降雨量及溝（河）道流量進行實時監(jiān)測。通過降雨-流量實時監(jiān)測，獲得代表性降雨產(chǎn)流成果，為數(shù)值模擬提供實測驗證資料。監(jiān)測儀器為雨量計和超聲波流量計。

2.2 監(jiān)測站點的布設(shè)原則

監(jiān)測站點選取和布設(shè)原則主要遵從以下3 點：①監(jiān)測點代表性；②監(jiān)測點信息傳輸可靠性；③監(jiān)測點施工及維護性。

為選擇具有代表性的監(jiān)測點，主要考慮以下因素：從監(jiān)測點所在的降雨分區(qū)、流域面積大小、鐵路跨河工程類型（橋梁、涵洞）。監(jiān)測設(shè)備平面布設(shè)如圖4所示。

圖4 監(jiān)測設(shè)備布置平面圖（m）

本研究從小流域站點的分布、鐵路跨河工程類型、下墊面植被類型，分別選取了2～6 m 香嫩涵洞、3×32 m 班那迷中橋兩種不同跨河類型的鐵路構(gòu)筑物，2 個站點流域特征如表1所示。

表1 磨萬鐵路典型小流域山洪監(jiān)測站點特征值表

雨量計應(yīng)選擇在地形開闊且無高喬木遮擋的位置，兼顧現(xiàn)場維護。因此將雨量計設(shè)在監(jiān)測流域靠近鐵路的開闊地帶。

流量計需對溝道內(nèi)水流的水位及流速進行實時監(jiān)測，雷達流量計應(yīng)位于河道水流正上方，確保下方河道無遮擋，無漂浮物阻礙雷達信號發(fā)射和傳輸。

3 小流域暴雨洪水數(shù)值模擬成果驗證

選擇實測典型降雨-產(chǎn)流實測資料進行數(shù)值模擬模型率定，以實際降雨資料模擬暴雨洪水形成過程。流域地形數(shù)據(jù)均由現(xiàn)場勘測DEM 獲得，假設(shè)降雨強度在流域范圍內(nèi)保持一致，流域降雨強度依據(jù)現(xiàn)場儀器監(jiān)測數(shù)據(jù)進行計算，最后將數(shù)值模擬成果與監(jiān)測數(shù)據(jù)進行對比，如圖5、圖6所示。

圖5 香嫩典型暴雨洪水數(shù)值模擬與實測數(shù)據(jù)對比圖

圖6 班那迷典型暴雨洪水數(shù)值模擬與實測數(shù)據(jù)對比圖

隨著降雨的持續(xù)，坡面開始產(chǎn)流。受地形影響，坡面流開始匯集，此時流量較小。值得注意的是，在流域內(nèi)坡面徑流形成過程中，由于地形影響，局部出現(xiàn)雨水積蓄現(xiàn)象。當蓄水高度超過洼地地形高度時，坡面流會繼續(xù)向下傳播。隨著坡面過流量的增加，坡面流逐漸趨于穩(wěn)定。計算監(jiān)測點最大過流量為28 m3／s，與實際測量值20 m3／s 較為接近，洪峰量線型與實測序列吻合。

以上2 個案例表明，所建立的耦合物理模型能夠較好地描述降雨誘發(fā)的坡面形成及傳播過程特征。

4 最大雨強條件下監(jiān)測流域暴雨洪水過程模擬

采用實測資料驗證的數(shù)值模擬模型，對磨萬鐵路2 個監(jiān)測流域開展最大雨強條件下暴雨誘發(fā)洪水形成過程的模擬分析，如表2所示。假設(shè)降雨強度在流域范圍內(nèi)保持一致，依據(jù)最大雨強進行計算，降雨時長到坡面產(chǎn)流達到穩(wěn)定狀態(tài)為止。

表2 研究流域5年至100年一遇雨強值表

由模型計算所得的各個流域溝道最大過流量與與二院法理論公式［5-6］的結(jié)果進行對比，如表3所示。兩者的結(jié)果較為接近，表明本文所使用的方法可較好應(yīng)用于暴雨誘發(fā)洪水形成過程描述及定量評價，可為降雨誘發(fā)流域坡面匯流及洪水危險性定量評價提供技術(shù)支撐［7-9］。

表3 計算最大流量與二院法公式計算結(jié)果對比表

5 結(jié)論

小流域暴雨洪水形成過程十分復(fù)雜，涉及多個演化過程，如入滲、坡面匯流、洪水溝道運動等。因此，精確模擬這一過程需要考慮小流域暴雨洪水形成過程中的關(guān)鍵環(huán)節(jié)。本文將小流域暴雨洪水形成過程劃分為入滲、匯流、運動3 個環(huán)節(jié)。針對每個環(huán)節(jié)開展物理建模，通過尋找各個環(huán)節(jié)之間的銜接要素，對不同環(huán)節(jié)對應(yīng)物理模型進行耦合得到了小流域暴雨洪水全過程物理模型，實現(xiàn)了小流域暴雨洪水形成過程的模擬計算。

利用所提出的小流域暴雨洪水耦合物理模型與數(shù)值計算方法，以中老鐵路磨萬段小流域工點為對象，以現(xiàn)場監(jiān)測數(shù)據(jù)為依照，開展了小流域暴雨洪水形成過程的模擬與對比。模擬結(jié)果表明，2 個小流域暴雨洪水的形成過程模擬符合實際情況，且模擬所得暴雨洪水流量過程線與實際測量數(shù)據(jù)較為吻合，反映了實際小流域暴雨洪水的過程演化特征，進一步驗證了小流域暴雨洪水耦合物理模型與數(shù)值計算方法實際應(yīng)用的可行性。

將所提出的小流域暴雨洪水模擬方法與理論計算法進行了對比驗證，以證明所提出小流域暴雨洪水模擬方法的可靠性。在相同邊界與初始條件（地形、最大雨強等）下，針對中老鐵路磨萬段小流域開展了暴雨洪水形成過程模擬。模擬結(jié)果表明，在不同流域條件下，通過本項目所提出小流域暴雨洪水模擬方法得到的洪峰流量值與二院法計算得到的洪峰流量值均十分接近，驗證了小流域暴雨洪水模擬方法的可靠性。