范嚴(yán)偉，趙文舉，畢貴權(quán)(蘭州理工大學(xué)能源與動力工程學(xué)院，蘭州 730050)

非飽和土壤導(dǎo)水率、擴散率以及土壤水分特征曲線等是數(shù)值模擬非飽和土壤水分運動必不可少的重要參數(shù)，準(zhǔn)確獲取這些參數(shù)是定量研究土壤水分和溶質(zhì)運動的基礎(chǔ)[1]。多年來，國內(nèi)外學(xué)者建立了許多經(jīng)驗?zāi)Ｐ蛠砻枋鐾寥浪痔卣髑€和非飽和導(dǎo)水率曲線，其中，Van Genuchten模型應(yīng)用最為廣泛[2]?；赗ichards方程開發(fā)的HYDRUS-1D/2D/3D軟件[3]，均采用VG模型擬合。為此，許多土壤物理學(xué)家致力于確定VG模型中的參數(shù)。邵明安等[4]提出推求非飽和土壤水分運動參數(shù)的簡單入滲法，用以確定VG模型中的參數(shù)α和n；李春友[5]、馬英杰[6]、陳大春[7]、廖林仙[8]、郭向紅[9]等分別采用單純形調(diào)優(yōu)法、阻尼最小二乘法、隨機粒子群算法、智能算法、混合遺傳算法擬合土壤水力參數(shù)；魏義長等[10]運用Matlab軟件的非線性擬合和非線性回歸函數(shù)獲得VG模型參數(shù)值；王薇等[11]利用RETC軟件依據(jù)土壤的顆粒組成和容重轉(zhuǎn)換出VG模型參數(shù)值；查元源等[12]根據(jù)瞬時剖面法，推求VG模型參數(shù)Ks、α和n。

模型參數(shù)的準(zhǔn)確程度是數(shù)值計算精度高低的關(guān)鍵，而模型中各參數(shù)對模擬結(jié)果的影響程度因其物理意義的不同而有所差異。因此，有必要了解VG模型參數(shù)變化對土壤入滲特性的影響規(guī)律，以有效識別關(guān)鍵參數(shù)，減少率定參數(shù)的數(shù)量，提高模型的運行效率。李法虎等[13]對土壤水分運動參數(shù)的靈敏性能進行了理論分析；王志濤等[14]分析了VG模型參數(shù)變化對粉壤土入滲特性的影響。目前，尚未開展VG模型參數(shù)變化對不同質(zhì)地土壤入滲特性的影響。為此，基于HYDRUS-1D模型，對VG模型中的θr、θs、α、n和Ks參數(shù)進行單因素擾動，模擬獲得土壤濕潤鋒運移距離和累積入滲量，分析VG模型參數(shù)變化對入滲特性的影響程度，以便為優(yōu)化和改進模型結(jié)構(gòu)提供依據(jù)。

1 數(shù)學(xué)模型建立

1.1 基本方程

HYDRUS-1D模型對土壤水分運動的模擬采用Richards方程[15]。其表達式為：

(1)

式中：z為垂向坐標(biāo)，規(guī)定z向上為正；θ為體積含水量，cm3/cm3；h在飽和區(qū)和非飽和區(qū)分別為壓力水頭和基質(zhì)勢，cm；t為入滲時間，min；K(h)為非飽和導(dǎo)水率，cm/min。

Richards方程中涉及的土壤水分特征曲線θ(h)和非飽和土壤導(dǎo)水率K(h)采用VG模型擬合。即

(2)

K(h)=KsSle[1-(1-S1/me)m]2

(3)

式中：Se=(θ-θr)/(θs-θr)；θs為殘余含水量，cm3/cm3；θs為飽和含水量，cm3/cm3；α、n和m為土壤物理特性有關(guān)的擬合參數(shù)，α單位為cm-1；n>1；m=1-1/n；l=0.5；Ks為飽和導(dǎo)水率，cm/min。

1.2 定解條件

