劉嘉男，李發(fā)文※，馮 平，王現(xiàn)領(lǐng)

(1.天津大學(xué)水利工程仿真與安全國家重點實驗室，天津 300072； 2.天津市水利科學(xué)研究院，天津 300061)

0 引言

21世紀(jì)以來，全球發(fā)生洪澇災(zāi)害的頻率高于過去[1]。隨著全球環(huán)境的變化，這種趨勢越加明顯。研究表明， 較長時間尺度上，氣候變化對水文水資源的影響更加明顯，但是短期內(nèi)，土地利用/覆被變化(LUCC)是水文水資源變化的主要驅(qū)動要素之一[2]。隨著經(jīng)濟(jì)的發(fā)展，農(nóng)業(yè)科學(xué)技術(shù)水平的提高，設(shè)施農(nóng)業(yè)迅速發(fā)展，對農(nóng)業(yè)產(chǎn)業(yè)結(jié)構(gòu)調(diào)整和農(nóng)民增收起到巨大作用。但是設(shè)施農(nóng)業(yè)用地增加導(dǎo)致了不透水面積增大，使地表徑流系數(shù)增大，進(jìn)而導(dǎo)致農(nóng)田區(qū)域產(chǎn)水量增大，對農(nóng)田區(qū)域的排澇模數(shù)造成一定影響。而排澇模數(shù)對合理確定排澇系統(tǒng)規(guī)模、保障農(nóng)田安全具有重要現(xiàn)實意義。因此，在設(shè)施農(nóng)業(yè)得到快速發(fā)展的情況下，有必要開展設(shè)施農(nóng)業(yè)對農(nóng)田排澇模數(shù)影響的研究。

排澇模數(shù)的確定主要依賴于設(shè)計徑流深，因此建立合理的降雨徑流模型是進(jìn)行排澇模數(shù)估算的基礎(chǔ)。Stanford模型、Tank模型、Boughton模型、新安江模型、陸渾模型、SWAT模型等都是流域產(chǎn)流計算比較成熟的模型。產(chǎn)流模型主要應(yīng)用于大尺度流域，對山區(qū)模擬效果較好，而對平原區(qū)產(chǎn)流模擬精度一般較差。平原區(qū)產(chǎn)流機(jī)制更為復(fù)雜，下墊面多樣且變化，超滲產(chǎn)流與蓄滿產(chǎn)流模式并存，采用單一產(chǎn)流模型時，往往不能更全面、更真實地反映平原區(qū)的產(chǎn)流機(jī)制，導(dǎo)致模型精度達(dá)不到使用要求。文章選用包為民[3]提出的超滲—蓄滿垂向耦合產(chǎn)流模型，該模型很好地解決了平原區(qū)超滲與蓄滿產(chǎn)流模式并存問題，在平原區(qū)得到了廣泛應(yīng)用[4]?？紤]到農(nóng)田區(qū)域?qū)崪y徑流資料往往需要反演推算，資料的準(zhǔn)確度不高，而農(nóng)田可以獲得實測土壤含水量和地下水位資料，因此在模型合理性評價中，除了采用實測徑流相對誤差作為模型的評價標(biāo)準(zhǔn)外，還采用實測土壤含水量以及實測地下水位的效率系數(shù)和相對誤差作為模型的評價標(biāo)準(zhǔn)，以有效提高模型模擬的精度。

模型參數(shù)敏感性分析和參數(shù)優(yōu)化對模型模擬結(jié)果有重要影響[5-6]。該文選用最具代表性的全局參數(shù)敏感性分析方法Sobol′法[7]，該方法已經(jīng)廣泛應(yīng)用于眾多水文模型參數(shù)敏感性分析中[8-10]。在模型參數(shù)優(yōu)化選擇中，SCE-UA算法也得到了廣泛應(yīng)用[11]，王維[4]首次將SCE-UA用于垂向耦合模型中，并與新安江模型進(jìn)行對比，證明了SCE-UA算法在兩種模型中的有效性。

該文選擇具有一定設(shè)施農(nóng)業(yè)規(guī)模的封閉小區(qū)為研究區(qū)，構(gòu)建超滲—蓄滿垂向耦合模型，基于Latin超立方抽樣，進(jìn)行Sobol′全局敏感性分析以及SCE-UA全局參數(shù)優(yōu)化，同時考慮土壤含水量與地下水位變化量，與人工率定相結(jié)合，進(jìn)行模型參數(shù)的優(yōu)選。針對不同重現(xiàn)期的設(shè)計暴雨，通過已經(jīng)構(gòu)建好的模型求得農(nóng)田小區(qū)的產(chǎn)流量，用平均排除法計算不同重現(xiàn)期下的設(shè)計排澇模數(shù)，與2003年《天津市平原地區(qū)農(nóng)田除澇水文手冊》(簡稱《除澇手冊》)進(jìn)行對比，探討設(shè)施農(nóng)業(yè)下對原除澇手冊修訂的必要性。

