基于MIKE耦合模型的平原區(qū)流域洪澇過(guò)程模擬

趙華青，周 璐，趙然杭，李占華，齊 真

（1.山東大學(xué)土建與水利學(xué)院，濟(jì)南 250061；2.北京水務(wù)咨詢(xún)有限公司，北京 100038；3.山東省海河淮河小清河流域水利管理服務(wù)中心，濟(jì)南 250000）

0 引 言

我國(guó)平原流域多屬于溫帶季風(fēng)區(qū)，汛期降雨強(qiáng)度較大，加之河流比降較低，過(guò)流能力有限，常因行洪不暢，引發(fā)洪澇災(zāi)害，容易造成嚴(yán)重的財(cái)產(chǎn)損失甚至人員傷亡，對(duì)一、二產(chǎn)業(yè)及生態(tài)環(huán)境造成了難以計(jì)量的損失［1］。洪澇災(zāi)害防治是目前我國(guó)水利建設(shè)的薄弱環(huán)節(jié)之一，需同時(shí)采用洪澇災(zāi)害防治工程措施與非工程措施，減少災(zāi)害損失［2］。其中，基于水文模型的洪澇過(guò)程模擬的預(yù)報(bào)預(yù)警，是洪澇災(zāi)害防治重要非工程措施之一。因此，實(shí)現(xiàn)平原區(qū)流域洪澇過(guò)程的精準(zhǔn)模擬，對(duì)于防洪減災(zāi)和洪水管理十分必要。

水文模型是以水文學(xué)、流體力學(xué)及水動(dòng)力學(xué)等理論為基礎(chǔ)，依托數(shù)學(xué)及物理模型，模擬洪水過(guò)程。國(guó)內(nèi)外水文模擬大部分使用概念模型與集總模型［3］，分布式水文模型因考慮了氣象及下墊面的空間分布，相較于概念模型及集總式模型，能更精確地模擬流域水文變化及分布過(guò)程［4］，從而能更好地分析流域洪澇災(zāi)害分布情況。在我國(guó)應(yīng)用較為廣泛的具有代表性的分布式水文模型主要有SWAT 模型、VIC 模型、SHE 模型等，其中丹麥水利科學(xué)研究所研制的MIKE SHE 模型以標(biāo)準(zhǔn)化網(wǎng)格為輸入單元，多模塊耦合的方式來(lái)表達(dá)水文循環(huán)過(guò)程［5］，在模擬坡面洪水流速、水深等方面具有顯著優(yōu)勢(shì)。

暴雨導(dǎo)致的超出流域承載力極限的洪水量是導(dǎo)致平原區(qū)洪澇災(zāi)害的主要因素，洪澇過(guò)程模擬是進(jìn)行洪澇災(zāi)害評(píng)價(jià)工作的基礎(chǔ)［6］。相較于集總式水文模型，分布式水文模型充分考慮了流域內(nèi)降水與下墊面要素的空間變化，對(duì)下墊面復(fù)雜的平原流域具有很好的適用性。因此，建立適合平原區(qū)特點(diǎn)的洪澇模型，對(duì)平原區(qū)流域洪澇過(guò)程進(jìn)行分布式模擬，對(duì)于流域防洪減災(zāi)和洪水管理十分必要。

本文以海河流域的徒駭河流域?qū)m家攔河閘上游為研究區(qū)域，利用水文站監(jiān)測(cè)數(shù)據(jù)，以GIS技術(shù)為支撐，建立分布式水文-水動(dòng)力耦合模型，模擬洪澇過(guò)程，包括河流外流域水深、流速分布及河流徑流過(guò)程，探討耦合模型在平原區(qū)洪澇模擬中的適用性，試圖為研究平原流域洪水模擬、洪澇災(zāi)害評(píng)估等提供科學(xué)的方法。

1 研究區(qū)概況和數(shù)據(jù)來(lái)源

1.1 研究區(qū)概況

