王子龍，常廣義，姜秋香，付 強(qiáng)，陳偉杰，林百健

（東北農(nóng)業(yè)大學(xué)水利與土木工程學(xué)院，哈爾濱 150030）

土壤水分特征曲線是模擬土壤水分和溶質(zhì)運(yùn)移重要輸入?yún)?shù)，廣泛應(yīng)用于地下水補(bǔ)給、農(nóng)業(yè)與土壤化學(xué)、水土保持等領(lǐng)域，可反映土壤孔隙分布狀況[1-3]。

測定土壤水分特征曲線費(fèi)時(shí)長、成本高、耗費(fèi)人力，常用土壤水分特征曲線模型有Brooks-Corey（BC）[4]，Gardner[5]、van Genuchten[6]（VG）及Kosugi[7]、Fredlund 和 Xing[8]、Dual-porosity模型[9]等。土壤傳遞函數(shù)模型分為點(diǎn)預(yù)測模型和參數(shù)預(yù)測模型[10]。點(diǎn)預(yù)測模型利用田間含水量33 kPa、凋萎含水量1 500 kPa和土壤水力性質(zhì)構(gòu)建函數(shù)關(guān)系；參數(shù)預(yù)測模型基于模型參數(shù)與物理性質(zhì)（如土壤質(zhì)地、有機(jī)質(zhì)、顆粒組成等）構(gòu)建函數(shù)關(guān)系。土壤傳遞函數(shù)模型研究主要有線性和非線性擬合、人工神經(jīng)網(wǎng)絡(luò)[11]等方法。土壤傳遞函數(shù)模型成果有Vereecken[12]模型，利用土壤中顆粒組成、有機(jī)碳等與VG模型參數(shù)建立關(guān)系式，Cosby[13]和Saxton[14]模型利用顆粒組成與實(shí)測土壤水分特征曲線參數(shù)作回歸分析，姚姣轉(zhuǎn)等在VG模型基礎(chǔ)上建立科爾沁沙地土壤傳遞函數(shù)等[15]。

比較分析土壤水分特征曲線傳遞函數(shù)模型和經(jīng)驗(yàn)?zāi)Ｐ皖A(yù)測能力，如Woston等比較21種土壤傳遞函數(shù)，結(jié)果表明點(diǎn)預(yù)測模型模擬精度優(yōu)于參數(shù)預(yù)測模型[16]；Kern分析對比6種土壤傳遞函數(shù)模型，發(fā)現(xiàn)Rawls模型預(yù)測精度最高[17]；Tietje等評價(jià)13種土壤傳遞函數(shù)預(yù)測效果，結(jié)果表明Vereecken模型對土壤水分特征曲線預(yù)測值與實(shí)測值偏差最小，模擬精度最高[18]；Matlan等對比4種經(jīng)驗(yàn)?zāi)Ｐ?，結(jié)果表明MG模型對不同質(zhì)地土壤土壤水分特征曲線模擬效果最好[19]；丁新原等利用VG和Gardner模型模擬塔里木沙漠公路防護(hù)林SWCC，發(fā)現(xiàn)VG模型模擬精度高于Gardner模型[20]。目前，土壤水分特征曲線模型對比研究較多，但模型普適性和精確度依賴土壤類型，對比分析特定土壤適應(yīng)性尤為必要。

本文以松嫩平原黑土區(qū)黑土為例，利用VG、BC、LND、MG 4種土壤水分特征曲線模型模擬不同質(zhì)地黑土土壤水分特征曲線，對比分析4種經(jīng)驗(yàn)?zāi)Ｐ湍M精度，選取最優(yōu)模型。研究結(jié)果為松嫩平原黑土區(qū)土壤水分特征曲線預(yù)測提供可靠模型和技術(shù)支持。

1 材料與方法

1.1 研究區(qū)概況

研究區(qū)處于松嫩平原黑土區(qū)南部黑龍江省哈爾濱市（見圖1），位于東經(jīng)125°42'～130°10'，北緯44°04'～46°40'，全年平均降水量569.1 mm，冬長夏短，屬于溫帶大陸性季風(fēng)氣候。哈爾濱市土壤類型較多，以黑土為主，包括壤土、黏土、砂壤土和砂黏壤土，少量砂土分布江河兩岸。黑土土壤養(yǎng)分含量較豐富，是主要耕作土壤，土壤有機(jī)質(zhì)含量高，適宜農(nóng)作物生長。

圖1 研究區(qū)地理位置及采樣點(diǎn)分布Fig.1 Geographical location and sampling point distribution in the study area

1.2 供試土樣采集與測定

選取壤土、砂黏壤土、砂壤土、黏土、砂土5種土壤質(zhì)地作為供試土樣，采樣點(diǎn)分布如圖1所示。