圖1 研究區(qū)概況及站點布設(shè)

1 材料與方法

1.1 數(shù)據(jù)獲取與處理

研究區(qū)域位于天津市西青區(qū)張家窩鎮(zhèn)，屬于半濕潤半干旱區(qū)。建立產(chǎn)流模型所需要的數(shù)據(jù)均為現(xiàn)場試驗實測資料，在試驗區(qū)內(nèi)設(shè)置雨量計、土壤墑情測量儀、地下水位觀測儀、泵站排水流量計等觀測設(shè)備，可得到降雨量、土壤含水量、地下水位變化量以及泵站排水量數(shù)據(jù)系列，為模型建立提供了可靠的數(shù)據(jù)來源。設(shè)備具體型號為：自動雨量計：LC-YL1、土壤水分傳感器：FDS120、超聲波流量計：SCL-80、壓阻式水位計：Unisens-WL20。研究區(qū)域概況及站點布設(shè)見圖1。

通過2016～2017年的現(xiàn)場試驗，采集到8場降雨過程。對觀測的降雨量、土壤含水量、地下水位相對變化量、泵站排水量數(shù)據(jù)以1h為計算時段進(jìn)行整理，并將土壤含水量與地下水位相對變化量轉(zhuǎn)化為水層深度(mm)，方便建模使用。

經(jīng)調(diào)研，該農(nóng)田小區(qū)在2008年之前全部為經(jīng)濟(jì)林木，隨著天津“菜籃子”工程的實施，設(shè)施農(nóng)業(yè)迅猛發(fā)展[12]，占地面積為3.709 8萬m2， 占總面積的7%； 2011年底，設(shè)施農(nóng)業(yè)占地面積進(jìn)一步擴(kuò)大，占地面積為7.753 1萬m2，占總面積的15%； 2012年底至今，設(shè)施農(nóng)業(yè)面積所占比例達(dá)到41%，研究區(qū)域設(shè)施農(nóng)業(yè)發(fā)展概況見圖2。

圖2 研究區(qū)域設(shè)施農(nóng)業(yè)發(fā)展概況

1.2 超滲—蓄滿耦合產(chǎn)流模型建立

1.2.1 透水面(農(nóng)田)

降雨P(guān)扣除蒸散發(fā)E后的凈雨量PE到達(dá)地面，首先通過空間分布的下滲曲線，判斷是否發(fā)生超滲產(chǎn)流，當(dāng)發(fā)生超滲產(chǎn)流時，采用具有流域分布特征的格林—安普特下滲曲線計算超滲產(chǎn)流量和下滲量[13]，然后以超滲產(chǎn)流的下滲量作為蓄滿產(chǎn)流中的降雨輸入量。在土壤缺水量大的面積上，補(bǔ)充土壤含水量不產(chǎn)流，在缺水量小的面積上，補(bǔ)足土壤缺水量后，發(fā)生蓄滿產(chǎn)流[3]。

超滲產(chǎn)流量計算過程為[3]：

PE=P-E-I

(1)

(2)

FMM=FM×(1+B1)

(3)

(4)

R1=PE-FA

(5)

圖3 模型垂向耦合形式

其中，P為降雨量，mm；E為蒸散發(fā)量，mm；I為植物截留量，mm；PE為凈雨量，mm；FM為農(nóng)田小區(qū)平均下滲能力，mm/h；FC1為穩(wěn)定下滲率，mm/h；K2為土壤缺水量對下滲率影響的靈敏度系數(shù)；FMM為平均下滲率為FM時上層最大點的下滲率，mm/h；FA為實際下滲量，mm；WM為土壤平均蓄水容量，mm；W為土壤含水量，mm；B1為下滲能力分配曲線的指數(shù)，反映了下滲能力分布的不均勻性；R1為超滲產(chǎn)流量，mm。

地面以下的徑流，采用蓄滿產(chǎn)流結(jié)構(gòu)，計算公式為：

(6)

WMM=WM×(1+B2)

(7)

(8)

式中，A1為前期影響雨量，mm；WMM為土壤內(nèi)最大點蓄水容量，mm；B2為蓄水容量—面積分配曲線的指數(shù)，反映了蓄水容量分布的不均勻性；R2為蓄滿產(chǎn)流量，mm； 式中其他符號含義同前。

總產(chǎn)流量R為超滲產(chǎn)流量R1與蓄滿產(chǎn)流量R2之和，即：

