改進的暴雨洪水淹沒模擬算法及其在淮河流域的應用研究

盧燕宇 鄧汗青, 田紅, 何冬燕, 戴娟,

（1安徽省氣象科學研究所，合肥 230031；2 安徽省大氣科學與衛星遙感重點實驗室，合肥 230031；3 安徽省氣候中心，合肥 230031）

0 引言

由局地強降水造成的中小河流突發性洪水頻繁發生，易造成人員傷亡。由于大部分中小河流站網偏稀，缺少必要的應急監測手段，預報方案和能力有待健全與完善，加上中小河流源短流急，洪水具有暴雨強度大、歷時短、難預報、難預防的特點，針對暴雨誘發中小河流洪水和山洪地質災害的預報和防御成為目前防洪減災工作中突出的難點。充分利用暴雨后的短暫時間及時快速地做出洪水或洪澇預警決策，做好人員轉移，是防范洪澇災害保證人民生命安全的有效途徑。因此基于降水監測結果，開展暴雨洪水演進模擬可以提前識別可能遭受災害影響的承災體信息，有助于做出合理預判和決策，能夠更加有效地發揮預警實效。基于此，本研究將在現有暴雨洪水淹沒算法的基礎上進行改進，使之與現實更為貼近，將有助于發揮預警信息的最大效益，推動暴雨災害風險實時動態評估的業務化進程，促進氣象災害風險管理理念和方法的發展與創新。

IPCC《管理極端事件和災害風險，推進氣候變化適應》特別報告（SREX）中指出構成災害風險有3個要素，即致災危險性、暴露度和脆弱性，只有將三要素有機結合才能合理表征災害風險。就暴雨災害風險而言，致災危險性可由暴雨致災臨界閾值來反映，脆弱性體現承災體對暴雨災害強度的響應，暴露度則是連接危險性和脆弱性的關鍵環節，只有將空中降水轉化為地面徑流、洪水，才能判斷出可能受災害影響的承災體，因而需要借助暴雨洪水淹沒模擬來識別暴雨災害的風險暴露度[1,2]。暴雨災害的危害與洪水的淹沒范圍和水深直接相關，20世紀90年代以來，基于GIS空間分析技術，進行給定淹沒高程下的洪水淹沒分析已有不少研究[3-6]。但是目前大多研究仍然以靜態分析為主，尚不能快速提供不同的入流量、不同時相的淹沒狀況和水深分布。由于水動力模型能夠實現洪水演進的動態模擬，可以比較準確反映淹沒范圍、淹沒深度及其歷時特征，因而成為當前研究的一個熱點方向[7-11]。而水動力學模型研究的重要問題是如何在充分考慮洪水演進的物理機制的同時高效快速地實現洪水演進的動態模擬，以及如何將洪水模擬結果與社會經濟等數據相匹配結合，以評估災害影響。當前，將水動力學模型與GIS技術相結合[1,7]，為這一問題的解決提供了思路，但仍然處于探索階段，并且目前關于洪水淹沒的模擬多直接以河流水文特征量為起始值，忽略了與災害發生的根源——降水之間聯系的分析，在這些方面仍有待進一步的研究。

針對風險預警評估業務服務的迫切需求，以及現有暴雨洪水淹沒算法存在的問題，本研究將融合降水產流模型和二維水動力模型，改進淹沒模擬算法，使之反映前期降水和下墊面要素空間變化對暴雨洪水形成及其演進的影響；通過優化降水空間插值，開展滾動化疊加運算，改進算法的輸入和運行模式，使之體現降水時空變化對暴雨災害風險分布的影響。選擇淮河流域典型地區完成算法的檢驗和應用分析，從而為實時氣象防災減災提供技術支撐。

1 資料與方法

1.1 研究思路與方案

開展暴雨洪水淹沒模擬的關鍵是要抓住降水→徑流→淹沒三個環節，在這三個環節中需解決以下關鍵問題：1）降水精細化時空分布，2）降水向地表徑流的轉化，3）洪水的演進動態。在以往的暴雨洪水研究中這三方面存在一定的脫節，在開展洪水淹沒模擬時往往僅考慮地形因素對暴雨洪水演進的影響，對降水產流過程以及降水時空分布變化的影響考慮不足，因而影響了暴雨洪水淹沒模擬的仿真效果。因此本研究將著眼于現狀問題與迫切需求，通過產流系數將降雨產流模型（SCS）與二維水動力模型（Floodarea）相融合（圖1），使暴雨洪水淹沒算法能夠體現前期降水和下墊面特性的影響，提升仿真模擬效果；改進降水格網化方案和模型運算方式，實現滾動疊加運算，更好地體現降水時空分布變化對暴雨災害風險分布的影響；最后基于改進的模型算法，在淮河流域典型區內對模型算法和參數進行檢驗和優化，實現暴雨洪水淹沒的分布式動態模擬。