山東省海河流域位于華北平原，山東省西北方，流域內(nèi)有徒駭河、馬頰河、德惠新河三條主干河流，流域面積約為2.97 萬(wàn)km2，流域地勢(shì)低平，洪澇災(zāi)害頻繁，防洪非工程措施不完善，洪澇災(zāi)害形勢(shì)較為嚴(yán)峻。研究區(qū)域位于山東省海河流域西南部，面積約為6 700 km2，主要河流為徒駭河上游段，河流方向?yàn)槲髂?東北［7］，有新金線(xiàn)河、趙牛新河等八條主要支流（見(jiàn)圖1）。流域形狀狹長(zhǎng)，坡度較小。流域降水具有顯著的季節(jié)性，豐水年與枯水年周期性交替，旱澇急轉(zhuǎn)與連旱連澇現(xiàn)象發(fā)生頻繁。研究區(qū)域內(nèi)河流均屬季節(jié)性河流，具有坡度小，匯流時(shí)間長(zhǎng)，洪水持續(xù)時(shí)間長(zhǎng)，水量年內(nèi)變化大等特點(diǎn)。

圖1 徒駭河宮家攔河閘上游流域圖Fig.1 Upper reaches of Gongjia Sluice of Tuhai River

1.2 數(shù)據(jù)來(lái)源

洪澇模擬需要的基礎(chǔ)數(shù)據(jù)包括DEM、歷史雨水情數(shù)據(jù)、土地利用類(lèi)型數(shù)據(jù)、土壤類(lèi)型數(shù)據(jù)、河流斷面數(shù)據(jù)等多種類(lèi)型數(shù)據(jù)。其描述和來(lái)源如表1所示。

表1 徒駭河宮家攔河閘上游流域MIKE耦合模型數(shù)據(jù)Tab.1 MIKE-A-R model data of Gongjia Sluice of Tuhai River

2 MIKE-A-R耦合模型構(gòu)建

2.1 MIKE耦合模型

（1）MIKE SHE 模型。MIKE SHE 是一種具有確定性以及基于物理機(jī)制的流域分布式水文模型。通過(guò)劃分流域網(wǎng)格并建立相鄰網(wǎng)格間的拓?fù)潢P(guān)系，表達(dá)參數(shù)的空間差異性：平面網(wǎng)格劃分形式為矩形網(wǎng)格，每個(gè)網(wǎng)格獨(dú)立作為一個(gè)輸入單元，表達(dá)流域下墊面的空間差異性；在豎直方向劃分水平層，模擬不同土層的洪水運(yùn)移過(guò)程［8］。MIKE SHE 模型通過(guò)蒸散發(fā)、坡面流、河流與湖泊、非飽和帶以及飽和帶6個(gè)模塊相互耦合，模擬完整的流域產(chǎn)匯流過(guò)程。

（2）MIKE 11 模型。MIKE 11 采用一維圣維南方程組，模擬河流洪水演進(jìn)過(guò)程：

式中：A為過(guò)水?dāng)嗝婷娣e；Q為過(guò)水?dāng)嗝媪髁?；t為時(shí)刻坐標(biāo)；x為沿河流的距離；q為河流旁側(cè)入流流量；h為水位；C為謝才系數(shù)；R為過(guò)水?dāng)嗝媪髁颗c面積的比值；α為動(dòng)量修正系數(shù)。

2.2 研究區(qū)域柵格與邊界提取

根據(jù)目前DEM 分辨率，基于GIS 技術(shù)提取河網(wǎng)的方法僅適用于地表比降較大，河網(wǎng)層級(jí)清晰的流域，針對(duì)地勢(shì)平緩、地貌起伏不顯著的研究區(qū)域，流域響應(yīng)單元?jiǎng)澐州^為困難［9］。為與流域?qū)嶋H相一致，需在研究區(qū)官方水系圖的基礎(chǔ)上進(jìn)行地理比對(duì)，進(jìn)而確定流域范圍。

在確定流域范圍過(guò)程中，為保證流域所有區(qū)域洪水均能匯入河流，對(duì)DEM采取以下處理：

降低河流對(duì)應(yīng)的DEM 柵格高程值，使其顯著低于兩岸高程，同時(shí)河流外柵格高程保持不變，然后通過(guò)ArcGIS 水文分析工具確定流域邊界，最后通過(guò)對(duì)比實(shí)際地形地貌及地圖影像做人工修正。

2.3 MIKE-A-R耦合模型

（1）河流外洪澇模型（MIKE-A）的構(gòu)建。根據(jù)流域的特點(diǎn)與資料情況，選擇MIKE SHE 模型各子模塊的計(jì)算方法（詳見(jiàn)表2）。結(jié)合DEM、流域邊界、土地利用類(lèi)型、土壤類(lèi)型及雨水情等數(shù)據(jù)，構(gòu)建河流外洪澇模型（MIKE-A）。