R=R1+R2

(9)

蓄滿產(chǎn)流計算出的產(chǎn)流量R2包括地面徑流R2S和地下徑流R2G兩部分，可以用穩(wěn)定下滲率進(jìn)一步劃分。計算過程為：

(10)

其中，F(xiàn)C2為下層土壤穩(wěn)定下滲能力，mm/h。

蒸散發(fā)量與土壤含水量的計算采用新安江模型中描述的三層蒸散發(fā)模式[14]，植物截留量對1次降雨徑流過程的影響量很小，但對年降水量而言，森林及茂密的植被，1年中截留的水量可達(dá)年降水量的25%～30%，這些截留的水量最終消耗于蒸散發(fā)，因此該文中將植物的截留作用合并于降雨的蒸散發(fā)損失中。這樣就構(gòu)成了完整的超滲—蓄滿垂向耦合產(chǎn)流模型。

1.2.2 不透水面(設(shè)施農(nóng)業(yè))

在不透水面積上的產(chǎn)流量為直接徑流DRS，其計算式為：

DRS=P×IMP×α0

(11)

α0=0.95，RS>0

(12)

其中，IMP為不透水面所占面積比例，經(jīng)測量取值為0.41；α0為不透水面產(chǎn)流系數(shù)；RS為透水面地面總徑流，即R1與R2S之和。

1.2.3 總產(chǎn)流計算

透水面積與不透水面積上的產(chǎn)流一起匯入河道，經(jīng)泵站排出，泵站出水口處產(chǎn)流深計算公式為：

R總=(1-IMP)×R+α0×IMP×P

(13)

其中，R總為區(qū)域內(nèi)總徑流模擬值，mm。

綜上所述，在研究區(qū)域透水面上采用超滲—蓄滿垂向耦合結(jié)構(gòu)，地面徑流R1按超滲產(chǎn)流計算，地面以下徑流按蓄滿產(chǎn)流計算，并依據(jù)下層土壤下滲能力分為地面徑流R2S和地下徑流R2G，其中，R1與R2S以及不透水面上的產(chǎn)流量DRS匯入河道，經(jīng)泵站排出。蓄滿產(chǎn)流中的地下徑流R2G補(bǔ)充地下水位。超滲—蓄滿垂向耦合模型流程圖見圖4。

圖4 超滲—蓄滿耦合產(chǎn)流模型流程圖

1.3 模型參數(shù)敏感性分析及率定

通過Latin超立方抽樣隨機(jī)生成符合要求的參數(shù)組，選用Sobol′敏感度分析方法進(jìn)行全局參數(shù)的一階敏感度與總敏感度計算，采用SCE-UA算法進(jìn)行全局參數(shù)自動率定，同時結(jié)合敏感性分析的成果與模型適用性評價標(biāo)準(zhǔn)進(jìn)行人工率定，最終得到模型的參數(shù)。

1.3.1 Latin超立方抽樣

Latin超立方抽樣作為Monte-Carlo方法的改進(jìn)，提供了1個有效而實用的小樣本采樣技術(shù)，在具有隨機(jī)輸入變量的復(fù)雜模型分析中應(yīng)用廣泛。該抽樣方法的本質(zhì)是將N維單位立方體中的每一維坐標(biāo)區(qū)間[0， 1]等分成n個等間距不重疊的子區(qū)間，對每個子區(qū)間分別進(jìn)行獨立的等概率抽樣，這樣每個樣本點都屬于特定的分層區(qū)間，然后將抽到的每一分量進(jìn)行隨機(jī)組合，從而可以使樣本點更均勻地分布在參數(shù)的可行域內(nèi)，通過較少的采樣得到相同的效果[15]，根據(jù)Tang等[16]的研究成果，選用Latin超立方技術(shù)對參數(shù)進(jìn)行抽樣，實現(xiàn)對模型參數(shù)的全局敏感性分析以及參數(shù)優(yōu)化。

1.3.2 Sobol′全局敏感度分析法

Sobol′方法是一種基于方差分解的全局模型敏感度分析方法，可以用于非線性和非單調(diào)函數(shù)和模型[17]，該方法的核心是將所研究的模型分解為單參數(shù)及多參數(shù)相互組合的函數(shù)，分別求出每項對應(yīng)于模型輸出的敏感度系數(shù)。假設(shè)產(chǎn)流模型為Y=f(X)，X=(x1，x2，…，xm)，其中，Y為模型輸出的目標(biāo)函數(shù)，X=(x1，x2，…，xm)為影響產(chǎn)流模型的m個待分析參數(shù)，X服從[0， 1]上的均勻分布，且f2(x)可積，定義Θm為模型輸入?yún)?shù)的空間域，表示為Θm=(X|0≤xi≤1，i=1， 2，…，m)，將模型分解為2m項之和[18]：

