基于Sobol方法的MIKE NAM模型參數(shù)敏感性分析

王復生，田 娟

(安徽省水利水電勘測設(shè)計研究總院有限公司,安徽 合肥 230088)

降雨徑流模擬是防洪排澇計算、流域洪水演算以及水資源保護規(guī)劃的重要研究方法。隨著數(shù)字化技術(shù)的發(fā)展，基于模型對流域內(nèi)的降雨徑流過程進行模擬計算已經(jīng)成為了一種主要趨勢。具有物理機制的降雨徑流水文模型可以綜合考慮降雨、土壤特性、下墊面條件等多方面因素，反演不同邊界條件下的降雨產(chǎn)匯流規(guī)律及空間分布特征[1- 2]，對洪水演進理論發(fā)展和工程應(yīng)用均具有重要價值。

MIKE NAM模型是一種集總式水文模型，可以用于一個流域或由多個子分區(qū)組成的較大河流及由河流、渠道構(gòu)成的復雜河系的降雨徑流模擬過程，在國內(nèi)外得到了廣泛的研究與應(yīng)用，并取得了良好的模擬效果[3- 5]。NAM模型的參數(shù)較多，模擬過程中存在模型參數(shù)的不確定性與敏感性，參數(shù)和邊界條件之間的相互影響等問題[6- 7]。參數(shù)取值的合理性以及對參數(shù)含義的準確理解對模型精度會產(chǎn)生很大影響，因此對模型參數(shù)進行全局敏感性分析，篩選出高敏感性參數(shù)進行重點率定十分有必要。

本文以巢湖流域桃溪水文站控制子流域為研究對象，建立MIKE NAM模型，綜合考慮NAM模型中各類參數(shù)的取值情況，運用Sobol方法對模型參數(shù)的敏感性進行定量分析，評估各參數(shù)對模型輸出結(jié)果方差的影響程度，將高敏感性參數(shù)作為主要率定參數(shù)，以期提高模型率定效率和模擬精度。

1 模型與方法

1.1 MIKE NAM模型

MIKE是由丹麥水資源與水環(huán)境研究所(DHI)開發(fā)的專業(yè)工程軟件，適用于河流、灌溉渠道及其它水體的一二維水動力、水質(zhì)、泥沙運輸模擬等[8- 9]。MIKE NAM模型作為一種新興的水文模擬工具，在流域降雨產(chǎn)匯流模擬中具有良好的適用性，被廣泛應(yīng)用到世界各地不同氣象水文條件的流域，在國內(nèi)長江、黃河、松花江流域均有應(yīng)用案例，是一個經(jīng)過大量工程實踐驗證的模型工具。NAM模型將土壤含水層分成積雪儲水層(Snow Storage)、地表儲水層(Surface Storage)、淺層儲水層(Lower Zone Storage)和地下水儲水層(Ground Water Storage)4個部分，模型分別進行連續(xù)計算以模擬流域中各種相應(yīng)的水文過程。

1.2 Sobol參數(shù)敏感性分析方法

Sobol方法是一種基于方差分解的定量全局敏感性分析算法，相比于定性敏感性分析方法，該方法能夠通過敏感性量化指標的計算直接給出模型參數(shù)的敏感性大小，其核心是將目標函數(shù)的總方差分解為單個參數(shù)的方差和多個參數(shù)間相互作用的方差，目前已被廣泛應(yīng)用于經(jīng)濟、環(huán)境等多個領(lǐng)域[10]。

Sobol方法假設(shè)模型結(jié)構(gòu)可以表示為函數(shù)形式u=f(x)，其中模型參數(shù)x=x1，x2，…，xn可視為n維離散點，u為標量輸出。

假設(shè)函數(shù)f(x)可積，xi在[0，1]上服從均勻分布，且f(x)滿足：

(1)

式中，k=i1，…，is。

則函數(shù)f(x)的方差分解表達式可以表示為：

(2)