在HYDRUS-1D模型中，土壤質(zhì)地、容重是通過VG模型中的參數(shù)體現(xiàn)的，壓力水頭和土壤初始含水量對累積入滲量和濕潤鋒運移距離影響較小。初始條件和邊界條件不會影響VG模型參數(shù)，故在模擬時以定水頭作為HYDRUS-1D模型的上邊界條件，上邊界水頭值取為2 cm，入滲過程中，濕潤鋒未到達下界面。土壤初始含水量根據(jù)土壤質(zhì)地取為一恒值，其中粉壤土為0.12 cm3/cm3，壤土為0.10 cm3/cm3，砂壤土為0.08 cm3/cm3，砂土為0.06 cm3/cm3。其定解條件為：

θ=θ0,0≤z≤L,t=0

(4)

h=h0,z=0,t>0

(5)

θ=θ0,z=L,t>0

(6)

式中：θ0為土壤的初始含水量，cm3/cm3；h0為初始壓力水頭，cm。

1.3 模型參數(shù)

各類土壤及其Van Genuchten模型參數(shù)均采用HYDRUS-1D模型中典型土壤的設(shè)定值，結(jié)果見表1。模擬過程中，參數(shù)擾動幅度設(shè)置為5個水平，分別為-20%、-10%、0%、10%和20%。

表1 不同土質(zhì)VG模型水力特性參數(shù)Tab.1 Hydraulic parameters in VG model of different soils

2 VG模型參數(shù)變化對入滲特性的影響分析

2.1 分析方法

為定量分析VG模型中各參數(shù)變化對土壤入滲特性的影響，采用靈敏度系數(shù)來表示參數(shù)變化對入滲過程的影響程度，即：

(7)

式中：RSi為靈敏度系數(shù)；ΔPi為某參數(shù)Pi的變化量；ΔQi為某參數(shù)Pi的變化而引起目標(biāo)函數(shù)Qi的變化量。

2.2 參數(shù)θr變化對入滲特性的影響

定量分析參數(shù)θr對不同質(zhì)地土壤濕潤鋒運移距離Z的影響，將不同擾動變幅下的θr值及T=100 min時的濕潤鋒運移距離值代入式(7)，得參數(shù)θr的靈敏度系數(shù)RSθr。如表2所示。

表2 參數(shù)θr對不同土質(zhì)濕潤鋒運移距離的影響Tab.2 Influence of parameter θr to transport distance of wetting front for different soils

由表2可知：θr的擾動對濕潤鋒運移距離無影響。

同理可得出參數(shù)θr對不同擾動變幅下累積入滲量I的靈敏度。如表3所示。

表3 參數(shù)θr對不同土質(zhì)累積入滲量的影響Tab.3 Influence of parameter θr to the cumulative infiltration for different soils

由表3可知：θr的擾動對土壤累積入滲量影響非常微弱。

綜合分析說明：土壤殘余含水量θr的變化對土壤入滲特性影響很小。考慮到土壤殘余含水量θr的獲取比較困難，因此，在沒有確定的或?qū)嶒灉y得的θr值的情況下，將穩(wěn)定枯竭點上的基質(zhì)勢(φ=-15 000 cm)所對應(yīng)的含水量作為θr的估算值是可行的[16]。

2.3 參數(shù)θs變化對入滲特性的影響

模擬參數(shù)θs擾動下的濕潤鋒運移距離、累積入滲量與時間t的關(guān)系曲線圖，如圖1所示。

圖1 參數(shù)θs對不同土質(zhì)入滲特性的影響Fig.1 Effect of parameter θs on infiltration characteristics for different soils

定量分析參數(shù)θs對不同質(zhì)地土壤濕潤鋒運移距離的影響，將不同擾動變幅下的θs值及T=100 min時的濕潤鋒運移距離值代入式(7)，得參數(shù)θs的靈敏度系數(shù)RSθs。如表4所示。

表4 參數(shù)θs對不同土質(zhì)濕潤鋒運移距離的影響Tab.4 Influence of parameter θs on transport distance of wetting front for different soils