每個(gè)土樣采集0～30 cm土層，分別用100 cm3環(huán)刀取原狀土，自封袋取散土。風(fēng)干散土土樣、過篩2 mm，采用MS2000激光粒度儀分析土壤顆粒，按國際制分類，分為砂粒（粒徑0.02～2 mm）、粉粒（粒徑0.002～0.02 mm）和黏粒（粒徑＜0.002 mm）。采用高速離心機(jī)法測定土壤水分特征曲線，分別測定土樣1、3、5、10、30、50、100、300、400、500、700、900、1 000、1 200 kPa吸力下體積含水率。干容重采用環(huán)刀法測定，土壤樣品質(zhì)地分類如表1所示。

表1 土壤樣品質(zhì)地分類Table1 Textureclassification of soil samples

1.3 土壤水分特征曲線經(jīng)驗(yàn)?zāi)Ｐ秃喗?/h3>

土壤水分特征曲線經(jīng)驗(yàn)?zāi)Ｐ透鶕?jù)參數(shù)數(shù)量分為三參數(shù)模型、四參數(shù)模型和五參數(shù)模型。研究表明，模型參數(shù)越多，實(shí)測數(shù)據(jù)擬合精度更高。

選取4種土壤水分特征曲線經(jīng)驗(yàn)?zāi)Ｐ蛌an Genuchten（VG）、Brooks and Corey（BC）、Modified Gardner（MG）、Log-Normal Distribution（LND）模型模擬黑土區(qū)土壤水分特征曲線。其中VG模型是van Genuchten在1981年提出應(yīng)用最廣泛模型；BC模型是由Brooks和Corey較早提出的土壤水分特征曲線模型，該模型簡單，應(yīng)用較廣；MG模型是由Matlan提出的指數(shù)型經(jīng)驗(yàn)?zāi)Ｐ停Ｐ蛥?shù)有明確物理意義，應(yīng)用范圍較廣；LND模型是Kosugi提出的對數(shù)正態(tài)分布模型。4種模型均為四參數(shù)模型，具有簡單性和普遍使用性等特點(diǎn)。

①VG模型表達(dá)式如下：

式中，θ—體積含水率；θs—土壤飽和體積含水率；θr—滯留土壤體積含水率；h—土壤水吸力（cm）；α、n、m為表征模型形狀參數(shù)，其中α數(shù)值上等于進(jìn)氣值倒數(shù)。

②BC模型表達(dá)式如下：

式中，Se—飽和度。

λ—土壤孔隙尺寸分布參數(shù)

③MG模型表達(dá)式如下：

式中，b、t為模型參數(shù)。

④LND模型表達(dá)式如下：

式中，q為模型參數(shù)。

4個(gè)模型除θs和θr參數(shù)外，還包含兩個(gè)模型參數(shù)，分別用P1和P2表示，其中P1單位為kPa、P2為分形參數(shù)，如表2所示。

表2 土壤水分特征曲線模型參數(shù)Table 2 Parameters of soil water characteristic curve models

1.4 模型參數(shù)優(yōu)化方法

1.41 粒子群優(yōu)化算法

采用粒子群優(yōu)化算法（Particle Swarm Optimization，PSO）優(yōu)化模型參數(shù)，模仿鳥群捕食行為提出粒子群優(yōu)化算法。首先初始化1個(gè)種群N，把每個(gè)種群中個(gè)體假設(shè)d維搜索空間中1個(gè)粒子，第i個(gè)粒子位置和速度分別為Xi=(xi1,xi2,…,xid)和Vi=(vi1,vi2,…,vid)，i=1,2,3…N，每次迭代中，粒子通過自我更新找出最優(yōu)解，即個(gè)體最優(yōu)解pbest，Pi=(pi1,pi2…,pid)；gbest是整個(gè)種群N種中最優(yōu)解，兩個(gè)最優(yōu)解確定后，粒子速度更新和位置更新方程如下：

其中，w為慣性權(quán)因子，c1和c2為學(xué)習(xí)因子，rand1和rand2為0～1均勻分布隨機(jī)數(shù)。

1.4.2 參數(shù)模型構(gòu)建

為精確獲得4種模型模型參數(shù)，其中飽和含水率θs采用實(shí)測數(shù)據(jù)，以模型擬合體積含水率值與土壤實(shí)測體積含水率殘差平方和最小為目標(biāo)函數(shù)，求解該非線性最優(yōu)化問題，即表達(dá)式如下：

式中，θobs為土壤體積含水率實(shí)測值，θmod為模型擬合值，N為樣本數(shù)。

1.5 模型精度評價(jià)指標(biāo)

采用均方根誤差（RMSE）和Pearson相關(guān)系數(shù)R作為定量評價(jià)指標(biāo)。

