胡 丹，伍靖?jìng)?，黃介生

(武漢大學(xué)水資源與水電工程科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，武漢 430072)

0 引 言

土壤水和地下水運(yùn)動(dòng)是水文循環(huán)的重要環(huán)節(jié)，能夠準(zhǔn)確模擬土壤水、地下水運(yùn)動(dòng)對(duì)水資源高效利用具有重要意義[1,2]。地下水與土壤水之間有著直接的水力聯(lián)系和水通量交換，為了更好地模擬區(qū)域尺度土壤水、地下水的運(yùn)動(dòng)，近年來(lái)發(fā)展了一些飽和-非飽和耦合模型，這些模型的主要區(qū)別在于對(duì)土壤水運(yùn)動(dòng)的模擬，其主要方式有3種：第一種是水均衡模型，如MODFLOW中的REC和ET模塊[3]，以土壤水質(zhì)量平衡為原理，獨(dú)立考慮蒸發(fā)、根系吸水等水文過(guò)程，由于水均衡模型不涉及水頭，無(wú)法與地下水的計(jì)算直接聯(lián)系，過(guò)分簡(jiǎn)化了土壤水運(yùn)動(dòng)，故會(huì)帶來(lái)較大誤差，但求解相對(duì)簡(jiǎn)單；第二種是基于動(dòng)力波方程，如MODFLOW-UZF1[4]，模擬土壤水運(yùn)動(dòng)時(shí)只考慮重力作用，而忽略毛管力作用；第三種是基于理查茲方程，如三維非飽和-三維飽和的耦合模型MODFLOW-VSF[5]、一維非飽和-三維飽和的耦合模型MODFLOW-HYDRUS[6]，理查茲方程是由達(dá)西定律和質(zhì)量守恒定律推導(dǎo)而得，對(duì)土壤水運(yùn)動(dòng)的描述更加準(zhǔn)確，但由于具有高度非線性，計(jì)算成本也更大。飽和-非飽和耦合模型能否做到模擬精度和計(jì)算成本的平衡直接關(guān)系到模型的運(yùn)用推廣。Twarakavi等[7]對(duì)MODFLOW-UZF1、MODFLOW-VSF和MODFLOW-HYDRUS三種耦合模型進(jìn)行了比較，發(fā)現(xiàn)進(jìn)行區(qū)域模擬時(shí)，MODFLOW-HYDRUS和MODFLOW-VSF的模擬精度比MODFLOW-UZF1更高，在計(jì)算要求方面，MODFLOW-HYDRUS和MODFLOW-UZF1的計(jì)算量比MODFLOW-VSF更小，尤其對(duì)于淺地下水、蒸發(fā)強(qiáng)烈的地區(qū)，MODFLOW-HYDRUS進(jìn)行土壤水-地下水區(qū)域模擬時(shí)，相比MODFLOW-UZF1和MODFLOW-VSF，在模擬精度和計(jì)算量方面能做到更好的平衡。

飽和-非飽和耦合模型進(jìn)行區(qū)域水分運(yùn)動(dòng)模擬時(shí)，往往需要輸入較多的參數(shù)，且水文地質(zhì)參數(shù)和土壤水力參數(shù)空間變異性較大，導(dǎo)致模擬結(jié)果存在較大不確定性。模型參數(shù)的獲得目前主要有實(shí)測(cè)法和參數(shù)率定法，其中實(shí)驗(yàn)法測(cè)得的數(shù)據(jù)由于尺度效應(yīng)的存在，很難直接運(yùn)用到模型中。參數(shù)率定法包括試錯(cuò)法、參數(shù)自動(dòng)優(yōu)化方法以及試錯(cuò)法與參數(shù)自動(dòng)優(yōu)化法聯(lián)合優(yōu)化。傳統(tǒng)的試錯(cuò)法依賴使用者對(duì)模型的了解及經(jīng)驗(yàn)，主觀性太強(qiáng)，花費(fèi)時(shí)間較長(zhǎng)，且難以達(dá)到參數(shù)最優(yōu)化。參數(shù)自動(dòng)優(yōu)化方法能夠較快得到優(yōu)化解，但是對(duì)初值的依賴較大，不同的初值可能導(dǎo)致差別較大的率定結(jié)果[8]。試錯(cuò)法與參數(shù)自動(dòng)優(yōu)化法聯(lián)合優(yōu)化既能發(fā)揮模型使用者的知識(shí)和經(jīng)驗(yàn)，又能利用先進(jìn)的優(yōu)化技術(shù)[9]。