表2 子模塊計(jì)算方法Tab.2 Submodule calculation method

（2）河流水文模型（MIKE-R）構(gòu)建。將河流數(shù)據(jù)、斷面數(shù)據(jù)、流域邊界、模型參數(shù)導(dǎo)入MIKE 11 水動(dòng)力模型，構(gòu)建河流水文模型（MIKE-R）。

（3）流域洪澇模型（MIKE-A-R）的構(gòu)建。在MIKE 11 的MIKE SHE LINKS 依次構(gòu)建干流與支流模型，完善相關(guān)信息，并設(shè)定河床透水系數(shù)，在MIKE SHE 的河流與湖泊模塊添加MIKE-R 模型文件，調(diào)整模擬時(shí)段使兩者保持一致。通過(guò)模塊設(shè)定完成MIKE-A 模型與MIKE-R 模型的耦合［12］，建立流域洪澇MIKE-A-R耦合模擬模型。

（4）模型邊界條件與初始條件確定。確定模型上、下游邊界條件分別為聊城水文站的實(shí)測(cè)洪水過(guò)程以及宮家攔河閘水文站的實(shí)測(cè)洪水過(guò)程。

以初始水深作為模型初始條件。由于降雨之前，流域未發(fā)生洪澇災(zāi)害，因此，河流外初始水深設(shè)置為0。河流水深根據(jù)河流水文站（包括上游聊城水文站、中游劉橋攔河閘水文站以及下游宮家攔河閘水文站）實(shí)測(cè)水深確定。

2.4 模型參數(shù)敏感性及率定

在模擬過(guò)程中，可將模型參數(shù)劃分為兩類(lèi)：

（1）由實(shí)測(cè)數(shù)據(jù)或遙感數(shù)據(jù)推算得到的參數(shù)。是基于實(shí)際情況的真實(shí)描述，無(wú)需率定，主要包括：降雨量、蒸散發(fā)量、土地利用類(lèi)型、土壤類(lèi)型、植被參數(shù)等，模擬過(guò)程中直接采用實(shí)測(cè)值或通過(guò)遙感數(shù)據(jù)提取值。

（2）通過(guò)非測(cè)量和計(jì)算得到的參數(shù)。需根據(jù)模型給出的理論參考范圍和經(jīng)驗(yàn)進(jìn)行率定，主要包括土壤水動(dòng)力學(xué)參數(shù)、坡面流參數(shù)、飽和帶參數(shù)、河流參數(shù)等。

對(duì)第2類(lèi)參數(shù)進(jìn)行敏感性分析，即在參數(shù)取值范圍內(nèi)，通過(guò)分析單個(gè)參數(shù)及其相應(yīng)組合的變化對(duì)輸出結(jié)果的影響，根據(jù)敏感性對(duì)第2類(lèi)參數(shù)進(jìn)一步分類(lèi)，敏感性較低的參數(shù)為：土壤水動(dòng)力學(xué)參數(shù)及坡面匯流參數(shù)，敏感性較高的參數(shù)為：河流曼寧系數(shù)、河床透水系數(shù)及壤中流水庫(kù)模塊的相關(guān)各項(xiàng)參數(shù)。

（1）土壤水動(dòng)力學(xué)參數(shù)確定。通過(guò)SPAW（土壤水文特性）軟件計(jì)算流域土壤水動(dòng)力學(xué)參數(shù)，主要包括田間含水率、飽和導(dǎo)水率、飽和含水率以及凋零含水率，詳見(jiàn)表3。

表3 土壤類(lèi)型及其水動(dòng)力學(xué)參數(shù)表Tab.3 Soil types and hydrodynamic parameters

（2）坡面流參數(shù)確定。用曼寧系數(shù)表示坡面流參數(shù)有效糙率，根據(jù)流域下墊面土地利用類(lèi)型以及參考相關(guān)文獻(xiàn)［13］進(jìn)行取值，詳見(jiàn)表4。

表4 坡面流有效糙率參數(shù)表Tab.4 Effective roughness of slope flow

將敏感性較高的參數(shù)作為“自由”參數(shù)，進(jìn)行參數(shù)率定，率定洪水場(chǎng)次為2010-08-10、2010-08-14、2010-08-22 及2010-08-25。