式中，X為體積含水率實(shí)測值，Xi為經(jīng)驗(yàn)?zāi)Ｐ蛿M合值。RSME評估模型整體誤差；R值反映實(shí)測值與模擬值趨勢。

2 結(jié)果分析

2.1 不同土壤質(zhì)地土壤水分特征曲線變化趨勢

土壤水分特征曲線反映土壤持水能力。通過離心機(jī)測定各吸力段下對應(yīng)土壤體積含水率，繪制不同質(zhì)地吸力段下體積含水率與土壤水基質(zhì)吸力變化趨勢如圖2所示。

圖2 不同質(zhì)地土壤水分特征曲線Fig.2 Soil moisture characteristic curve of different soil textures

由圖2可知，吸力值小于100 kPa時(shí)的低吸力段，土壤水分特征曲線急劇下降，變化率大，因?yàn)樵撾A段土壤大孔隙排水占主導(dǎo)；當(dāng)土壤吸力大于100 kPa時(shí)，隨著吸力增大，曲線變化率逐漸趨于0，吸力值大于100 kPa時(shí)，中高吸力段土壤大孔隙水分在重力作用下排空，僅小孔隙存留水分，毛管力對水分產(chǎn)生作用力，土壤具有較好持水能力，土壤水分特征曲線緩慢降低并趨于平穩(wěn)。在相同吸力條件下，體積含水率黏土＞壤土＞砂黏壤土＞砂壤土＞砂土，說明黏土持水能力最強(qiáng)，壤土次之，砂土最弱。結(jié)合表1不同質(zhì)地顆粒組成百分含量可知，黏粒含量越高，持水能力越強(qiáng)，由于黏粒粒徑較小，黏粒含量越高，小孔隙越多，比表面積越大，對水作用力越強(qiáng)；而砂土砂粒含量較高，粒徑較大，比表面積較小，水吸附能力弱。因此降水后砂土大孔隙較多易飽和，形成徑流，不利于水土保持。

2.2 土壤水分特征曲線擬合

本文采用VG、BC、MG、LND 4種經(jīng)驗(yàn)?zāi)Ｐ蛿M合5種土壤質(zhì)地實(shí)測土壤含水率數(shù)據(jù)曲線，基于粒子群優(yōu)化算法對4種模型參數(shù)優(yōu)化值如表3所示，繪制5種土壤類型土壤水分特征曲線實(shí)測值與采用4種模型擬合值關(guān)系如圖3所示。由圖3可知，VG模型對5種土壤擬合值與實(shí)測值接近，其中實(shí)測值在VG模型曲線上。滯留土壤含水率為0.034～0.092，符合實(shí)際土壤滯留含水率。5種土壤質(zhì)地中黏土飽和含水率和滯留含水率相對較高；砂土飽和含水率和滯留含水率最低，由土壤性質(zhì)決定。參數(shù)n值反映土壤水分特征曲線坡度，土壤水分特征曲線坡度越緩，n值越大，反之越小[23]。砂土坡度最陡，n值不是最小，故認(rèn)為n為表征曲線形狀系數(shù)[24]。VG模型參數(shù)P1值為0.0375～0.323，則土壤進(jìn)氣值為3.096～26.667，與飽和含水率對應(yīng)進(jìn)氣值相符。BC模型對砂土和砂黏壤土擬合趨勢較好，壤土擬合趨勢較差，可較好擬合5種土壤質(zhì)地曲線。BC模型優(yōu)化參數(shù)中滯留含水率為1.20～5.68×10-4，可知BC模型低估5種土壤質(zhì)地滯留含水率；模型參數(shù)P1值為0.2051～0.9001，土壤進(jìn)氣值為1.1109～4.8756。可知BC模型低估飽和含水率對應(yīng)進(jìn)氣值。MG模型對5種類型土壤擬合效果較好，實(shí)測值在MG模型曲線上。其中除砂土外，對其他4種土壤擬合值更接近實(shí)測值，砂土實(shí)測值擬合偏差較大。MG模型擬合參數(shù)滯留含水率為0.0124～0.1469，與土壤中滯留含水率相近。LND模型優(yōu)化參數(shù)中滯留含水率為0.00012～0.0736。可知LND模型低估壤土和砂黏壤土滯留含水率。LND模型對5種類型土壤擬合值與實(shí)測值偏差較小，表明LND模型可較好擬合5種土壤質(zhì)地土壤水分特征曲線。

表3 不同土壤質(zhì)地模型參數(shù)Table 3 Model parameters of different soil textures

2.3 土壤水分特征曲線模型比較及優(yōu)選

比較Brooks-Corey（BC模型）、Van Genuchten（VG模型）、Modified Garder（MG模型）和Log-Normal Distribution（LND模型）適用性，對比分析模擬結(jié)果。