PEST[10]是獨(dú)立于計(jì)算模型的參數(shù)優(yōu)化方法，已經(jīng)分別在地下水、土壤水運(yùn)動(dòng)模型中得到廣泛運(yùn)用[11-13]，但是在耦合模型中運(yùn)用效果和效率如何尚未有研究。

本文以內(nèi)蒙古河套灌區(qū)永聯(lián)試驗(yàn)區(qū)地下水觀測(cè)井?dāng)?shù)據(jù)和土壤含水率實(shí)測(cè)數(shù)據(jù)為基礎(chǔ)，采用試錯(cuò)法和PEST相結(jié)合的方法，對(duì)MODFLOW-HYDRUS模型的水文地質(zhì)參數(shù)和土壤水力參數(shù)進(jìn)行優(yōu)化，以提升區(qū)域尺度上飽和-非飽和帶中水分運(yùn)動(dòng)的模擬效果。

1 材料與方法

1.1 研究區(qū)域及數(shù)據(jù)

研究區(qū)域?yàn)閮?nèi)蒙古河套灌區(qū)五原縣永聯(lián)村境內(nèi)的永聯(lián)試驗(yàn)區(qū)(見(jiàn)圖1)，東經(jīng)108°00′～108°01′，北緯41°04′～41°05′，南北長(zhǎng)約13 km，東西寬約5 km，總面積約2 875 hm2，其中灌溉面積約1 453 hm2，土地利用類型主要是耕地、裸地(包括荒地、渠道、排水溝等)和村莊。試驗(yàn)區(qū)地勢(shì)自南向北稍有傾斜，中間較高，東西兩側(cè)略低，整體較為平坦，地面高程1 028.6～1 025.6 m。試驗(yàn)區(qū)邊界條件較為清楚，北邊是六排干排水邊界，南邊是皂火渠，東西兩邊分別是永什分干溝和乃永分干溝(見(jiàn)圖1)。

試驗(yàn)區(qū)內(nèi)布置有11口地下水觀測(cè)井(見(jiàn)圖1)，每5 d觀測(cè)1次地下水位，每個(gè)地下水觀測(cè)井附近布置有剖面土壤含水率監(jiān)測(cè)點(diǎn)，每10 d進(jìn)行1次取樣測(cè)定，取樣深度為0～100 cm，分為0～10，10～30，30～50，50～70，70～100 cm 5個(gè)層次。地下水位觀測(cè)采用測(cè)鐘法，土壤含水率采用烘干法測(cè)定。

圖1 研究區(qū)示意圖Fig.1 Sketch map of study area

1.2 飽和-非飽和耦合模型

非飽和帶水流運(yùn)動(dòng)考慮為一維垂向運(yùn)動(dòng)，用一維Richards方程來(lái)描述：

(1)

式中：θ為土壤的體積含水量，cm3/cm3；t為時(shí)間，d；h為壓力水頭，m；z為空間坐標(biāo)，向上為正，m；S為匯項(xiàng)，m3/(m3·d)；K(h)為非飽和水力傳導(dǎo)度，m/d。

飽和帶水流運(yùn)動(dòng)為三維運(yùn)動(dòng)，其控制偏微分方程為：

(2)

式中：h為水頭，m；x，y，z為空間坐標(biāo)，m；Kxx，Kyy，Kzz為滲透系數(shù)在x，y，z方向上的分量，這里假定滲透系數(shù)主軸方向與坐標(biāo)軸方向一致，m/d；W為匯項(xiàng)，m3/(m3·d)；t為時(shí)間，d；SS為給水度，無(wú)量綱。

非飽和帶與飽和帶的耦合過(guò)程參考Seo等[5]的介紹。

本文采用2007年5月1日至11月30日(共214 d)期間實(shí)測(cè)的地下水位和土壤含水率數(shù)據(jù)對(duì)模型參數(shù)進(jìn)行率定。