參數(shù)率定目標(biāo)為洪峰流量、洪量誤差盡可能小，流量過(guò)程與觀(guān)測(cè)流量過(guò)程盡可能吻合。因此，選取洪峰流量、時(shí)段洪水總量的相對(duì)誤差、相關(guān)性系數(shù)以及納什效率系數(shù)，共同檢驗(yàn)?zāi)Ｐ途取?/p>

（1）相關(guān)性系數(shù)。

（2）納什效率系數(shù)。

式中：Qobs，i為第i斷面處的實(shí)測(cè)流量值；Qsim，i為第i斷面處的模擬流量值；為實(shí)測(cè)流量平均值；為模擬流量平均值。

率定參數(shù)結(jié)果詳見(jiàn)表5。

表5 率定參數(shù)結(jié)果表Tab.5 Calibration parameter results

3 流域洪澇過(guò)程分布式模擬結(jié)果分析與討論

3.1 河流徑流模擬結(jié)果

將率定期（2010年8月）與驗(yàn)證期（2013年7-8月、2020年8月）的河流洪水徑流過(guò)程模擬值與實(shí)測(cè)值進(jìn)行對(duì)比驗(yàn)證，實(shí)測(cè)值取自劉橋攔河閘水文站。如表6與圖2所示。

圖2 劉橋攔河閘水文站處場(chǎng)次洪水模擬值與實(shí)測(cè)值對(duì)比圖Fig.2 Comparison between simulated and measured flood values at Luqiao Sluice hydrological station

表6 模型模擬結(jié)果分析Tab.6 Analysis of model results

模擬結(jié)果與實(shí)測(cè)值的對(duì)比分析表明，模擬流量過(guò)程線(xiàn)與實(shí)測(cè)流量過(guò)程線(xiàn)趨勢(shì)基本一致，二者擬合程度較高，確定性系數(shù)介于0.803～0.936，效率系數(shù)介于在0.666～0.866，洪峰流量、洪水總量的模擬值與實(shí)測(cè)值較為一致，其相對(duì)誤差分別為-6.72%～8.10%、-17.8%～18.82%，均保持在允許誤差范圍內(nèi)（20%以下）。

3.2 淹沒(méi)水深模擬結(jié)果

調(diào)查并收集研究區(qū)域內(nèi)山東省聊城市2010-08-09 場(chǎng)次暴雨和2013-07-26 場(chǎng)次暴雨的洪澇相關(guān)資料（見(jiàn)表5），對(duì)模型淹沒(méi)水深進(jìn)行驗(yàn)證。模型淹沒(méi)水深模擬結(jié)果如圖3、4 及表7所示。

圖3 研究區(qū)域2010-08-09場(chǎng)次暴雨不同時(shí)刻的淹沒(méi)水深Fig.3 Inundation depth at different times of 2010-08-09 rainstorms in the study area

圖4 研究區(qū)域2013-07-26場(chǎng)次暴雨不同時(shí)刻的淹沒(méi)水深Fig.4 Inundation depth at different times of 2013-07-26 rainstorms in the study area

表7中調(diào)查情況摘錄資料為實(shí)地調(diào)研的最大或較大水深，模擬水深為流域1 km×1 km 柵格的平均水深。同時(shí)，經(jīng)向熟悉情況的相關(guān)技術(shù)人員查證，模擬水深與實(shí)際水深較為吻合。

表7 研究區(qū)域2010-08-09與2013-07-26場(chǎng)次暴雨洪澇淹沒(méi)情況模擬結(jié)果對(duì)比分析表Tab.7 Comparison between simulation and survey data of 20130726 rainstorm in the study area

3.3 討 論

（1）數(shù)據(jù)來(lái)源分析。MIKE-A-R 耦合模型DEM 為2017年谷歌高程，土地利用數(shù)據(jù)為國(guó)家地球系統(tǒng)科學(xué)數(shù)據(jù)中心2015年土地利用數(shù)據(jù)，土壤類(lèi)型數(shù)據(jù)為世界統(tǒng)一土壤數(shù)據(jù)庫(kù)1995年土壤類(lèi)型數(shù)據(jù)，雨水情數(shù)據(jù)為研究區(qū)域2010-2020年實(shí)測(cè)雨水情數(shù)據(jù)。其中，研究區(qū)域內(nèi)建筑用地面積僅占流域總面積的10%，近年流域城鎮(zhèn)化對(duì)土地利用變化的影響較?。煌寥李?lèi)型數(shù)據(jù)隨時(shí)間變化極小，年份對(duì)數(shù)據(jù)的準(zhǔn)確性影響不大。綜上，DEM、土地利用類(lèi)型、土壤類(lèi)型數(shù)據(jù)與雨水情數(shù)據(jù)統(tǒng)一。