由圖1及表4可知：θs的擾動與垂直濕潤鋒運移距離成負(fù)相關(guān)，對其影響較大，且隨土質(zhì)由細(xì)變粗而影響逐漸增大；θs的負(fù)擾動對濕潤鋒運移距離的影響程度明顯強于正擾動。

同理可得出參數(shù)θs對不同擾動變幅下累積入滲量I的靈敏度。如表5所示。

由圖1及表5可知：θs的擾動與土壤累積入滲量成正相關(guān)，對其影響較大，且隨土質(zhì)由細(xì)變粗而影響逐漸減小。

綜合分析說明：土壤飽和含水量θs的變化對土壤入滲特性影響較大。但θs可通過實驗準(zhǔn)確測定。因此，在求解VG模型參數(shù)時，可將θs取為實測值，避免VG模型參數(shù)確定過程中θs的誤差。

表5 參數(shù)θs對不同土質(zhì)累積入滲量的影響Tab.5 Influence of parameter θs to the cumulative infiltration for different soils

2.4 參數(shù)α變化對入滲特性的影響

模擬參數(shù)α擾動下的濕潤鋒運移距離、累積入滲量與時間t的關(guān)系曲線圖，如圖2所示。

圖2 參數(shù)α對不同土質(zhì)入滲特性的影響Fig.2 Effect of parameters α on infiltration characteristics for different soils

定量分析參數(shù)α對不同質(zhì)地土壤濕潤鋒運移距離的影響，將不同擾動變幅下的α值及T=100 min時的濕潤鋒運移距離值代入式(7)，得參數(shù)α的靈敏度系數(shù)RSα。如表6所示。

表6 參數(shù)α對不同土質(zhì)濕潤鋒運移距離的影響Tab.6 Influence of parameter α on transport distance of wetting front for different soils

由圖2及表6可知：參數(shù)α的擾動與濕潤鋒運移距離成負(fù)相關(guān)，對其影響較小，且隨土質(zhì)由細(xì)變粗而影響逐漸減小。

同理可得出參數(shù)α對不同擾動變幅下累積入滲量I的靈敏度。如表7所示。

表7 參數(shù)α對不同土質(zhì)累積入滲量I的影響Tab.7 Influence of parameter α to the cumulative infiltration for different soils

由圖2及表7可知：參數(shù)α的擾動與土壤累積入滲量成負(fù)相關(guān)，對其影響較小，且隨土質(zhì)由細(xì)變粗而影響逐漸減小。

綜合分析說明：參數(shù)α的變化對土壤入滲特性影響較小，且隨土質(zhì)由細(xì)變粗而影響逐漸減小。在實際應(yīng)用時，VG模型參數(shù)α一般是通過多次迭代擬合確定，但擬合過程中往往造成參數(shù)的不唯一性。因此，在擬合求解VG模型參數(shù)時，可給出參數(shù)α的取值范圍，以解決參數(shù)α不唯一問題，或采用簡單入滲法來直接推求，推求過程中應(yīng)保證吸滲率的準(zhǔn)確性[4]。

2.5 參數(shù)n變化對入滲特性的影響

模擬參數(shù)n擾動下的濕潤鋒運移距離、累積入滲量與時間t的關(guān)系曲線圖，如圖3所示。

圖3 參數(shù)n對不同土質(zhì)入滲特性的影響Fig.3 Effect of parameter n on infiltration characteristics for different soils

定量分析參數(shù)n對不同質(zhì)地土壤濕潤鋒運移距離的影響，將不同擾動變幅下的n值及T=100 min時的濕潤鋒運移距離值代入式(7)，得參數(shù)n的靈敏度系數(shù)RSn。如表8所示。

由圖3及表8可知：參數(shù)n的擾動與濕潤鋒運移距離成正相關(guān)，對粉壤土、壤土影響較大，隨土質(zhì)由細(xì)變粗而影響迅速減小，表現(xiàn)為對砂壤土、砂土影響較小。