(1)網(wǎng)格劃分與土壤剖面設(shè)置。將研究區(qū)域劃分為25×50×2網(wǎng)格，其中有效網(wǎng)格為742個(gè)，水平方向上網(wǎng)格大小為160 m×264 m，垂直方向上，第一、二層底部高程分別為1 020.29 m和973.75 m。根據(jù)地面高程、地下水位、土地利用類型等情況，采用k均值聚類分析方法，在研究區(qū)內(nèi)設(shè)置25個(gè)土壤剖面(見(jiàn)圖2)，其中1～7號(hào)土壤剖面為耕地，8～14號(hào)為排水溝，15～16號(hào)為人民支渠，17號(hào)為皂火渠，18～21號(hào)為荒地，22～25號(hào)為村莊。土壤剖面的高度為地面到1 020.29 m高程處。

圖2 網(wǎng)格與土壤剖面Fig.2 Grids and soil profiles

(2)初始輸入數(shù)據(jù)。以2007年5月1日實(shí)測(cè)的地下水位數(shù)據(jù)插值得到整個(gè)區(qū)域的地下水位作為初始地下水位[見(jiàn)圖3(b)]；土壤剖面的初始?jí)毫λ^分布按照Seo等[5]的方法計(jì)算得到。

圖3 地面高程和初始地下水位Fig. 3 Ground surface elevation and initial groundwater table elevation

(3)邊界條件。上邊界條件為變化的大氣邊界和灌溉補(bǔ)給。從五原縣氣象站獲得該研究區(qū)域氣象資料，包括降雨量、氣溫、濕度、風(fēng)速、日照時(shí)數(shù)等，采用FAO-56[14]介紹的方法計(jì)算每日的潛在蒸發(fā)量和蒸騰量，將人民支渠口部的凈引水量作為區(qū)域內(nèi)的灌溉水量，而灌溉滲漏量概化為人民支渠上的線性補(bǔ)給。耕地(1～7號(hào)土壤剖面)的上邊界條件包括降雨、蒸發(fā)、蒸騰和灌溉[見(jiàn)圖4(a)]，荒地(18～21號(hào)土壤剖面)、排水溝(8～14號(hào)土壤剖面)和皂火渠(17號(hào)土壤剖面)的上邊界條件包括降雨和蒸發(fā)[見(jiàn)圖4(b)]，人民支渠(15～16號(hào)土壤剖面)的上邊界條件包括降雨、蒸發(fā)和灌溉滲漏量[見(jiàn)圖4(c)]，村莊(22～25號(hào)土壤剖面)的上邊界條件只有降雨[見(jiàn)圖4(d)]。第二含水層下是約1 m的黏土層，故將下邊界條件視為零通量邊界。將六排干、永什分干和乃永分干三條排水溝設(shè)置為排水邊界，皂火渠灌溉渠道設(shè)置為河流邊界。

圖4 不同土壤剖面的上邊界條件Fig.4 Upper boundary conditions of different soil profiles

(4)模型參數(shù)。研究區(qū)域兩個(gè)含水層分別是沙壤土和砂土，水文地質(zhì)參數(shù)的初始值根據(jù)五原縣永聯(lián)村鉆孔水文地質(zhì)報(bào)告取值(見(jiàn)表1)，van Genuchten[15]模型中土壤水力參數(shù)的初始值根據(jù)Carsel and Parrish[16]取值(見(jiàn)表1)。然后采用試錯(cuò)法對(duì)模型各參數(shù)進(jìn)行調(diào)整(見(jiàn)表1)，最后進(jìn)行參數(shù)自動(dòng)優(yōu)化。

表1 土壤水力參數(shù)和水文地質(zhì)參數(shù)取值Tab.1 Values of soil hydraulic parameters and hydrogeological parameters

1.3 PEST參數(shù)優(yōu)化算法

PEST算法采用高斯牛頓法與梯度下降法結(jié)合快速收斂性，又具備梯度下降法的全局搜索性，它能夠使目標(biāo)函數(shù)快速有效地收斂，進(jìn)行參數(shù)優(yōu)化時(shí)模型運(yùn)行的次數(shù)比其他算法更少[17]。

PEST算法的核心是對(duì)目標(biāo)函數(shù)求解最小值，目標(biāo)函數(shù)是實(shí)驗(yàn)觀測(cè)值與計(jì)算值的帶權(quán)重差異函數(shù)，根據(jù)參數(shù)優(yōu)化的需求和實(shí)驗(yàn)觀測(cè)值的種類，可將觀測(cè)值分組并賦予相應(yīng)的權(quán)重，以達(dá)到模擬值與實(shí)測(cè)值的最佳匹配。本研究包含地下水位、土壤含水率2組目標(biāo)函數(shù)，其公式如下：