式中，1≤i1<…is≤n(1≤s≤n)，則被加數(shù)共有2n項。

類似的，模型的總方差也可以分解為單個參數(shù)和多個參數(shù)相互作用的組合：

(3)

式中，D—模型的總方差；Di—參數(shù)xi產(chǎn)生的方差；Dij—參數(shù)xi和xj相互作用產(chǎn)生的方差，D1，2，…，n—n個參數(shù)共同作用產(chǎn)生的方差。

將上式歸一化即可得到各參數(shù)和參數(shù)之間相互作用的敏感性：

(4)

則模型的敏感性指數(shù)可以表示為：

式中，Si—參數(shù)xi單獨作用時的敏感度；Sij—參數(shù)xi和xj相互作用的敏感度；STi—參數(shù)xi和其它參數(shù)共同作用的敏感度；D～i—除參數(shù)xi外其它參數(shù)共同作用產(chǎn)生的方差。

2 研究區(qū)域與模型構(gòu)建

2.1 研究區(qū)域概況

本研究對巢湖流域上游豐樂河桃溪水文站控制流域進行建模分析。

豐樂河發(fā)源于大別山余脈，流經(jīng)雙河、桃溪、豐樂至三河鎮(zhèn)下大潭灣與杭埠河匯合后流入巢湖，河道全長117.5km，流域面積2124km2，是杭埠河最大的支流。豐樂河干流上的桃溪水文站設(shè)立于1951年7月，控制集水面積1510km2，約占豐樂河流域總面積的75%，有長系列的實測流量資料，且精度可靠，是巢湖流域水文分析計算重要的參證站點。豐樂河桃溪站控制流域地形圖如圖1所示。

圖1 桃溪站控制流域地形圖

根據(jù)研究區(qū)域內(nèi)實測降雨資料及雨量站的分布情況，選取9個代表站點進行降雨徑流計算，利用泰森多邊形法分析得到各個雨量站的計算權(quán)重值見表1。

表1 桃溪站控制流域面雨量計算各雨量站權(quán)重

2.2 NAM模型相關(guān)參數(shù)

NAM模型所需數(shù)據(jù)資料主要為流域降雨及蒸發(fā)數(shù)據(jù)，模型的初始條件包括開始時刻流域地表儲水層和根區(qū)儲水層的土壤含水率，以及坡面流、壤中流和基流的初始值[11]。NAM模型的輸出結(jié)果主要包括地表坡面徑流、壤中流和基流過程等。模型主要參數(shù)及取值范圍詳見表2。作為集總式水文模型，NAM模型中每個子流域的參數(shù)均為本區(qū)域的參數(shù)平均值。

表2 NAM模型主要參數(shù)

3 計算結(jié)果與分析

3.1 參數(shù)敏感性分析

由于Umax與Lmax相關(guān)性較高，可視為一個參數(shù)，本文以Umax、CQOF、CKIF、CK1、2、TOF、TIF、TG、CKBF共8個模型參數(shù)作為自變量函數(shù)，用蒙特卡羅法(Monte Carlo Simulation)對自變量進行抽樣，每個參數(shù)得到10組樣本值，根據(jù)Sobol參數(shù)敏感性分析方法產(chǎn)生10×(8+2)=100個參數(shù)組合，將這100組參數(shù)分別帶入桃溪站控制流域NAM模型中，對2020年6月8日—8月8日的區(qū)域降雨進行模擬。分別以洪峰流量、峰現(xiàn)時間和最大7d洪量作為NAM模型輸出結(jié)果的目標函數(shù)，用一階敏感性指數(shù)Si和全局敏感性指數(shù)STi來量化輸出結(jié)果對各個參數(shù)的敏感度高低。模型參數(shù)敏感性分析成果如圖2所示。

圖2 參數(shù)敏感性分析成果圖

由敏感性分析成果可知：

