基于MIKE耦合模型的西苕溪洪水演進及風險分析

李小天，戚 藍，胡琳琳，滕 暉，施 征，甄億位

(1.天津大學水利工程仿真與安全國家重點實驗室，天津300072；2.浙江省水利河口研究院，浙江杭州310020)

0 引 言

洪水災害目前仍是世界上最為頻發的自然災害之一，每年往往會造成大量的人員傷亡和巨大的財產損失。隨著城市化進展的不斷加快，下墊面原有的組成受到破壞；從而導致防洪形式發生新的變化。人們逐漸意識到單純依靠工程措施已經不能滿足防洪要求；同時，超標準洪水風險已經對城市發展構成極大的威脅[1]。正是在這種背景下，防洪戰略從“控制洪水”轉為“洪水管理”顯得尤為必要。洪水管理是人類按照可持續發展的原則，以協調人與洪水關系為目的，理性規范洪水調控行為與增強適應能力等一系列活動的總稱[2]。

洪水風險管理包括洪災風險分析、洪災風險評價和洪災風險管理與決策。洪災風險分析是洪災風險評價的前提，而洪災風險評價又是洪災風險管理和決策的依據。城市洪災風險分析是指借助洪水模擬技術對城市不同強度暴雨洪水可能形成的淹沒特征進行研究。目前，采用的模擬方法包括水文學和水力學。隨著計算機及數值計算技術的提高，城市洪澇模擬的水力學方法得以迅速發展并在城市雨洪模擬中取得了許多成功。Mark等[3]采用隱式差分求解圣維南方程組，建立了城市街道與地下管網聯動的一維水動力學模型，成功模擬了地表一維非恒定流過程。Schmitt等[4]采用有限體積法求解淺水波方程組，用于模擬城市暴雨地表二維非恒定流過程。張新華等[5]建立的任意多邊形網格2D FVM模型在時空上均具有二階精度，應用于城市的二維非恒定流模擬取得了較好的計算結果。

圖1 計算范圍及阻水建筑物概化示意

本文以西苕溪干流長潭以下至長興界之間的區域為例，構建了一、二維水動力耦合模型，模擬了多種計算方案下洪水淹沒演進過程并進行災害風險分析，為研究區防洪減災和洪水風險管理等提供重要參考。

1 研究區域概化及模型構建

1.1 區域概化

西苕溪干流長潭以下至長興界河道蜿蜒扭轉，水流不暢，為洪水多發區；所以，將此區域作為此次洪水風險分析的研究區域。該區對應河長約48.0 km，區域面積約180 km2，此處河段彎曲支流上溯到大溪蔣塘里、滸溪白水灣及遞鋪溪鳳凰水庫壩下。根據洪水來源、河道走向、鄉鎮道路等阻水建筑物的分布情況，將計算邊界概化為：將西溪賦石水庫、南溪老石坎水庫、大溪蔣塘里、滸溪白水灣、遞鋪港鳳凰水庫、渾泥港、昆銅港等7處作為上邊界，采用洪水流量過程線；將港口水文站作為下邊界，采用水位過程線；研究區涉及支流概化為9條集中入流和12段沿程入流，概化水工建筑物84個，其中堰壩53座、水閘4座、阻水橋梁27座；最終得到的洪水風險分析計算邊界及主要阻水建筑物(見圖1)。

1.2 模型方案構建

考慮安吉縣西苕溪流域洪水特性、基礎資料等因素，本次采用河道一維和區域二維水動力耦合模型模擬洪水要素。基于MIKE FLOOD對MIKE 11和MIKE 21進行耦合，分析河道洪水演進及區域洪水淹沒過程。

(1)一維水動力模型(MIKE11)。采用MIKE11河流模型建立適合西苕溪流域的一維河道非恒定流數學模型。該數學模型對分析區域的河道及現有工程做了合理概化處理，以便較好地反映流域的水流特性。連續性方程

(1)

運動方程

(2)

式中，Q為斷面洪水流量，m3/s；A為過水斷面面積，m2；q為均勻旁側入流的單寬流量，m2/s；x為沿程距離，m；t為時間，s；R為水力半徑，m；C為謝才系數，s/m1/3；Z為斷面水位，m；α為動量修正系數。

(2)二維水動力模型(MIKE 21)。堤防外平面水流模擬，采用MIKE 21中的FM模塊。該模塊的計算采用非結構化網格劃分系統，能有效控制劃分數量和單元數量。模型依據的是描述平面水流運動的二維非恒定流方程組，共包括三個方程：水流連續性方程、水流沿x方向的動量方程及水流沿y方向的動量方程。即

(3)

(4)

(5)

式中，Z為水深，m；H為水位，m；q為連續方程中的源匯項，m/s；M為垂直于x方向的平均單寬流量，m2/s；N為垂直于y方向的平均單寬流量，m2/s；u為x方向的平均流速，m/s ；v為y方向的平均流速，m/s ；n為糙率；g為重力加速度，m/s2。

(3)一、二維模塊耦合(MIKE FLOOD)。MIKE FLOOD模型耦合連接一、二維模塊計算，采用側向連接方式耦合兩個模塊。側向連接允許MIKE 21的網絡從側面連接到MIKE 11的部分河道或整個河道，利用結構物能量公式來計算通過側向連接的水流。潰堤部分計算采用能量方程，一維模塊提供潰口位置河道側水位，二維模塊提供潰口位置堤后水位，兩側水位按照模型時間步長進行實時更新，互相影響，并維持模型的總水量守恒。能量方程

(6)