(14)

其中，f0為常數(shù)項，保證模型分解的形式具有唯一性，要求其余子項對其所包含的變量的積分為0，即滿足公式：

(15)

則Y的方差D(Y)可以作分解：

(16)

其中，Di為參數(shù)ii產(chǎn)生的方差；Dij為兩個參數(shù)和j相互作用產(chǎn)生的方差；Dijk為3個參數(shù)ijk相互作用產(chǎn)生的方差；D12…m為m個參數(shù)共同作用產(chǎn)生的方差。將式(16)歸一化后得到各參數(shù)和參數(shù)相互作用的敏感性：

(17)

產(chǎn)流模型參數(shù)i對目標(biāo)函數(shù)的敏感度值可表示為：

(18)

(19)

其中，一階敏感度Si為參數(shù)i作用的敏感度； 總敏感度STi為參數(shù)i單獨作用以及與其他參數(shù)相互作用的敏感度；D～i為除了參數(shù)i之外的參數(shù)的方差。

由于水文模型具有很大的非線性和復(fù)雜性，運用Sobol′法進(jìn)行敏感度計算時，對產(chǎn)流模型中待分析的變量進(jìn)行兩次獨立抽樣，得到ON×m和PN×m兩個矩陣，其中，N為抽樣次數(shù)，m為參數(shù)個數(shù)，即：

(20)

(21)

(22)

(23)

1.3.3 SCE-UA算法

SCE-UA(Shuffled Complex Evolution Algorithm)算法是由美國亞利桑那州大學(xué)Duan等[11]于1992年提出來的，其綜合了隨機(jī)搜索算法、單純形法、聚類分析法及生物競爭演化等方法的優(yōu)點，收斂速度較快、魯棒性強(qiáng)、穩(wěn)定性好，可以高效解決高維非線性參數(shù)的全局優(yōu)化問題，且不需要顯式的目標(biāo)函數(shù)或者目標(biāo)函數(shù)的偏導(dǎo)數(shù)。該算法的詳細(xì)步驟如圖5。

圖5 SCE-UA算法流程圖

1.4 模型適用性評價

產(chǎn)流計算的準(zhǔn)確率與可信度是進(jìn)行排澇模數(shù)計算的保障，根據(jù)《水文情報預(yù)報規(guī)范》(SD 138-85)[19]，水量平衡評價指標(biāo)選取徑流深相對誤差，地下水位變化量與土壤含水量變化量選取相對誤差Re和Nash效率系數(shù)Ens評價模型模擬精度。計算公式為：

(24)

(25)

為了進(jìn)一步驗證上述產(chǎn)流計算模型，分別采用單一的超滲產(chǎn)流模型、蓄滿產(chǎn)流模型以及該文中建立的超滲—蓄滿垂向耦合模型計算，計算模型相對有效性[3]。計算公式：

(26)

(27)

其中，CEFc代表耦合產(chǎn)流相對于超滲產(chǎn)流的有效性；CEFx代表耦合產(chǎn)流相對于蓄滿產(chǎn)流的有效性；R0為實測降雨徑流深，mm；Rc、Rx、RH分別為用超滲、蓄滿和耦合模型計算的降雨徑流深，mm。

1.5 排澇模數(shù)計算

采用平均排除法計算排澇模數(shù)，即將相應(yīng)不同重現(xiàn)期下的降雨產(chǎn)流量在規(guī)定的排澇天數(shù)內(nèi)均勻排除[20]。該方法概念明確、所需資料較少，計算方便。計算公式為：

(28)

其中，M為設(shè)計排澇模數(shù)，m3/(s·km2)；R為設(shè)計徑流深，mm；td為排澇天數(shù)，d。天津市的排澇標(biāo)準(zhǔn)為1日暴雨2日排除[21]，即td=2。

表1 超滲—蓄滿垂向耦合模型參數(shù)含義及取值范圍

序號參數(shù)參數(shù)含義上限下限1K1陸面蒸散發(fā)能力折算系數(shù)050952C深層蒸散發(fā)系數(shù)01023FC1上層穩(wěn)定下滲率4654K2土壤缺水量靈敏度系數(shù)15355WUM上層土壤蓄水容量20506WLM下層土壤蓄水容量40707WDM深層土壤蓄水容量65908B1土壤下滲能力分配曲線指數(shù)03089B2蓄水容量—面積分配曲線指數(shù)020610FC2下層土壤穩(wěn)定下滲率254

2 實例應(yīng)用

2.1 研究區(qū)域水文模型構(gòu)建