(3)

式中：Φ為目標(biāo)函數(shù)；Oi和Si分別表示第i個(gè)實(shí)驗(yàn)觀測(cè)值和第i個(gè)模擬值；p表示實(shí)驗(yàn)觀測(cè)值的個(gè)數(shù)；ω表示權(quán)重系數(shù)；下標(biāo)wt和sm分別對(duì)應(yīng)地下水位和土壤含水率。

1.4 評(píng)價(jià)方法

采用決定系數(shù)(R2)、納什系數(shù)(NSE)、均方根誤差(RMSE)、標(biāo)準(zhǔn)化均方根誤差(NRMSE)4種指標(biāo)來(lái)評(píng)價(jià)模型的模擬效果：

(7)

決定系數(shù)(R2)用來(lái)衡量模擬值和實(shí)測(cè)值的相關(guān)密切程度，R2越接近1表示相關(guān)性越好；納什系數(shù)(NSE)可用來(lái)評(píng)價(jià)模型的模擬質(zhì)量，NSE越接近1說(shuō)明模擬質(zhì)量越好；均方根誤差(RMSE)用來(lái)衡量模擬值與實(shí)測(cè)值之間的偏差，RMSE越接近0，表示模擬值與觀測(cè)值的差異越小，模擬效果越好；標(biāo)準(zhǔn)化均方根誤差(NRMSE)也是用來(lái)衡量模擬誤差的大小，NRMSE值越小，說(shuō)明模擬效果越好。

2 結(jié)果與分析

參數(shù)自動(dòng)優(yōu)化共選取19個(gè)參數(shù)作為優(yōu)化對(duì)象，包括15個(gè)土壤水力參數(shù)和4個(gè)水文地質(zhì)參數(shù)，參數(shù)優(yōu)化初值為試錯(cuò)法調(diào)整后的結(jié)果，各個(gè)參數(shù)的取值范圍和優(yōu)化結(jié)果見(jiàn)表1。

將試錯(cuò)法和PEST參數(shù)優(yōu)化后的地下水位模擬值和觀測(cè)井實(shí)測(cè)值進(jìn)行比較(見(jiàn)圖5)，同時(shí)將模擬開(kāi)始后第46、86、128、168天土壤剖面含水率的模擬值和實(shí)測(cè)值進(jìn)行比較(見(jiàn)圖6)。將所有的地下水位觀測(cè)值與對(duì)應(yīng)時(shí)刻的模擬值繪制在一張圖中(見(jiàn)圖7)，將所有的土壤含水率實(shí)測(cè)值與對(duì)應(yīng)時(shí)刻的模擬值繪制在一張圖中(見(jiàn)圖7)，并計(jì)算評(píng)價(jià)模型效果的4個(gè)指標(biāo)(見(jiàn)表2)。

采用2008年5月1日至11月31日期間的地下水位和土壤含水率實(shí)測(cè)值對(duì)模型進(jìn)行驗(yàn)證，將地下水位觀測(cè)值與對(duì)應(yīng)時(shí)刻的模擬值繪制在一張圖中(見(jiàn)圖8)，將土壤含水率實(shí)測(cè)值與對(duì)應(yīng)時(shí)刻的模擬值繪制在一張圖中(見(jiàn)圖8)，并計(jì)算評(píng)價(jià)模型效果的4個(gè)指標(biāo)(見(jiàn)表2)。

圖5 地下水位模擬值與實(shí)測(cè)值對(duì)比Fig.5 Computed vs. observed groundwater table

圖6 不同時(shí)刻土壤剖面含水率模擬值與實(shí)測(cè)值對(duì)比Fig.6 Computed vs. observed soil moisture of different time

模擬期觀測(cè)變量試錯(cuò)法R2NSERMSENRMSEPESTR2NSERMSENRMSE率定期地下水位0．7650．7310．5210．1040．8080．7860．4640．093土壤含水率0．4160．4080．0640．1370．5290．5280．0570．122驗(yàn)證期地下水位0．8040．6030．5800．1290．8240．6410．5520．123土壤含水率0．3250．0320．0660．1740．4380．3950．0530．140

圖7 率定期地下水位與土壤含水率散點(diǎn)圖Fig. 7 Scatter plots of groundwater table and soil moisture for calibration period

