設(shè)施農(nóng)業(yè)對(duì)農(nóng)田排澇模數(shù)影響*

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

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

0 引言

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

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

模型參數(shù)敏感性分析和參數(shù)優(yōu)化對(duì)模型模擬結(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)行對(duì)比，證明了SCE-UA算法在兩種模型中的有效性。

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

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

1 材料與方法

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

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

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

經(jīng)調(diào)研，該農(nóng)田小區(qū)在2008年之前全部為經(jīng)濟(jì)林木，隨著天津“菜籃子”工程的實(shí)施，設(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)流時(shí)，采用具有流域分布特征的格林—安普特下滲曲線計(jì)算超滲產(chǎn)流量和下滲量[13]，然后以超滲產(chǎn)流的下滲量作為蓄滿產(chǎn)流中的降雨輸入量。在土壤缺水量大的面積上，補(bǔ)充土壤含水量不產(chǎn)流，在缺水量小的面積上，補(bǔ)足土壤缺水量后，發(fā)生蓄滿產(chǎn)流[3]。

超滲產(chǎn)流量計(jì)算過程為[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為土壤缺水量對(duì)下滲率影響的靈敏度系數(shù)；FMM為平均下滲率為FM時(shí)上層最大點(diǎn)的下滲率，mm/h；FA為實(shí)際下滲量，mm；WM為土壤平均蓄水容量，mm；W為土壤含水量，mm；B1為下滲能力分配曲線的指數(shù)，反映了下滲能力分布的不均勻性；R1為超滲產(chǎn)流量，mm。

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

(6)

WMM=WM×(1+B2)

(7)

(8)

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

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

R=R1+R2

(9)

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

(10)

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

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

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

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

DRS=P×IMP×α0

(11)

α0=0.95，RS>0

(12)

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

1.2.3 總產(chǎn)流計(jì)算

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

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

(13)

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

綜上所述，在研究區(qū)域透水面上采用超滲—蓄滿垂向耦合結(jié)構(gòu)，地面徑流R1按超滲產(chǎn)流計(jì)算，地面以下徑流按蓄滿產(chǎn)流計(jì)算，并依據(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ù)的一階敏感度與總敏感度計(jì)算，采用SCE-UA算法進(jìn)行全局參數(shù)自動(dòng)率定，同時(shí)結(jié)合敏感性分析的成果與模型適用性評(píng)價(jià)標(biāo)準(zhǔn)進(jìn)行人工率定，最終得到模型的參數(shù)。

1.3.1 Latin超立方抽樣

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

1.3.2 Sobol′全局敏感度分析法

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

(14)

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

(15)

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

(16)

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

(17)

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

(18)

(19)

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

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

(20)

(21)

(22)

(23)

1.3.3 SCE-UA算法

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

圖5 SCE-UA算法流程圖

1.4 模型適用性評(píng)價(jià)

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

(24)

(25)

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

(26)

(27)

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

1.5 排澇模數(shù)計(jì)算

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

(28)

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

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

序號(hào)參數(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 實(shí)例應(yīng)用

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

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

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

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

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

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

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

通過分析圖6，發(fā)現(xiàn)當(dāng)以徑流深和徑流深相對(duì)誤差為評(píng)價(jià)指標(biāo)時(shí)，參數(shù)在相同等級(jí)的降雨中所呈現(xiàn)的規(guī)律是一致的。在特大暴雨時(shí)，一階敏感度最大的參數(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的自動(dòng)率定

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

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

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

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

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

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

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

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

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

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

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

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

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

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

2.5 模型適用性評(píng)價(jià)

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

表7 模型評(píng)價(jià)指標(biāo)分析

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

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

降雨場(chǎng)次降雨量(mm)實(shí)測(cè)流量(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可知，率定期與驗(yàn)證期地下水變化量與土壤含水量變化量的Nash效率系數(shù)均大于0.6，各項(xiàng)相對(duì)誤差均控制在20%之內(nèi)，可以認(rèn)為模型是合理的。

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

2.6 排澇模數(shù)計(jì)算

2.6.1 降雨趨勢(shì)分析

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

圖7 降雨趨勢(shì)分析

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

2.6.2 設(shè)計(jì)暴雨計(jì)算

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

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

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

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

2.6.3 排澇模數(shù)計(jì)算

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

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

3 結(jié)論

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

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

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

[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] 胡華浪， 李偉方，易湘生，等.土地利用/覆蓋變化對(duì)永定河流域水文過程的可能影響模擬.中國(guó)農(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)用研究.中國(guó)農(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ù)敏感性分析方法評(píng)述.水利水電科技進(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ā)展模式與對(duì)策研究.中國(guó)農(nóng)業(yè)資源與區(qū)劃， 2012, 33(6): 79～84

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

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

[15]方磊. 基于LHS抽樣的不確定性分析方法在概率安全評(píng)價(jià)中的應(yīng)用研究.合肥：中國(guó)科學(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]中華人民共和國(guó)水利部. 水文情報(bào)預(yù)報(bào)規(guī)范SD 138-85.北京：中國(guó)標(biāo)準(zhǔn)出版社,2008

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

[21]水利部. 灌溉排水工程設(shè)計(jì)規(guī)范GB50288-99.北京：中國(guó)計(jì)劃出版社,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é)報(bào)(地)， 2011, 41(2): 459～464

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

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