按照1.2所述過程構(gòu)建了該農(nóng)田小區(qū)超滲—蓄滿垂向耦合模型，模型共有10個參數(shù)。利用兩年現(xiàn)場試驗所經(jīng)歷的8場降雨事件，進(jìn)行模型參數(shù)的敏感性分析以及參數(shù)率定與驗證，以得到優(yōu)選的模型參數(shù)。

2.2 基于Sobol′的參數(shù)敏感性分析

因為沒有關(guān)于參數(shù)的先驗信息，進(jìn)行敏感性分析輸入的參數(shù)值是從均勻分布中抽取的，不同的參數(shù)范圍在0～1之間進(jìn)行線性變換[22]。產(chǎn)流模型中參數(shù)的取值上下限見表1。利用Latin超立方抽樣，在參數(shù)的取值范圍內(nèi)隨機(jī)生成包含2 000組參數(shù)的兩個樣本[3， 14， 23]。

模型參數(shù)的敏感性與目標(biāo)函數(shù)相關(guān)聯(lián)，該文選取徑流深與徑流深相對誤差兩個指標(biāo)作為目標(biāo)函數(shù)。根據(jù)雨量等級將率定期5作為目標(biāo)函數(shù)。根據(jù)雨量等級將率定期5場降雨進(jìn)行劃分，最終選取3場降雨(其中160720場次降雨221.8mm，為特大暴雨； 160807場次降雨92mm，為暴雨； 160724場次降雨41.2mm，為大雨)進(jìn)行各參數(shù)一階敏感度與總敏感度的分析。

3種典型降雨下10個參數(shù)對2個目標(biāo)函數(shù)的一階敏感度與總敏感度值分布情況見圖6。每一列代表同一目標(biāo)函數(shù)下模型參數(shù)在3種典型降雨中的表現(xiàn)。該文中指定總敏感度超過閾值10%(圖6中虛線)的參數(shù)為敏感參數(shù)。

圖6 3種典型降雨下10個參數(shù)對不同目標(biāo)函數(shù)的一階敏感度與總敏感度圖注:每幅圖中的左邊圖柱為一階敏感度；右邊圖柱為總敏感度

通過分析圖6，發(fā)現(xiàn)當(dāng)以徑流深和徑流深相對誤差為評價指標(biāo)時，參數(shù)在相同等級的降雨中所呈現(xiàn)的規(guī)律是一致的。在特大暴雨時，一階敏感度最大的參數(shù)是FC2，較敏感的參數(shù)是FC1，總敏感度最大的是FC2，其次為K2； 在暴雨中，一階敏感度與總敏感度最大的參數(shù)均是WUM，WLM、WDM、FC2、K2、FC1均為敏感參數(shù)，K1為較敏感參數(shù)； 在大雨中，一階敏感度與總敏感度最大的參數(shù)依然均是WUM，WLM、WDM、FC2，F(xiàn)C1為敏感參數(shù)。

2.3 基于SCE-UA的自動率定

在進(jìn)行參數(shù)優(yōu)化時，將率定期5場降雨綜合考慮，以徑流深相對誤差為核心，采用表2目標(biāo)函數(shù)，通過SCE-UA算法，自動尋找到使目標(biāo)函數(shù)值達(dá)到最小，且滿足單場次降雨徑流深相對誤差在20%范圍內(nèi)的最優(yōu)參數(shù)組。

SCE-UA算法自身包含多個參數(shù)，根據(jù)前人研究成果[11， 24]，該文中算法參數(shù)取值見表3，參數(shù)優(yōu)化的可行域與表1相同。

表2 SCE-UA算法目標(biāo)函數(shù)

評價目標(biāo)目標(biāo)函數(shù)平均相對誤差MAPEMAPE=1N∑Ni=1Ri，obs-Ri，simRi，obs最大相對誤差MMPEMMPE=maxRi，obs-Ri，simRi，obs()平均相對對數(shù)絕對誤差MLGMLG=1N∑ni=1logRi，obs-Ri，simRi，obs+1[]總目標(biāo)函數(shù)Fmin(F)=∑MLG+MAPE+MMPE() 注：表2中，Ri，obs為場次降雨i的實際徑流深，Ri，sim是場次降雨i的模擬徑流深；i=1，2，…，N；N為場次降雨總數(shù)

表3 SCE-UA算法參數(shù)設(shè)置

參數(shù)名稱參數(shù)意義參數(shù)取值n待優(yōu)化參數(shù)的個數(shù)10p復(fù)合形個數(shù)2m每個復(fù)合形的頂點數(shù)21s樣本點數(shù)42q每個子復(fù)合形的頂點數(shù)11α父代產(chǎn)生子代的個數(shù)1β父代產(chǎn)生子代的代數(shù)21T最大循環(huán)次數(shù)100/200/1000/2000