表8 參數(shù)n對不同土質(zhì)濕潤鋒運移距離的影響Tab.8 Influence of parameter n on transport distance of wetting front for different soils

同理可得出參數(shù)n對不同擾動變幅下累積入滲量I的靈敏度。如表9所示。

表9 參數(shù)n對不同土質(zhì)累積入滲量的影響Tab.9 Influence of parameter s to the cumulative infiltration for different soils

由圖3及表9可知：參數(shù) 的擾動與累積入滲量成正相關(guān)，對粉壤土、壤土影響較大，隨土質(zhì)由細(xì)變粗而影響迅速減小，表現(xiàn)為對砂壤土、砂土影響較小。

綜合分析說明：參數(shù)n的變化對細(xì)質(zhì)土壤入滲特性影響較大，對粗質(zhì)土壤入滲特性影響較小。參數(shù)n與α相同，在實際應(yīng)用時，一般是通過多次迭代擬合來確定，但擬合過程中容易造成參數(shù)n的不唯一性。因此，在擬合求解VG模型參數(shù)時，需給出參數(shù)n的取值范圍，以解決參數(shù)n不唯一問題?；虿捎煤唵稳霛B法來直接推求，推求過程中應(yīng)保證吸滲率的準(zhǔn)確性[4]。

2.6 參數(shù)Ks變化對入滲特性的影響

模擬參數(shù)Ks擾動下的濕潤鋒運移距離、累積入滲量與時間t的關(guān)系曲線圖，如圖4所示。

圖4 參數(shù)Ks對不同土質(zhì)入滲特性的影響Fig.4 Effect of parameter Ks on infiltration characteristics for different soils

定量分析參數(shù)Ks對不同質(zhì)地土壤濕潤鋒運移距離的影響，將不同擾動變幅下的Ks值及T=100 min時的濕潤鋒運移距離值代入式(7)，得參數(shù)Ks的靈敏度系數(shù)RSKs。如表10所示。

表10 參數(shù)Ks對不同土質(zhì)濕潤鋒運移距離的影響Tab.10 Influence of parameter Ks on transport distance of wetting front for different soils

由圖4及表10可知：參數(shù)Ks的擾動與濕潤鋒運移距離成正相關(guān)，對其影響較大，且隨土壤質(zhì)地由細(xì)變粗而影響逐漸變大。

同理可得出參數(shù)Ks對不同擾動變幅下累積入滲量I的靈敏度。如表11所示。

表11 參數(shù)Ks對不同土質(zhì)累積入滲量的影響Tab.11 Influence of parameter Ks to the cumulative infiltration for different soils

由圖4及表11可知：參數(shù)Ks的擾動與累積入滲量成正相關(guān)，對其影響較大，且隨土壤質(zhì)地由細(xì)變粗而影響逐漸變大。

綜合分析說明：飽和導(dǎo)水率Ks的變化對土壤入滲特性影響較大，且隨土壤質(zhì)地由細(xì)變粗而影響逐漸變大。Ks可通過實驗測定，但空間變異較大，應(yīng)考慮其尺度效應(yīng)。或在求解VG模型參數(shù)時，同時測定土壤水分特征曲線θ(h)和非飽和土壤水?dāng)U散率D(θ)，采用同步擬合方法，得出VG模型參數(shù)θr、θs、α、n和Ks值，以消除單個參數(shù)引起的總體誤差。

3 結(jié) 語

基于HYDRUS-1D模型，采用單因素影響法，分析VG模型參數(shù)變化對土壤入滲特性的影響，結(jié)果表明：

(1)參數(shù)θr對濕潤鋒運移距離和累積入滲量影響很小，基本無影響。

(2)參數(shù)θs對濕潤鋒運移距離影響較大，成負(fù)相關(guān)關(guān)系，隨土質(zhì)由細(xì)變粗而影響增大，負(fù)擾動強于正擾動；同時，θs與累積入滲量成正相關(guān)關(guān)系，對粉壤土、壤土影響較大，對砂壤土、砂土影響較小。