圖1 技術路線和算法流程Fig. 1 Technical route and algorithm flow

1.2 降水產流模型

構建降水量與地表徑流之間關系的關鍵問題是建立方便有效的降水產流模型。在對土壤入滲產流機理研究的基礎上，許多學者提出了不同的估算坡面產流過程模型，其中美國農業部根據美國氣候特征及多年水文徑流資料所研發的SCS-CN模型由于其結構簡單、所需參數較少、模擬結果準確度較高，而被廣泛應用于場次降雨地表產流及其過程的預測中。國內外一系列分布式—半分布式生態水文模型，如SWAT，EPIC，CREAMS等，均采用 SCS-CN模型預測徑流量。

SCS-CN模型包含一個水分平衡方程和兩個基本假設，表達如下：

水分平衡方程：

比例相等假設：

降雨初損量與潛在滯蓄量正比關系假設：

式中，P為降雨量（mm）；Q為徑流深（mm）；F為累計入滲量（mm）；S為流域最大蓄水能力（mm）；Ia為初損量（mm）；λ為初損率，無量綱。

聯立方程（1）～（3）可以得到計算徑流深的表達式：

S可由CN系數計算得到，二者有以下統計關系：

該模型的主要參數徑流曲線數CN是反映流域前期土壤濕潤程度（Antecedent moisture condition，AMC）、坡度、土壤類型和土地利用現狀等綜合特性的參數。SCS模型對CN的敏感性很高，有研究指出，CN取值±10%的變化可導致計算徑流量-45%～55%的變化?？梢?，CN值的確定對降雨徑流量的準確預測非常重要。除CN系數外，初損率λ以往通常取標準值0.2，但近年來研究發現初損率隨不同區域取值變異性也較大，其值變化對模擬效果有較為明顯的影響。

1.3 二維水動力模型

洪澇災害的危害與洪水的淹沒范圍和水深直接相關，因此，確定災害范圍和程度可通過模擬洪水演進及其水文特征來實現。

本研究引入了德國Geomer公司研制的基于GIS的水動力模型——FloodArea模型。其原理是充分利用GIS柵格數據在水文-水動力學建模上的優勢，實現GIS與水文-水動力學模型的數據融合。模型以柵格為基本單元，淹沒模擬基于二維非恒定流水動力學模型，用Manning-Strickler公式計算每個柵格與周圍柵格之間的水量交換。用Manning-Stricker公式計算每個柵格單元與周圍8個單元之間的洪水流量，FloodArea模型計算洪水匯流示意圖見圖2，相鄰單元的水流寬度被認為是相等的，位于對角線的單元，以不同的長度算法來計算；圖中R為相鄰單元的柵格距離，R為對角線單元的柵格距離，陰影部分指柵格面積，箭頭指水流方向。水流的淹沒深度為淹沒水位高程和地面高程之間差值，由下式表示：

淹沒過程中的水流方向由地形坡向所決定，地形坡向反映了斜坡所面對的方向，坡向指地表面上一點的切平面的法線矢量在水平面的投影與過該點的正北方向的夾角，表征該點高程值改變量的最大變化方向，計算公式如下：

其中，α為地形坡度。FloodArea模型有3種基本的淹沒情景：漫頂、潰口以及暴雨，能夠以給定水位，給定流量或給定面雨量3種方式進入模型，并可根據水文過程線進行實時調整，可視化表達流向、流速和淹沒水深等水文要素的時空物理場，為洪水淹沒風險動態模擬提供了有效工具。

圖2 FloodArea匯流計算原理示意[12]Fig. 2 Schematic illustration of the confluence calculation in FloodArea model

1.4 研究區概況及資料

以淮河流域中的史河上游為例開展算法改進實驗和應用研究。研究樣區如圖3所示，位于大別山腹地金寨縣境內，以黃泥莊水文站作為研究區流域的出水口斷面，獲取了2008—2010年多場洪水過程的水文數據，用于計算實際徑流深。研究河段屬山區型河道，水淺流急，支流較多；流域面積約為800 km2，海拔高度懸殊大，研究區內高差超過1200 m，土壤類型以黃棕壤為主，植被覆蓋類型主要為林地。研究區及周邊布有多個氣象觀測站，可獲取逐時雨量觀測值用來計算流域面雨量。