式中,Z1、Z2為上、下斷面的水位，m；hf、hj為上下斷面間的沿程、局部水頭損失，m；u1、u2為上、下斷面的流速，m/s；α1、α2為上、下斷面的動能修正系數；g為重力加速度，m/s2。

本次采用側向連接的方式耦合河道一維模型與二維保護區模型，以適用性較好的堰流公式計算通過側向連接的水流，模擬一維河道洪水以漫堤的方式向二維保護區演進過程。模型共設置側向連接45處。一、二維耦合模型如圖2所示。

圖2 一、二維耦合模型示意(單位:m)

2 模型參數率定與驗證

2.1 參數率定

結合研究區域的歷史洪水資料，選擇2009年作為參數率定年份，2013年作為驗證年份。在水動力模型中，需要率定的主要參數是糙率系數，糙率是反映水流流動時受阻程度的一個綜合系數。本研究對不同的河段分別賦予不同的糙率。參數率定時將已有的參數成果作為糙率初值，在此基礎上不斷進行調試。

表1 計算區域糙率取值分布

2.2 率定結果

將梅溪水文站與橫塘水文站實測資料與模型計算值進行比較，通過不斷調整參數，使模擬結果達到最優。2009年的率定誤差結果見表2。

表2 洪水率定計算成果

由表2可知，率定誤差滿足《洪水風險圖編制技術細則(試行)》規定的最大水位誤差絕對值小于等于0.2 m，計算值與實測值吻合較好。

2.3 參數驗證

為了檢驗西苕溪流域洪水演進計算模型的準確性及合理性，根據率定的模型參數，選取2013年“菲特”典型洪水進行參數驗證，驗證誤差結果見表3。

表3 洪水驗證計算成果

由表3可知，梅溪站與橫塘站最大水位誤差絕對值分別為0.05 m和0.08 m，遠小于《洪水風險圖編制技術細則(試行)》規定的最大水位誤差絕對值小于等于0.2 m，說明計算采用的網格概化及參數是合理有效的。

3 洪水演進計算及風險分析

3.1 計算方案設定

根據分析安吉西苕溪流域防洪標準和邊界洪水來源分析，本文擬定計算方案數量11個，包括5個漫溢方案和6個潰堤方案，如表4所示。

表4 計算方案

3.2 計算成果

根據計算結果可知，方案1～5的洪水淹沒面積分別為2.81、8.77、21.84、37.07、45.57 km2。5年一遇仍有少量的淹沒面積，表明西苕溪流域部分河段防洪能力較低，洪水淹沒概率較高。橫塘和梅溪站個重現期下設計水位見表5，主要漫堤點見表6。

以方案4為例，漫堤發生在泥河沙河匯合口～草蕩水庫左、右岸和西苕溪左岸安城。不同淹沒水深對應的淹沒面積、不同淹沒歷時對應的面積見表7和表8，淹沒水深、淹沒歷時見圖3和圖4。

表5 水位成果 m

表6 主要漫堤點位統計 m

表7 淹沒水深—淹沒面積關系

表8 淹沒歷時—淹沒面積關系

由表7和表8可知，西苕溪流域遭遇50年一遇洪水時最大淹沒面積為37.07 km2，以大于3 m水深為主，占總淹沒面積的53.96%。淹沒面積隨著淹沒歷時先增加后減少，以24～72 h為主，占總淹沒面積的80.32%。

圖3 西苕溪流域50年一遇漫堤淹沒水深

3.3 合理性分析

(1)水量平衡。根據模型運算結果統計出計算區域上游來水量、下游出水量、河道槽蓄量與淹沒區域積水量，各變量滿足洪水淹沒區域積水量=上游來水量-下游出水量-河道槽蓄量，說明構建模型洪水計算滿足水量平衡[6]。

(2)淹沒面積與重現期相關性[7]。統計出不同重現期洪水淹沒面積，并繪制其與重現期的相關圖，如圖5所示。由圖5可知，隨著設計洪水重現期的增加，淹沒面積隨之增大，但增速隨著重現期的增大而減小，這是由于設計洪量增速隨重現期的增大而減小，與實際相符，說明計算結果合理。

圖5 淹沒面積—重現期關系

(3)高防洪能力區域分析。針對西苕溪流域100年一遇淹沒水深分布圖中未受淹的區域，將本次淹沒圖層疊加到2 m×2 m精度DEM上，分析未受淹位置與周圍受淹區域之間的地形關系。分析可知，未受淹位置高程明顯高于周圍受淹區域高程，說明淹沒范圍分布合理。

(4)束窄斷面流場分析。以涉水建筑物橋梁為例，分析涉橋處模擬流場分布。由模擬結果可知，洪水進入束窄斷面后流場收縮，流速增大，水位抬高，符合實際情況。

4 結 論

(1)本文選取西苕溪干流長潭以下至長興界之間河段作為研究區域，基于實測河道斷面資料和3 m×3 m高精度DEM構建了MIKE一二維耦合模型，選用歷史洪水數據進行模型參數率定驗證。利用水文學方法給定各個計算方案對應的模型邊界條件，進行洪水演進及風險分析。

(2)本文從水量平衡、淹沒面積與重現期相關性、淹沒范圍與地形相關性以及流場分布對計算結果進行了合理性分析，研究表明，建模精度滿足精度要求、結果較為合理可靠。可為當地防洪減災制定策略提供參考。

(3)針對發生五年一遇洪水就淹沒的防洪能力較低的地區，應當對河道進行疏浚清淤，適度拓寬河道，加高加固堤防，提高河道過流能力，從而提高河道整體防洪能力。