(3)參數(shù)α對濕潤鋒運移距離和累積入滲量影響都較小，成負(fù)相關(guān)關(guān)系，隨土質(zhì)由細(xì)變粗而影響減小。

(4)參數(shù)n與濕潤鋒運移距離和累積入滲量成正相關(guān)，對粉壤土、壤土影響較大，對砂壤土、砂土影響較小。

(5)參數(shù)Ks對濕潤鋒運移距離和累積入滲量影響較大，成正相關(guān)關(guān)系，隨土質(zhì)由細(xì)變粗而影響增大。

□

[1] 邵明安,王全九,黃明斌. 土壤物理學(xué)[M]. 北京:高等教育出版社, 2006.

[2] Van Genuchten M Th. A closed form equation for predicting the hydraulic conductivity of unsaturated soils[J]. Soil Sci.Soc.Am.J.,1980,44:892-898.

[3] Simunek J, Van Genuchten M Th, Sejna M. HYDRUS: Model Use, Calibration, and Validation[J]. Transactions of the ASABE,2012,55(4):1 261-1 274.

[4] 邵明安,王全九,Robert Horton. 推求土壤水分運動參數(shù)的簡單入滲法I.理論分析[J]. 土壤學(xué)報,2000,37(1):1-6.

[5] 李春友,任 理,李保國. 利用優(yōu)化方法求算Van Genuehten方程參數(shù)[J]. 水科學(xué)進展,2001,12(4):473-478.

[6] 馬英杰,虎膽·吐馬爾拜,沈冰. 利用阻尼最小二乘法求解Van Genuchten方程參數(shù)[J]. 農(nóng)業(yè)工程學(xué)報, 2005,21(8):179-181.

[7] 陳大春,馬英杰. 基于隨機粒子群算法的Van Genuchten方程參數(shù)優(yōu)化求解[J]. 農(nóng)業(yè)工程學(xué)報,2006,22(12):82-85.

[8] 廖林仙,邵孝侯,徐俊增. 基于智能算法推求Van Genuchten方程的參數(shù)[J]. 水利學(xué)報,2007,38(S1):696-700.

[9] 郭向紅,孫西歡,馬娟娟. 基于混合遺傳算法估計Van Genuchten方程參數(shù)[J]. 水科學(xué)進展,2009,20(5):677-682.

[10] 魏義長,劉作新,康玲玲,等. 土壤持水曲線Van Genuchten模型求參的Matlab實現(xiàn)[J]. 土壤學(xué)報,2004,41(3):380-386.

[11] 王 薇,孟 杰,虎膽·吐馬爾白. RETC推求土壤水動力學(xué)參數(shù)的室內(nèi)試驗研究[J]. 河北農(nóng)業(yè)大學(xué)學(xué)報2008,31(1):99-102,106.

[12] 查元源,周發(fā)超,楊金忠. 一種由土壤剖面含水率估算土壤水力參數(shù)的方法[J]. 水利學(xué)報,2011,,42(8):883-891.

[13] 李法虎,傅建平,孫雪峰. 土壤水分運動參數(shù)的確定及其靈敏性能分析[J].灌溉排水學(xué)報,1993,12(3):6-14.

[14] 王志濤,繳錫云,韓紅亮,等. 土壤垂直一維入滲對VG模型參數(shù)的敏感性分析[J]. 河海大學(xué)學(xué)報(自然科學(xué)版),2013,41(1):80-84.

[15] Simunek J, Van Genuchten M T, Sejna M. The HYDRUS-1D software package for simulating the one-dimensional movement of water, heat, and multiple solutes in variably-saturated media[J]. University of California-Riverside Research Reports, 2005,3:1-240.

[16] 王金生,楊志峰, 陳家軍,等. 包氣帶土壤水分滯留特征研究[J]. 水利學(xué)報, 2000,31(2):1-6.