以洪峰流量作為目標函數(shù)時，坡面流匯流系數(shù)CQOF和坡面流匯流時間常數(shù)CK1、2的一階敏感性指數(shù)、全局敏感性指數(shù)遠高于其它參數(shù)，說明這兩個參數(shù)對洪峰流量大小起決定性作用。CQOF的一階敏感性指數(shù)大于0，CK1、2的一階敏感性指數(shù)小于0，說明坡面流匯流系數(shù)與流量呈正相關(guān)，坡面流匯流時間常數(shù)與流量呈負相關(guān)，這與參數(shù)代表的物理意義相符。

以峰現(xiàn)時間作為目標函數(shù)時，坡面流匯流時間常數(shù)CK1、2的敏感性指數(shù)高于其它參數(shù)，說明該參數(shù)對洪峰出現(xiàn)的時間有較大影響。

以最大7d洪量作為目標函數(shù)時，坡面流匯流系數(shù)CQOF的敏感性指數(shù)高于其它參數(shù)，說明該參數(shù)對洪量有較大影響。

3.2 模型率定

根據(jù)敏感性分析結(jié)論，NAM模型模擬結(jié)果對不同參數(shù)變化的敏感性是有差異的，參數(shù)CQOF和CK1、2對NAM模型精度的影響比其它參數(shù)大，率定時需重點考慮。因此本研究在對2020年6月8日—8月8日降雨進行模擬時，先采用NAM模型自帶的SCE全局尋優(yōu)算法對8個參數(shù)進行自動率定，迭代次數(shù)設(shè)置為1000次，在此基礎(chǔ)上對CQOF和CK1、2兩個重點參數(shù)進行人工率定，最終得到的參數(shù)率定成果如表3所示。桃溪站模擬洪水與實測洪水擬合效果見圖3，洪量的擬合效果見圖4。

表3 桃溪站參數(shù)率定成果

圖3 桃溪站模擬洪水與實測洪水擬合效果

圖4 桃溪站累計洪量擬合效果

2020年桃溪站實測洪峰流量為1239m3/s，模型計算得到的洪峰流量為1304m3/s，差值為5.2%，原因在于2020年桃溪站上游柏林圩等圩口發(fā)生了潰破，對豐樂河洪峰流量起到了削峰作用。2020年桃溪站實測最大7d洪量為4.26億m3/s，模型計算最大7d洪量為4.41億m3，誤差為3.5%，桃溪站實測洪峰與模型模擬洪峰的出現(xiàn)時間相差10h，誤差在合理范圍內(nèi)?？紤]破圩還原后，誤差進一步減小，計算誤差結(jié)果見表4。

表4 計算誤差統(tǒng)計表

4 結(jié)論

(1)基于Sobol方法對MIKE NAM模型中的主要參數(shù)進行全局敏感性分析，分析結(jié)果表明：參數(shù)CQOF和CK1、2為高敏感度參數(shù)，對模型精度影響較大，建議作為重點率定參數(shù)。

(2)對CQOF和CK1、2兩個參數(shù)采用人工率定結(jié)果，其它參數(shù)采用模型自動率定結(jié)果，以豐樂河桃溪水文站控制流域為研究區(qū)，模擬了2020年汛期降雨作用下的區(qū)域降雨徑流過程，并根據(jù)實測流量數(shù)據(jù)對模擬效果進行擬合檢驗，結(jié)果表明模型計算誤差在合理范圍內(nèi)。

(3)本文構(gòu)建的NAM模型中未考慮流域內(nèi)圩區(qū)的調(diào)蓄以及流域下游水位實時變化的作用，對模擬精度有一定的影響，并且參數(shù)敏感性分析成果缺少不同區(qū)域不同場次降雨的驗證。后續(xù)可以將NAM模型與MIKE11水動力模型進行耦合，針對不同降雨強度、不同流域特征進行建模和率定，進而實現(xiàn)巢湖全流域的降雨徑流模擬。本研究可為巢湖流域水文水動力模型的構(gòu)建和驗證提供技術(shù)支撐和應(yīng)用基礎(chǔ)。