圖3 研究樣區概況、地形、水文站和氣象站分布Fig. 3 Overview of the study area, topography and distribution of hydrological and meteorological stations

下墊面信息包括數字地面高程（DEM）、土地利用類型等，DEM數據主要來自于1:5萬地理信息數據，根據研究區范圍進行裁剪得到，數據分辨率為25 m，土地利用類型來自于GLC30土地類型數據庫，分辨率為30 m。氣象數據來自于全省自動氣象站的逐小時降水數據集。水文數據主要來自于安徽省水文遙測信息網。

2 結果與討論

2.1 降水空間插值方案的優選

為使氣象數據與模型輸入相匹配，獲取精細化的降水時空分布特征，采用空間插值的方法將離散的站點數據插值成細網格化的面狀數據[13]。插值方法將考慮反距離權重（IDW）、自然鄰域、克里金（Kriging）和薄盤樣條 等方法，同時根據降水空間分布特點，采用海拔高度、地形起伏度、坡度、坡向等地形因子作為協同變量，構建空間插值模型。采用交叉驗證（Cross-validation）的方法來比對不同插值方法和模型的實際效果，從中選擇最優化處理方案，作為臺站降水數據的空間拓展和網格化推算的方法。

選擇近年來多次降水過程進行插值方法對比分析，其中包括了熱帶氣旋、梅雨鋒等不同天氣系統所導致的強降水過程，具有一定的代表性。采用交叉驗證的平均相對誤差、均方根誤差來衡量模擬效果，由于研究重點關注暴雨，因而進一步引入納什效率系數來考察對于極值的擬合效果。從圖4可以看出就平均相對誤差來看，IDW和RBF 法擬合值總體比觀測值偏小，其他3種方法則略偏高，其中以Co-kriging方法誤差最?。粚τ诰礁`差來看，RBF法的誤差最大，而Kriging法的估計誤差最??；納什效率系數越接近于1說明擬合結果與觀測值越接近，可以看出同樣RBF法的估計偏差最大，而Kriging法則與觀測結果最為接近。綜合3種指標可以看出以Kriging法的擬合效果最佳。

圖4 不同降水插值方法交叉驗證的平均相對誤差（a）、均方根誤差（b）和納什效率系數（c）Fig. 4 Mean relative error (a), root mean square error(b) and Nash efficiency factor (c) for cross-validation of different precipitation interpolation methods

2.2 降水產流模型的參數分析與應用

2.2.1CN系數計算

SCS模型中CN系數反映了流域下墊面的產流能力，對后續降水產流及淹沒模擬計算具有重要意義，CN值與土壤類型、土地利用類型以及前期土壤濕潤程度密切相關。根據研究區土壤類型和土地利用類型數據，通過查表可以得到中等濕潤程度下CN值分布（圖5），再通過計算研究區前5 d累計面雨量即可獲取實時的CN系數。

由于SCS模型構建和參數化過程中主要依托于美國大量緩坡（坡度5%左右）降雨產流資料發展起來，而研究區地形復雜，海拔起伏大，在實際應用時還應進一步考慮坡度對CN系數的影響。根據已有研究結果，采用Huang等[14]的坡度修正公式（見式（8））對CN系數進行修正，使之適用于起伏山地的降雨產流模擬。

圖5 研究區土壤類型（a）、土地利用類型（b）和AMC II 下的CN值（c）分布Fig. 5 Distribution of soil type (a), land use type (b) and CN values (c) under AMCII in the study area

式中：CNII為中等濕潤程度下的CN系數，slp為平均坡度（%），CNIIs為坡度修正后的CN系數。

基于DEM計算了研究區的坡度分布，并采用公式（8）對CN系數進行了修正，確定了AMCII條件下的CN值分布結果（圖6），作為SCS模型應用的輸入參數。

圖6 研究區坡度（a）及其修正后的CN系數（b）分布Fig. 6 Distribution of slope (a) and its corrected CN coefficient (b) in the study area

2.2.2 初損率λ優化

當確定了AMC II條件下的CN系數后，只需要計算API判定實時AMC條件，根據SCS模型公式可知僅有一個未知參數即初損率λ，在以往的SCS模型應用中，一般取λ=0.2，但該值是根據美國土壤水文觀測資料推算而來，是否適用于本地區還需進一步研究，這里選擇2008年的多次典型降水過程分析λ的取值對模擬結果的影響，并確定優化后的參數值。

參數優化的依據為確定系數（R2）、相對誤差（PBIAS）和效率系數（NSE）。使用確定系數R2評價模型追蹤觀測值變化的準確程度，越接近于1表明模型越能準確反映變化趨勢；利用相對誤差PBIAS評價模擬結果高于或低于觀測值的平均趨勢，該值大于0表示模型低估產流量，小于0表示模型高估產流量，越接近于0表明偏差越??；采用NSE表示模擬結果與觀測值的契合程度，越接近于1表明模型越準確。

從表1可以看出模擬徑流深基本上隨著λ的增加而減少，這也可以用λ的物理意義來解釋，即λ增大，降水損失量增加，轉化為徑流的比例相應減少。對比歷次過程，模擬徑流深基本上是當初損率λ為0.05～0.15范圍內與實測結果較為接近。進一步通過R2、PBIAS和NSE來看，當λ=0.10時SCS模型模擬效果表現最佳（圖7），因此選擇0.10作為研究區的λ的優化取值。

表1 典型洪水過程實測與模擬徑流深Table 1 Observed and simulated runoff depths for typical flood processes

2.2.3 SCS 模型應用及檢驗

在確定了SCS模型的CN系數和λ后，進一步建立不同條件下CN系數的概化轉換公式（圖9 ），以中等濕潤條件B類土壤的CN系數作為模型輸入參數，然后結合API的分布結果進行CN系數的實時計算，從而實現降水產流量的實時分布式模擬（圖8）。

基于模型應用流程，采用率定好的模型參數集，以降雨量作為輸入實時模擬了2010年多次降水過程的產流量（表2），模擬與實測的相關系數在0.9以上，相對誤差約為-15%，模擬值略小于實測?？偟膩砜?，率定后的SCS模型能夠較好地適應研究區，可以實時反映降水產流關系，從而為后續暴雨洪水淹沒模擬提供可靠的輸入前提。

2.3 二維水動力模型的參數分析與應用

在二維水動力模型FloodArea中影響模型效果的主要控制參數為Manning-Sticker公式中的地面糙率，根據已有研究表明，糙率系數與土地利用類型有密切關系，參考相關研究及標準規范[15]，按不同土地利用類型設定了糙率系數（表3）。該套系數已在大別山區和淮河王家壩等地區經過了驗證[11,16]，能夠較好地反映實測淹沒水深的時空動態分布。

在以往暴雨洪水淹沒模擬中，多以過程模擬和事后評估為主，降水輸入模型時僅體現降水隨時間變化，不考慮降水空間分布的動態變化，即在一次過程模擬中，暴雨中心位置常常假定不變。本研究選擇史河上游2010年7月11日暴雨過程為例研究FloodArea模型在強降水致洪淹沒分析中的應用，基于實況觀測資料采用Kriging方法計算了該日逐小時面雨量，經統計24 h累計面雨量超過80 mm，流域出口斷面黃泥莊水文站水位上漲接近1 m?；?:5萬地理信息數據和糙率信息，將降水空間分布和過程雨量帶入到FloodArea模型進行洪水淹沒動態分析，同時為了考察降水輸入方式對模擬結果的影響，我們進一步改變了模型輸入方式，采用滾動迭代輸入將逐小時的降水空間分布代入模型，進行疊加運算。

圖7 不同初損率模擬結果的確定系數（a）、相對誤差（b）和效率系數（c）Fig. 7 Determination factor (a), relative error (b) and efficiency factor (c) for different initial loss rate simulations

從圖9可以看出，在整個過程中降水中心的改變十分明顯，即使是比較接近的時刻降水的空間分布仍有較大變化。而如果在模型輸入時采用固定的降水空間分布將不可避免地帶來較大的偏差，例如在本例中過程累計雨量主要集中在流域的南部，而實際上在過程中降水中心時常在南北擺動。

對比兩種輸入方式的模擬結果（圖略），可以看出淹沒水深的空間分布具有一定差異，主要表現為相比于默認輸入方式，滾動模擬結果流域北部淹沒范圍更大，而南部相對淹沒范圍較小，水深較淺。進一步采用黃泥莊水文站的觀測數據進行分析（圖10），可以看出采用默認的固定降水輸入模擬的水深變化略滯后于實測，而滾動模擬結果則與實測較為同步，這可能是由于過程累計降水中心距流域出口斷面較遠，采用固定的降水空間分布使得匯水時間增長，從而流域出口的水位變化與實際相比存在滯后。

總的來看，FloodArea模型對水深的動態變化模擬效果較好，能夠較為準確地反映水位的漲落情況，但兩種輸入方式下水深均顯著高于實際，這主要是由于FloodArea模型在模擬過程中假定降水全部轉換為徑流參與匯流和淹沒分析，而實際情況則是在降水產流過程中，有部分降水以下滲等方式損耗，產流量要低于降水量，因此在模擬中還需要進一步考慮產流系數，將降水轉化為產流降水才能與實際情況更為接近。

圖8 SCS模型實時應用流程示意Fig. 8 Schematic representation of the SCS model real-time application flow

2.4 模型耦合及算法驗證

2.4.1 耦合思路與算法改進

基于精細的格網化降水量，通過引入SCS降水產流模型，綜合前期降水量、土壤和土地利用類型、地形因子等參數來計算歷次降水的產流系數，利用產流系數來對實時降水進行修正，并計算每個格點上實際產生的地表徑流量，之后將實際產流量代入到FloodArea中運算模擬暴雨洪水的動態演進情況。即以“產流系數”這個關鍵環節來將SCS降水產流模型與FloodArea洪水淹沒模型相融合，形成分布式的暴雨洪水淹沒算法，能夠充分反映流域內前期降水和下墊面要素空間變化對暴雨洪水形成及其演進的影響。

此外，在以往暴雨洪水淹沒模擬中，多以均勻雨量輸入模型，不考慮降水時空變化的影響，本文進一步采用滾動迭代的方式進行逐小時模擬以體現降水空間分布隨時間的變化，以實現更好的仿真效果。

圖9 累計降水分布（a）和降水空間格局動態變化（b）Fig. 9 Cumulative precipitation distribution (a) and spatial pattern dynamics of precipitation (b)

圖10 不同降水輸入方式下的淹沒模擬結果與觀測值對比Fig. 10 Comparison of inundation simulation results with observations for different precipitation inputs

表2 典型暴雨過程的觀測與模擬徑流深Table 2 Observed and simulated runoff depths for typical storm events

表3 不同土地利用類型對應的糙率系數Table 3 Roughness coefficients corresponding to different land use types

2.4.2 改進后算法模擬效果分析

為進一步驗證改進算法的模擬效果，仍然采用2010年7月11日史河上游的模擬案例，從圖11可以看出，在考慮了降水產流和損耗后，模擬的黃泥莊水深變化與實際更為接近，動態趨勢吻合度更好，并且系統偏差也明顯縮小，改進后的算法對暴雨洪水淹沒模擬誤差比之前縮小約50%。

3 結論

針對現有暴雨洪水淹沒算法存在的問題，本文在降水空間插值、降水產流模型參數優化、二維水動力模型本地化、模型耦合及應用等方面開展了系統分析，融合SCS模型和FloodArea模型改進淹沒模擬算法，使之實現更好的仿真效果。主要結論包括：

1）交叉驗證結果表明以Kriging方法插值效果最佳，因此選用Kriging法作為降水空間插值和面雨量計算的優選方法。

2）基于下墊面信息實現了SCS模型的參數化，采用地形坡度對CN系數進行了修正；并通過典型個例模擬及率定，對SCS模型中初損率值進行了優化，率定后的SCS模型能夠較好地反映降水產流關系，從而為后續暴雨洪水淹沒模擬提供可靠的輸入前提。

3）根據土地利用類型確定了FloodArea模型的地表糙率系數，降水輸入方式對洪水動態模擬結果有一定影響，改進后的滾動模擬方式可以實現更好的仿真效果。

4）以降水產流過程為紐帶實現了SCS模型與FloodArea模型的耦合，形成分布式的暴雨洪水淹沒算法，來反映流域內前期降水和下墊面要素空間變化對暴雨洪水形成及其演進的影響，改進后算法對模擬效果有顯著提升作用，模擬誤差縮小約50%。

總的來說，通過參數優化、模型耦合等方法對現有暴雨洪水淹沒算法進行了改進，使之與現實更為貼近。由于本研究重點關注降水落到地面后的產匯流過程及其可能帶來的洪澇災害，因此主要以降水實況監測作為輸入來研究算法模擬效果。未來為進一步延長預見期，還有必要將預報降水與模擬算法相結合，探討誤差來源，優化匹配策略，以更好地與風險預警服務業務相銜接，在發布災害性天氣預報后，如果達到或超過臨界條件，進一步開展暴雨洪水淹沒模擬分析，指出災害高風險區及可能災損，為政府和群眾提供直觀的暴雨災害風險信息，有助于發揮預警信息的最大效益。

圖11 改進后算法的模擬效果對比Fig. 11 Comparison of the simulation effects of the improved algorithm