（2）MIKE-A-R 耦合模型適用于平原流域的洪澇過(guò)程模擬。平原區(qū)流域坡面情況較為復(fù)雜，比降較小，在遭遇暴雨時(shí)坡面流在較短時(shí)間內(nèi)不能全部匯入河流，部分會(huì)滯蓄在河流之外形成內(nèi)澇。河流外洪澇模型MIKE-A能結(jié)合分布式下墊面條件，基于圣維南方程的擴(kuò)散波逼近法，模擬平原地區(qū)的地表徑流及淹沒(méi)過(guò)程；河流水文模型MIKE-R，基于一維非恒定流圣維南方程組，能夠快速、準(zhǔn)確地模擬河流內(nèi)洪水運(yùn)動(dòng)過(guò)程。MIKE-A 模型作為MIKE-R 模型旁側(cè)入流的邊界條件，構(gòu)建的MIKE-A-R 耦合模型，能實(shí)現(xiàn)平原流域河流內(nèi)、外的洪澇過(guò)程耦合模擬。

（3）MIKE-A-R 耦合模型的適用性分析。模型構(gòu)建時(shí)，首先根據(jù)流域的特點(diǎn)與資料情況，在建模過(guò)程中，根據(jù)流域的特點(diǎn)與資料情況，基于流域柵格與邊界的正確提取，選擇MIKE SHE 中子模塊的計(jì)算方法，利用DEM 及流域相關(guān)數(shù)據(jù)，形成河網(wǎng)數(shù)據(jù)文件，在MIKE 11 中添加干支流相關(guān)信息，在MIKE SHE的河流與湖泊模塊添加MIKE11 模型文件，實(shí)現(xiàn)流域洪澇MIKE-A-R 耦合分布式模擬。因此，結(jié)合流域特點(diǎn)，MIKE-A-R耦合模型可適用于其他平原區(qū)流域。

4 結(jié) 論

基于MIKE SHE 與MIKE 11，構(gòu)建了平原區(qū)流域洪澇過(guò)程的MIKE-A-R 耦合模型，并對(duì)徒駭河宮家攔河閘上游流域的洪澇過(guò)程進(jìn)行模擬與分析，主要結(jié)論如下。

（1）采取改進(jìn)提取流域柵格與邊界的方法，使其空間集水能力加強(qiáng)，修正水流流向，河網(wǎng)得到優(yōu)化，坡面匯流與實(shí)際更相符合，從而提取出較為準(zhǔn)確的流域范圍。

（2）建立的MIKE-A-R 耦合模型，既保證了河流在流域出口處洪水過(guò)程的精準(zhǔn)度，也通過(guò)考慮流域氣候、下墊面以及模型計(jì)算參數(shù)的空間分布不均勻等影響因素對(duì)流域河流外產(chǎn)匯流過(guò)程的影響，較好地模擬平原流域的洪澇過(guò)程。

（3）建立的MIKE-A-R 耦合模型應(yīng)用于徒駭河宮家攔河閘上游流域的洪澇過(guò)程模擬，其結(jié)果具有較好的精度，作為率定與驗(yàn)證的模擬流量過(guò)程線(xiàn)與實(shí)測(cè)流量過(guò)程線(xiàn)趨勢(shì)基本一致，確定性系數(shù)介于0.803～0.936，效率系數(shù)介于0.666～0.866，洪峰流量、洪水總量的模擬值與實(shí)測(cè)值較為一致，相對(duì)誤差均在20%以下。水深分布模擬值與調(diào)查結(jié)果較一致，洪澇過(guò)程的模擬打破了行政區(qū)劃的限制，在空間分布格局上更加精細(xì)化、精準(zhǔn)化。MIKE-A-R 耦合模型能較好地模擬平原流域河流外各地區(qū)水深分布、時(shí)變過(guò)程以及河流洪水演進(jìn)過(guò)程，模擬結(jié)果可為流域洪澇災(zāi)害風(fēng)險(xiǎn)評(píng)估、預(yù)報(bào)及預(yù)警提供全面數(shù)據(jù)支撐。結(jié)合流域特點(diǎn)，MIKE-A-R耦合模型可適用于其他平原區(qū)流域。