通過Latin超立方抽取42組參數(shù)，以表2總目標(biāo)函數(shù)對率定期5場降雨綜合考慮，最大迭代次數(shù)T分別設(shè)置為100/200/1 000/2 000，并分別連續(xù)運行5次，初步率定出一組較優(yōu)參數(shù)，結(jié)果見表4。此時，總目標(biāo)函數(shù)值以及單場次徑流深相對誤差值見表5。各場次徑流深相對誤差均在20%內(nèi)，滿足要求。

表4 參數(shù)初步率定結(jié)果

參數(shù)K1CFC1K2WUMWLMWDMB1B2FC2參數(shù)值055013643307453568417801031021254

表5 目標(biāo)函數(shù)值及單場次徑流深相對誤差

表6 參數(shù)率定結(jié)果

2.4 自動率定與人工率定相結(jié)合

在SCE-UA參數(shù)優(yōu)選過程中，蒸散發(fā)計算與兩水源劃分參數(shù)，即K1、C、WUM、WLM、FC2波動比較大； 產(chǎn)流計算參數(shù)FC1、WDM相對比較收斂； 產(chǎn)流計算參數(shù)K2、B1、B2較穩(wěn)定。其中K1、C、B1、B2為不敏感參數(shù)，在率定過程中，可以對K1取多次尋優(yōu)結(jié)果的平均值0.75； C根據(jù)經(jīng)驗可以取0.17； B1、B2取自動率定結(jié)果值，分別為0.3和0.2。其他敏感參數(shù)的率定還需要結(jié)合產(chǎn)流區(qū)域?qū)嶋H情況進(jìn)一步分析。

鑒于產(chǎn)流區(qū)域?qū)崪y降雨場次較少，沒有實測流量過程線，缺少對產(chǎn)流過程的模擬，因此在參數(shù)率定過程中，不僅以徑流深相對誤差為評價指標(biāo)，還需考慮地下水位變化與土壤含水量變化過程，提高模型的模擬精度。地下水位主要受地下徑流的影響，著重調(diào)節(jié)FC1與FC2使地下水位變化過程滿足要求； K2、WUM、WLM、WDM與土壤含水量有密切關(guān)系，調(diào)節(jié)這4個參數(shù)使土壤含水量評價指標(biāo)滿足要求。自動率定與人工率定結(jié)合得到的參數(shù)數(shù)值見表6。

2.5 模型適用性評價

按照1.4計算式(24)(25)計算，率定期與驗證期各項評價指標(biāo)結(jié)果見表7，按照式(26)(27)選取8場降雨計算垂向耦合模型相較于單一超滲模型和單一蓄滿模型的有效性，結(jié)果見表8。

表7 模型評價指標(biāo)分析

模擬期降雨場次評價對象Nash系數(shù)相對誤差(%)率定期160720地下水位變化量0671564土壤含水量079491徑流深相對誤差(%)512160724地下水位變化量0741257土壤含水量061561徑流深相對誤差(%)562160801地下水位變化量0681437土壤含水量074605徑流深相對誤差(%)254160807地下水位變化量0781324土壤含水量080811徑流深相對誤差(%)1210160818地下水位變化量0751263土壤含水量073417徑流深相對誤差(%)1067驗證期170706地下水位變化量0721058土壤含水量081763徑流深相對誤差(%)271170709地下水位變化量0681239土壤含水量071834徑流深相對誤差(%)1593170722地下水位變化量0801161土壤含水量087637徑流深相對誤差(%)510

表8 產(chǎn)流結(jié)果計算成果

降雨場次降雨量(mm)實測流量(m3)R0(mm)RC(mm)R0?RCR0(%)RX(mm)R0?RXR0(%)RH(mm)R0?RHR0(%)1607202218962561899815438187423320-2275178905831607243823380667635479746-117963055516080141248809639610271199-2451940239160807922850056256323-1240501810796305-120916081835232806475641288732-1305578106617070649449809831016-3371333-35621010-27617070933448609591140-18851316-37201112-159317072247849609799433671305-3331929510匯總CEFc=0869CEFx=0912

由表7可知，率定期與驗證期地下水變化量與土壤含水量變化量的Nash效率系數(shù)均大于0.6，各項相對誤差均控制在20%之內(nèi)，可以認(rèn)為模型是合理的。

由表8可知，垂向耦合模型相對于超滲產(chǎn)流的有效性為0.869，相對于蓄滿產(chǎn)流的有效性為0.912，進(jìn)一步證明了在該農(nóng)田小區(qū)采用超滲蓄滿垂向耦合產(chǎn)流模型是合理的。