圖8 驗(yàn)證期地下水位與土壤含水率散點(diǎn)圖Fig. 8 Scatter plots of groundwater table and soil moisture for validation period

從圖5可以看出，地下水位模擬值與觀測(cè)值的變化趨勢(shì)基本一致；從圖6可以看出，模擬開(kāi)始后不同時(shí)刻土壤剖面含水率模擬值與實(shí)測(cè)值較為接近；從圖7、圖8可以發(fā)現(xiàn)，地下水位和土壤含水率散點(diǎn)都較好地集中在1:1直線附近。

表2的模型效果評(píng)價(jià)指標(biāo)顯示，PEST算法進(jìn)行參數(shù)優(yōu)化，提高了模型的模擬效果。

在率定期，地下水位的決定系數(shù)由0.765提高到0.808，提高了6%；納什系數(shù)由0.731提高到0.786，提高了8%；均方根誤差由0.521減小到0.464，減小了11%；標(biāo)準(zhǔn)化均方根誤差由0.104減小到0.093，減小了11%。

在率定期，土壤含水率的決定系數(shù)由0.416提高到0.529，提高了27%；納什系數(shù)由0.408提高到0.528，提高了29%；均方根誤差由0.064減小到0.057，減小了11%；標(biāo)準(zhǔn)化均方根誤差由0.137減小到0.122，減小了11%。

在驗(yàn)證期，地下水位的決定系數(shù)由0.804提高到0.824，提高了2%；納什系數(shù)由0.603提高到0.641，提高了6%；均方根誤差由0.580減小到0.552，減小了5%；標(biāo)準(zhǔn)化均方根誤差由0.129減小到0.123，減小了5%。

在驗(yàn)證期，土壤含水率的決定系數(shù)由0.325提高到0.438，提高了35%；納什系數(shù)由0.032提高到0.395，提高了11.3倍；均方根誤差由0.066減小到0.053，減小了20%；標(biāo)準(zhǔn)化均方根誤差由0.174減小到0.140，減小了20%。

同時(shí)從表2可以看出，耦合模型對(duì)地下水位的模擬效果要優(yōu)于土壤含水率。從試錯(cuò)法的結(jié)果看，地下水位在率定期和驗(yàn)證期決定系數(shù)分別為0.765和0.804，納什系數(shù)分別為0.731和0.603；而土壤含水率在率定期和驗(yàn)證期決定系數(shù)分別為0.416和0.325，納什系數(shù)分別為0.408和0.032。通過(guò)PEST參數(shù)優(yōu)化后，地下水位在率定期和驗(yàn)證期決定系數(shù)分別為0.808和0.824，納什系數(shù)分別為0.786和0.641；而土壤含水率在率定期和驗(yàn)證期決定系數(shù)分別為0.529和0.438，納什系數(shù)分別為0.528和0.395。土壤含水率的各項(xiàng)評(píng)價(jià)指標(biāo)都低于地下水位的相應(yīng)指標(biāo)。對(duì)土壤水運(yùn)動(dòng)模擬效果較差，可能的原因是本次模擬將整個(gè)土壤剖面視為均質(zhì)土壤，與實(shí)際土壤剖面為非均質(zhì)土壤的情況不符。

3 結(jié) 語(yǔ)

(1)本文采用試錯(cuò)法和PEST優(yōu)化算法相結(jié)合的方法對(duì)MODFLOW-HYDRUS耦合模型的水文地質(zhì)參數(shù)和土壤水力參數(shù)進(jìn)行優(yōu)化，提高了模擬精度。與試錯(cuò)法的結(jié)果相比，PEST參數(shù)優(yōu)化后決定系數(shù)提高了2%～35%，納什系數(shù)提高了6～11.3倍，均方根誤差減小了5%～20%，標(biāo)準(zhǔn)化均方根誤差減小了5%～20%，其中驗(yàn)證期土壤含水率的各項(xiàng)指標(biāo)提高或減小的幅度均為最大，決定系數(shù)、納什系數(shù)、均方根誤差和標(biāo)準(zhǔn)化均方根誤差變化幅度分別為35%、11.3倍、20%和20%。