由表4可知，VG模型對5種類型土壤擬合值與實(shí)測值相關(guān)系數(shù)R為0.996～0.998；BC模型R為0.9547～0.9869；MG模型R為0.994～0.998；LND模型R為0.9820～0.9960，可見4種模型擬合值與實(shí)測值變化趨勢相似。其中，除壤土外，MG模型對其他4種土壤擬合效果優(yōu)于VG、BC、LND模型；另外，VG模型對壤土模擬效果最佳。由圖3可知，VG、MG、LND 3種模型在低吸力段（0～100 kPa）壤土、砂黏壤土、砂壤土、黏土、砂土擬合值與實(shí)測值重合，而BC模型擬合值明顯偏離實(shí)測值；在中高吸力段（＞100 kPa）VG、MG、BC、LND 4種模型擬合值均在實(shí)測值上方。其中BC模型擬合值與實(shí)測值偏差較大，而VG、MG、LND 3種模型擬合值略高于實(shí)測值，可知VG、MG、LND3種模型對5種類型土壤擬合效果優(yōu)于BC模型。

由殘差平方和值可知，4種模型均可較好擬合不同土壤質(zhì)地含水量與水吸力間關(guān)系，VG、MG模型殘差平方和值顯著明顯低于BC和LND模型，其中BC模型殘差平方和值最高，可知MG和VG模型對5種類型土壤擬合效果最佳，BC模型擬合效果最差。從RMSE值來看，MG模型對5種土壤質(zhì)地RMSE值為0.0033～0.0107；VG模型RMSE變化為0.0039～0.0092；LND 模型 RMSE 變化為 0.0096～0.0114；BC模型RMSE變化為0.0140～0.0250。結(jié)果表明，MG和VG模型RSME值明顯小于LND和BC模型，其中BC模型RMSE值比其他3種模型RMSE值高一個(gè)數(shù)量級，BC模型對5種土壤類型模擬效果最差。綜合比較擬合誤差分析可知，VG和MG模型對黑土區(qū)土壤模擬精度最高。

表4 4種模型模擬土壤水分特征曲線R、RSME和殘差平方和Table4 R,RSMEand sum of squared residualsof soil water characteristic curvecompared to simulationsfor four models

3 討 論

本研究利用4種土壤水分特征曲線經(jīng)驗(yàn)?zāi)Ｐ蛿M合5種不同質(zhì)地土壤曲線，4種模型擬合結(jié)果較佳，其中BC模型在擬合過程中低估土壤滯留含水率和土壤進(jìn)氣值，由于實(shí)測飽和含水率存在誤差，導(dǎo)致出現(xiàn)低估土壤滯留含水率和土壤進(jìn)氣值現(xiàn)象。模型模擬結(jié)果中高吸力段擬合大于實(shí)測值，由于高速離心機(jī)測定過程中土壤結(jié)構(gòu)被破壞，容重改變，導(dǎo)致誤差[25]。

對比分析4種經(jīng)驗(yàn)?zāi)Ｐ湍M結(jié)果可知，VG和MG模型模擬5種不同質(zhì)地土壤效果最好。丁新原等比較VG和其他模型的土壤水分特征曲線模擬結(jié)果，VG模型模擬精度最高，與本研究結(jié)果一致[20,23,26]，因VG和MG模型中各參數(shù)均有明確物理意義[19,27]。Matlan等利用4種模型模擬分析黏土、砂土、粉土3種土壤質(zhì)地，表明BC模型對實(shí)測值擬合效果最差，與本文研究結(jié)果一致[19]。其他4種土壤質(zhì)地，BC模型對砂土模擬效果最優(yōu)，表明BC模型更適合模擬粗質(zhì)地土壤水分特征曲線[28]。

4 結(jié) 論

黑土區(qū)不同質(zhì)地土壤土壤體積含水率因土壤吸力變化持水能力存在差異，在相同吸力下，土壤黏粒含量越高，土壤體積含水率越高。在低吸力段，大孔隙排水，砂粒含量多，曲線變化率越快，釋水能力越強(qiáng)；在中高吸力段，土壤小孔隙比例較高，土壤對水吸附作用較強(qiáng)，因此土壤水分特征曲線逐漸趨于平穩(wěn)。

本文利用粒子群優(yōu)化算法優(yōu)化VG、BC、MG、LND 4種土壤水分特征曲線模型參數(shù)。比較4種經(jīng)驗(yàn)?zāi)Ｐ蛿M合值與實(shí)測值，VG、MG模型對黑土區(qū)土壤水分特征曲線模擬精度較高。預(yù)測黑土區(qū)土壤水分特征曲線推薦使用VG和MG模型。