2.6 排澇模數(shù)計算

2.6.1 降雨趨勢分析

利用5年滑動平均值對1990～2015年年降雨量及汛期(6～9月)降雨量進(jìn)行趨勢分析[25]，結(jié)果見圖7。

圖7 降雨趨勢分析

通過圖7(a)與(b)發(fā)現(xiàn)，1990～2015年，年降雨量以及汛期降雨量波動平穩(wěn)，沒有顯著線性變化趨勢。因此，分析設(shè)施農(nóng)業(yè)對農(nóng)田排澇模數(shù)影響時，可以不考慮降雨變化對其影響。

2.6.2 設(shè)計暴雨計算

采用2015年《天津市暴雨圖集》24h暴雨參數(shù)等值線圖，得研究區(qū)域中心年最大24h暴雨量均值為98.8mm，變差系數(shù)CV為0.45，采用CS=3.5CV。查P-Ⅲ型曲線KP表，即可算出研究區(qū)域不同重現(xiàn)期24h設(shè)計點雨量，因研究區(qū)域面積僅50.70萬m2，點面折算系數(shù)為1，故用設(shè)計點雨量代表面雨量計算，不同重現(xiàn)期24h設(shè)計面雨量成果見表9。

表9 不同標(biāo)準(zhǔn)24h設(shè)計面雨量成果

重現(xiàn)期3年5年10年20年Kp10813116018824h設(shè)計點雨量(mm)1066112900157941859224h設(shè)計面雨量(mm)10661129001579418592

表10 不同重現(xiàn)期下設(shè)計徑流深與排澇模數(shù)對比

2.6.3 排澇模數(shù)計算

根據(jù)研究區(qū)附近天津耳閘站24h設(shè)計暴雨的時程分配，將不同重現(xiàn)期設(shè)計面雨量進(jìn)行分配，代入超滲蓄滿垂向產(chǎn)流模型中，計算不同重現(xiàn)期下的設(shè)計徑流深，代公式(28)求得相應(yīng)的排澇模數(shù)，與2003年《除澇手冊》中給出的設(shè)計排澇模數(shù)進(jìn)行對比，結(jié)果見表10。

由表10可知，不同重現(xiàn)期下的設(shè)計徑流深和排澇模數(shù)與原手冊設(shè)計排澇模數(shù)相比發(fā)生了較大變化，與原設(shè)計排澇模數(shù)相比均有所增加，增大的幅度在16.44%～3.77%之間，且變化幅度隨重現(xiàn)期增大而減小。由降雨趨勢分析可知，基準(zhǔn)期與評價期的降雨沒有發(fā)生變化，因此設(shè)計排澇模數(shù)變化，主要是由于設(shè)施農(nóng)業(yè)占比增加造成的。

3 結(jié)論

(1)通過開展農(nóng)田小區(qū)現(xiàn)場試驗，布設(shè)水文要素測量站點，減小了模型輸入數(shù)據(jù)的誤差，同時考慮了地下水位與土壤含水量的變化，提高了模型的模擬精度。與單一超滲模型和蓄滿模型進(jìn)行對比，有效性分別提高了0.869和0.912，驗證了在該農(nóng)田小區(qū)建立的超滲—蓄滿垂向耦合模型的有效性。

(2)在SCE-UA全局參數(shù)優(yōu)化中結(jié)合Sobol′全局敏感性分析的成果，綜合考慮模型結(jié)構(gòu)，充分利用實測資料，高效優(yōu)化了模型參數(shù)，各項評價指標(biāo)的Nsah系數(shù)均高于0.6，相對誤差在20%以內(nèi)。

(3)隨著農(nóng)業(yè)現(xiàn)代化進(jìn)程的加快，試驗農(nóng)田小區(qū)中設(shè)施農(nóng)業(yè)所占比例逐年提高，對水文要素和產(chǎn)流過程有較大影響， 3年、5年、10年、20年一遇重現(xiàn)期下的設(shè)計排澇模數(shù)與原《除澇手冊》相比，分別增加16.44%、8.63%、7.50%、3.77%。通過降雨趨勢分析，基準(zhǔn)期與評價期的降雨沒有顯著線性變化趨勢，因此設(shè)計排澇模數(shù)的變化主要是由于設(shè)施農(nóng)業(yè)占比增加引起。由此可以看出設(shè)施農(nóng)業(yè)對設(shè)計排澇模數(shù)的影響較大，這為天津市近郊區(qū)設(shè)施農(nóng)業(yè)下農(nóng)田設(shè)計排澇模數(shù)的修訂提供了有效依據(jù)，對合理確定排澇系統(tǒng)規(guī)模、有效控制澇災(zāi)、保障農(nóng)業(yè)安全生產(chǎn)具有重要的現(xiàn)實意義。