(2)MODFLOW-HYDRUS耦合模型模擬區(qū)域尺度上地下水、土壤水運(yùn)動(dòng)時(shí)，在地下水運(yùn)動(dòng)模擬上的表現(xiàn)要優(yōu)于土壤水運(yùn)動(dòng)的模擬。從試錯(cuò)法與PEST參數(shù)優(yōu)化后結(jié)果來(lái)看，地下水位的決定系數(shù)在0.765～0.824之間，納什系數(shù)在0.603～0.786之間；而土壤含水率的決定系數(shù)在0.325～0.529之間，納什系數(shù)在0.032～0.528之間，明顯低于地下水位相對(duì)應(yīng)的指標(biāo)。導(dǎo)致土壤水運(yùn)動(dòng)模擬效果不佳的原因可能是將整個(gè)土壤剖面視為均質(zhì)土壤，與實(shí)際土壤剖面非均質(zhì)土壤的情況不符，在今后運(yùn)用模型時(shí)可以考慮按照實(shí)際的非均質(zhì)土壤剖面來(lái)設(shè)置。

□

[1] 楊欣坤, 王 宇, 趙蘭坡, 等. 土壤水動(dòng)力學(xué)參數(shù)及其影響因素研究進(jìn)展[J]. 中國(guó)農(nóng)學(xué)通報(bào), 2014,30(3):38-43.

[2] 薛禹群. 中國(guó)地下水?dāng)?shù)值模擬的現(xiàn)狀與展望[J]. 高校地質(zhì)學(xué)報(bào), 2010,16(1):1-6.

[3] Harbaugh A W, E R Banta, M C Hill, et al. MODFLOW-2000, the U.S. Geological Survey modular ground-water model user guide to modularization concepts and the ground-water fl ow process.USGS, Denver, CO,2000.

[4] Niswonger R G, D E Prudic, R S Regan. Documentation of the Unsaturated-Zone Flow (UZF1) package for modeling unsaturated flow between the land surface and the water table with MODFLOW-2005. Techniques and Methods 6-A19. Available at http:∥pubs.usgs.gov/tm/2006/tm6a19/ (verified 14 Apr. 2008). USGS, Denver, CO,2006.

[5] Thoms R B, R L Johnson, R W Healy. User’s guide to the Variably Saturated Flow (VSF) process for MODFLOW. Techniques and Methods 6-A18. Available at http://pubs.usgs.gov/tm/2006/ tm6a18/ (verified 15 Apr. 2008). USGS, Reston, VA, 2006.

[8] 舒曉娟, 陳洋波, 黃鋒華, 等. PEST在WetSpa分布式水文模型參數(shù)率定中的應(yīng)用[J]. 水文, 2009,29(5):45-49.

[9] 章四龍, 劉九夫. 通用模型參數(shù)率定技術(shù)研究[J]. 水文, 2005,25(1):9-12.

[10] Doherty J. PEST Model-Independent Parameter Estimation User Manual[M]. 5th Edition. Brisbane, Australia: Watermark Numerical Computing, 2004.

[11] Zyvoloski G, E Kwicklis, A A. Eddebbarh, et al. The site-scale saturated zone flow model for Yucca Mountain: calibration of different conceptual models and their impact on flow paths[J]. Journal of Contaminant Hydrology, 2003,62-63(2):731-750.

[12] 董艷輝, 李國(guó)敏, 郭永海, 等. 應(yīng)用并行PEST 算法優(yōu)化地下水模型參數(shù)[J]. 工程地質(zhì)學(xué)報(bào), 2010,18(1):140-145.

[13] Fang Q X, Green T R, Ma L, et al. Optimizing soil hydraulic parameters in RZWQM2 under fallow conditions[J]. Soil Science Society of America Journal, 2010,74(6):1 897-1 913.

[14] Allen R G, Pereira L S, Raes D, et al. Crop Evapotranspiration: Guidelines for Computing Crop Water Requirements[D]. Rome, Italy: United Nations Food and Agriculture Organization, 1998.

[15] van Genuchten, M.Th. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils[J]. Soil Science Society of America Journal, 1980,44(44):892-898.

[16] Carsel R F, Parrish R S. Developing joint probability distributions of soil water retention characteristics[J]. Water Resources Research, 1988,24(5):755-769.

[17] Bahremand A, De Smedt F. Distributed hydrological modeling and sensitivity analysis in Torysa Watershed, Slovakia[J]. Water Resources Management, 2008,22(3):393-408.