[1] Mwangi H M， Julich S，Patil S D，et al.Modelling the impact of agroforestry on hydrology of Mara River Basin in East Africa.Hydrological Processes， 2016, 30(18): 3139～3155

[2] 胡華浪， 李偉方，易湘生，等.土地利用/覆蓋變化對永定河流域水文過程的可能影響模擬.中國農(nóng)業(yè)資源與區(qū)劃， 2016, 37(3): 74～83

[3] 包為民， 王從良.垂向混合產(chǎn)流模型及應(yīng)用.水文, 1997，(3)： 18～21

[4] 王維， 馮忠倫，楊偉，等.基于SCE-UA算法的新安江模型與垂向混合產(chǎn)流模型參數(shù)優(yōu)化應(yīng)用研究.中國農(nóng)村水利水電， 2017，(3)： 26～30

[5] Song X， Zhang J，Zhan C，et al.Global sensitivity analysis in hydrological modeling：Review of concepts，methods，theoretical framework，and applications.Journal of Hydrology， 2015, 523(523): 739～757

[6] 宋曉猛， 張建云，占車生，等.水文模型參數(shù)敏感性分析方法評述.水利水電科技進(jìn)展， 2015, 35(6): 105～112

[7] Sobol I M. Sensitivity analysis for non-linear mathematical models.Mathematical Modeling & Computational Experiment， 1993， 1： 407～414

[8] 齊偉， 張弛，初京剛，等.Sobol′方法分析TOPMODEL水文模型參數(shù)敏感性.水文， 2014, 34(2): 49～54

[9] Zhang C， Chu J，F(xiàn)u G.Sobol′s sensitivity analysis for a distributed hydrological model of Yichun River Basin，China.Journal of Hydrology， 2013, 480：58～68

[10]張小麗， 彭勇，徐煒，等.基于 Sobol′ 方法的新安江模型參數(shù)敏感性分析.南水北調(diào)與水利科技， 2014，(2)： 20～24

[11]Duan Q Y， Gupta V K，Sorooshian S.Shuffled complex evolution approach for effective and efficient global minimization.Plenum Press. 1993

[12]黃學(xué)群， 李瑾，宋建輝，等.天津現(xiàn)代農(nóng)業(yè)園區(qū)發(fā)展模式與對策研究.中國農(nóng)業(yè)資源與區(qū)劃， 2012, 33(6): 79～84

[13]包為民. 格林—安普特下滲曲線的改進(jìn)和應(yīng)用.人民黃河， 1993,(9)： 1～3

[14]趙人俊. 流域水文模擬.北京：中國水利電力出版社，1984

[15]方磊. 基于LHS抽樣的不確定性分析方法在概率安全評價中的應(yīng)用研究.合肥：中國科學(xué)技術(shù)大學(xué)， 2015

[16]Tang Y， Reed P，Wagener T，et al.Comparing sensitivity analysis methods to advance lumped watershed model identification and evaluation.Hydrology & Earth System Sciences， 2007, 3(6): 793～817

[17]Sobol I M. On sensitivity estimation for nonlinear mathematical models.Keldysh Applied Mathematics Institute， 1990, 2(1): 112～118

[18]Sobol I M. Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates.Mathematics & Computers in Simulation， 2001， 55(1-3)： 271～280

[19]中華人民共和國水利部. 水文情報預(yù)報規(guī)范SD 138-85.北京：中國標(biāo)準(zhǔn)出版社,2008

[20]于暢， 郝曼秋，高成，等.平原城市自排區(qū)排澇模數(shù)計算方法研究.水資源與水工程學(xué)報， 2014,(6)： 184～186

[21]水利部. 灌溉排水工程設(shè)計規(guī)范GB50288-99.北京：中國計劃出版社,1999

[22]Nossent J， Elsen P，Bauwens W.Sobol′ sensitivity analysis of a complex environmental model.Elsevier Science Publishers B.V.2011, 26(12): 1515～1525

[23]李發(fā)文， 馮平，張超.天津北三河地區(qū)垂向耦合產(chǎn)流模型及應(yīng)用.吉林大學(xué)學(xué)報(地)， 2011, 41(2): 459～464

[24]馬海波， 董增川，張文明，等.SCE-UA算法在TOPMODEL參數(shù)優(yōu)化中的應(yīng)用.河海大學(xué)學(xué)報(自然科學(xué)版)， 2006, 34(4): 361～365

[25]傅寧， 史珺，任雨.1951～2012年天津地區(qū)汛期降水變化特征分析.干旱區(qū)資源與環(huán)境， 2014, 28(4): 83～88