基于加速遺傳算法的免耕農(nóng)田土壤水分特征曲線分析

袁宏偉，湯廣民，袁先江，趙 暉

(1.安徽省水利部淮委水利科學(xué)研究院 水利水資源安徽省重點(diǎn)實(shí)驗(yàn)室，合肥 230088；2.山東省臨沂市國(guó)土資源局土地整理中心，山東 臨沂 276000)

0 引 言

土壤物理性質(zhì)的測(cè)定是測(cè)量和確定土壤不發(fā)生化學(xué)變化就表現(xiàn)出來(lái)的性質(zhì)，土壤物理性質(zhì)與其水分、空氣、熱量狀況以及對(duì)農(nóng)田灌排的要求和耕作效果密切相關(guān)。因此，測(cè)定土壤物理性質(zhì)是調(diào)節(jié)水、肥、氣、熱矛盾，進(jìn)行土壤管理的基礎(chǔ)工作，具有非常重要的意義[1-3]。

安徽省淮北平原砂姜黑土區(qū)是全省主要的糧食產(chǎn)區(qū)，種植的主要糧食作物為冬小麥和夏玉米。隨著近幾年政府倡導(dǎo)秸稈禁燒政策之后，夏玉米的耕作方式發(fā)生了極大的變化，由焚燒冬小麥秸稈后翻耕播種變?yōu)榻斩掃€田免耕播種。耕作方式的變化也造成了夏玉米播種前的土壤物理性質(zhì)與以往不同[1]，因此有必要對(duì)夏玉米播種后的土壤基礎(chǔ)物理性質(zhì)進(jìn)行重新取樣分析。

土壤水分特征曲線作為土壤的最基本物理特性，是研究包氣帶土壤水分保持能力和運(yùn)移的基礎(chǔ)及與水分運(yùn)動(dòng)相關(guān)的熱量、溶質(zhì)運(yùn)動(dòng)的基礎(chǔ)[6]。土壤水分特征曲線表示土壤水吸力和含水率之間的關(guān)系，可用于分析土壤中孔隙大小及分布、土壤的持水性及土壤水的有效性，關(guān)系到農(nóng)作物灌溉制度的確定[4-10]。土壤水分特征曲線測(cè)定方法主要有經(jīng)驗(yàn)公式法及直接測(cè)定法，其中經(jīng)驗(yàn)公式法主要有Brooks-Corey模型、Van Genuchten模型(簡(jiǎn)稱VG模型)和Arya-Paris模型等，直接測(cè)定法則包括離心機(jī)法、張力計(jì)法、壓力膜儀法、平衡水汽壓法、砂芯漏斗法等[4-10]。因VG模型的線型與實(shí)測(cè)曲線的吻合程度較高且參數(shù)具有明確物理意義，故目前國(guó)內(nèi)外多采用該模型對(duì)土壤水分特征曲線進(jìn)行描述[4,5]。如何得到精確的模型參數(shù)成為相關(guān)研究的重點(diǎn)。VG模型參數(shù)的確定是典型的非線性擬合問題，傳統(tǒng)最小二乘法等擬合方法效果不理想，而隨著人工智能技術(shù)的發(fā)展，可以利用智能算法對(duì)模型參數(shù)進(jìn)行自動(dòng)優(yōu)化，粒子群算法、模擬退火算法、遺傳算法和螢火蟲算法等智能算法已經(jīng)被應(yīng)用于解決該問題[4,5]。上述算法各有優(yōu)缺點(diǎn)，但普遍存在計(jì)算量大，優(yōu)化速度慢等問題，因此本文利用金菊良[13]等提出的加速遺傳算法(accelerating genetic algorithm)，來(lái)優(yōu)化VG模型參數(shù)，為參數(shù)優(yōu)化提供更為高效且精確的方法，同時(shí)可以為探索加速遺傳算法在該領(lǐng)域的應(yīng)用提供參考。

本文對(duì)淮北平原區(qū)免耕農(nóng)田不同深度的原狀土壤進(jìn)行取樣，用激光法和壓力膜儀法測(cè)定了土壤的機(jī)械組成和水分特征曲線，分析其分布變化規(guī)律，并利用實(shí)測(cè)數(shù)據(jù)通過加速遺傳算法對(duì)VG模型的4項(xiàng)參數(shù)進(jìn)行了優(yōu)化率定，對(duì)土壤脫濕過程進(jìn)行了模擬。研究結(jié)果為淮北平原砂姜黑土土壤水分運(yùn)移以及合理設(shè)計(jì)農(nóng)田灌排工程、田間管理和灌溉制度提供重要的最新理論依據(jù)。

1 材料與方法

1.1 研究區(qū)概況

選取新馬橋農(nóng)水綜合試驗(yàn)站為試驗(yàn)點(diǎn)，試驗(yàn)站地處東經(jīng)117°22′，北緯33°09′，位于皖北平原中南部、安徽省蚌埠市北25 km的固鎮(zhèn)縣新馬橋鎮(zhèn)境內(nèi)，地面高程19.7 m，試驗(yàn)站的地理環(huán)境、自然條件和作物種植等均在淮北平原具有很好的代表性。據(jù)統(tǒng)計(jì)，試驗(yàn)區(qū)多年平均降雨量911 mm、蒸發(fā)量916 mm，每年降水多集中于6-9月，約占全年降水量的60%，降雨多以暴雨形式降落，時(shí)空分布不均、極易形成區(qū)域農(nóng)作物旱澇漬災(zāi)害。試驗(yàn)區(qū)土壤為砂姜黑土，土壤分類中劃分為青黑土，屬于典型的中低產(chǎn)田土壤，其理化性狀均屬不良，質(zhì)地黏重，脹縮率大，滲透性差，易澇易旱。

1.2 試驗(yàn)方法

夏玉米完成播種后(前茬為冬小麥)，在試驗(yàn)站大田內(nèi)選取一點(diǎn)開挖剖面。用厚度3 cm，內(nèi)徑6.1 cm的環(huán)刀取0～10、10～20、20～30、30～40、40～50、50～60 cm深度的6組原狀土，每組取3個(gè)原狀土樣(覆蓋10 cm土層)，以3個(gè)土樣的平均值代表該土層的物理性質(zhì)，同時(shí)每層取混合土樣用于土壤機(jī)械組成分析。

(1)土壤機(jī)械組成測(cè)定方法。風(fēng)干后的土樣過2.5 mm篩；稱取土樣0.5 g，用30%過氧化氫加熱去除有機(jī)質(zhì)；加入六偏磷酸鈉超聲60 s后用激光粒度儀Mastersizer 3000測(cè)量土壤粒徑的體積百分比[11，12]。

(2)土壤水分特征曲線實(shí)測(cè)方法。本文用壓力膜儀法測(cè)量土壤水分特征曲線，壓力膜儀為日本產(chǎn)DIK-3404。將環(huán)刀和原狀土樣下墊濾紙放在壓力膜板上，將其一起浸沒入純凈水中浸泡72 h使其達(dá)到飽和，后放入壓力膜儀中逐次加壓，本文選取10個(gè)壓強(qiáng)點(diǎn)，包括10.00、15.85、39.82、100.03、316.32、631.13、1 000.28、1 585.33、4 080.00、10 200.00(以上單位均為cm水柱/cm)，每次加壓后以排水管不再出水為達(dá)到平衡狀態(tài)標(biāo)志，每次平衡后均需稱出環(huán)刀加濕土的質(zhì)量，實(shí)驗(yàn)最后將土樣及環(huán)刀放入烘箱在105 ℃下烘干，得出干土質(zhì)量。

(3)Van Genuchten 模型。Mualem在1976年提出可以利用土壤水分特征曲線和飽和導(dǎo)水率推導(dǎo)非飽和導(dǎo)水率的模型。1980年Van Genuchten在Mualem等的模型基礎(chǔ)上提出了描述土壤水分特征曲線的新模型[4,7-9]，其表達(dá)形式為：

(1)

式中：θh為土壤體積含水率，cm3/cm3；h為土壤水吸力，cm；θs為土壤飽和含水率，cm3/cm3；θr為土壤殘余含水率，cm3/cm3，實(shí)際應(yīng)用中一般取凋萎含水率[4]，cm3/cm3；α、n為表征土壤水分特征曲線的形狀參數(shù)。

(4)基于加速遺傳算法的參數(shù)率定。以VG模型的土壤飽和含水率θs、土壤殘余含水率θr及形狀參數(shù)α、n為優(yōu)化變量，根據(jù)試驗(yàn)站30多年試驗(yàn)經(jīng)驗(yàn)和前人相關(guān)研究，對(duì)4個(gè)參數(shù)設(shè)定合理的上下限，以擬合值與實(shí)測(cè)值的絕對(duì)誤差和最小為目標(biāo)函數(shù)，采用加速遺傳算法進(jìn)行優(yōu)化求解[13]，為：

(2)

(3)

(5)評(píng)價(jià)指標(biāo)。采用均方根誤差(root mean square error，RMSE)、平均絕對(duì)誤差(mean absolute error，MAE)和平均相對(duì)誤差(mean relative error，MRE)對(duì)基于VG模型和加速遺傳算法的土壤水分特征曲線擬合方法進(jìn)行適用性、精確性評(píng)估[4,7,10,14]。

(4)

(5)

(6)

式中：m為土壤水吸力取值樣本數(shù)。RMSE、MAE以及MRE越小，表明擬合值與實(shí)測(cè)值的偏差越小，估算方法的適用性和精確性越強(qiáng)。

2 結(jié)果與分析

2.1 不同土層干容重和粒徑組成特征分析

郝振純[6]等研究表明國(guó)際制分類比美國(guó)制能更好地反映砂姜黑土土壤剖面質(zhì)地的變異性，因此本研究對(duì)土壤粒徑及質(zhì)地分類均采用國(guó)際制。依照國(guó)際制土壤粒徑分級(jí)標(biāo)準(zhǔn)對(duì)0～60 cm土層的土壤粒徑進(jìn)行分級(jí)，每層的土壤干容重和粒徑分級(jí)數(shù)據(jù)見表1。

表1 不同土層深度土壤干容重及粒徑組成

根據(jù)試驗(yàn)站30多年的試驗(yàn)經(jīng)驗(yàn)和郝振純等的研究表明[6]，翻耕耕作條件下農(nóng)田表層0～10 cm土層土壤干容重最小，其值一般在1.20～1.30 g/cm3之間；而表1中免耕條件下0～10 cm土層土壤干容重為1.57 g/cm3，明顯大于翻耕農(nóng)田土壤干容重，其值與10～30 cm土層土壤干容重相差很小，表明農(nóng)田0～30 cm土層土壤經(jīng)過冬小麥210 d左右的生長(zhǎng)期后，因自然沉降及人類活動(dòng)使得土壤容重明顯增加，土壤變緊實(shí)。0～60 cm土層土壤干容重呈現(xiàn)先逐漸增大，在30～40 cm土層突然明顯變小，然后又逐漸增大的變化規(guī)律，其與30～50 cm土層存在大量砂姜有關(guān)。以上分析表明免耕會(huì)增大土壤表層容重，不利于苗期玉米根系向深層土壤的生長(zhǎng)，增大玉米倒伏風(fēng)險(xiǎn)。

由表1可知，砂姜黑土0～60 cm土層黏粒含量在23.33%～30.57%之間，10～20 cm土層黏粒含量最大，為30.57%。根據(jù)國(guó)際制土壤質(zhì)地分類，0～60 cm土層土壤分類包括粉質(zhì)黏土和粉質(zhì)黏壤土2種，20～40 cm土壤存在明顯分層。本文分析認(rèn)為，由于受人類活動(dòng)影響較重，土壤的翻耕及作物秸稈還田等，導(dǎo)致0～30 cm土層的平均砂粒含量要小于其他土層，進(jìn)而使得其粉粒和黏粒含量亦較高于其他土層。

2.2 不同土層實(shí)測(cè)土壤水分特征曲線規(guī)律分析

土壤水分特征曲線斜率的倒數(shù)即單位基質(zhì)勢(shì)(吸力)的變化引起的含水量變化，稱為比水容量(C)，是分析土壤水分運(yùn)動(dòng)的一個(gè)重要參數(shù)。同一吸力下，土壤的含水率越大，表示其持水性越好，斜率越大，C值越小[6]。由圖1可知，各土層的土壤水分特征曲線變化規(guī)律基本一致，在10～631.13 cm土壤水吸力范圍內(nèi)隨著壓力的增加土壤含水率顯著降低，當(dāng)壓力繼續(xù)增加時(shí)土壤含水率隨著壓力的升高變化較為平緩。由圖1還可看出，10～50 cm土層在相同吸力下土壤含水率逐漸增大，說(shuō)明其持水性隨深度增加而逐漸增強(qiáng)，到50～60 cm土層時(shí)土壤持水性明顯下降；0～10 cm與20～30 cm土層的土壤水分特征曲線基本重合，免耕條件下農(nóng)田表層0～10 cm土壤持水性并非最小；在低吸力時(shí)，0～30 cm土層土壤水分特征曲線基本重合，表明低吸力時(shí)該土層的土壤持水性基本相同。

圖1 0～60 cm分層土壤實(shí)測(cè)水分特征曲線Fig.1 0～60 cm layered soil measured moisture characteristic curve

2.3 不同土層土壤水分特征曲線擬合結(jié)果分析

通過壓力膜儀法測(cè)得0～60 cm土層的分層土壤水分特征曲線，基于實(shí)測(cè)數(shù)據(jù)利用加速遺傳算法對(duì)VG模型進(jìn)行優(yōu)化率定，優(yōu)化參數(shù)值見表2。將優(yōu)化的參數(shù)值代入VG模型，得到不同土壤水吸力下對(duì)應(yīng)的土壤體積含水率，作為擬合含水率，將其與實(shí)測(cè)體積含水率進(jìn)行對(duì)比，結(jié)果見表2和圖2。

表2 加速遺傳算法優(yōu)化VG模型參數(shù)結(jié)果

圖2 0～60 cm分層土壤實(shí)測(cè)與擬合水分特征曲線Fig.2 0～60 cm layered soil measured and fitted water characteristic curve

由表2可知，實(shí)測(cè)含水率和擬合含水率的MAE、MRE和RMSE都很小，各土層的MAE和RMSE均小于0.01 cm3/cm3，MRE均小于1.6%，而付強(qiáng)[4]等利用改進(jìn)螢火蟲算法得到的0～20 cm土層MRE為0.91%，RMSE為1.11 cm3/cm3；由圖2可知，各土層絕大部分實(shí)測(cè)點(diǎn)落在擬合曲線上，近乎重合，綜合以上結(jié)果表明利用加速遺傳算法優(yōu)化率定VG模型后的擬合結(jié)果誤差小精度高，擬合效果好，適用于砂姜黑土土壤水分特征曲線的擬合，所以用優(yōu)化模型參數(shù)的物理含義來(lái)描述實(shí)際的土壤水分特征曲線的特征是可行可靠的。0～20 cm土層在土壤水吸力為10 200.00 cm時(shí)存在較大誤差，這可能與砂姜黑土的脹縮率較大有關(guān)，在不斷加壓失水過程中砂姜黑土體積會(huì)有不同程度縮小，但本文由于試驗(yàn)條件所限未考慮該因素的影響，可能會(huì)導(dǎo)致測(cè)量的土壤體積含水率比實(shí)際值偏大，這也表明0～20 cm土層土壤的脹縮率較其他土層更高。

由表1和表2中數(shù)據(jù)對(duì)比分析可知，土壤干容重越大，飽和含水率越小，凋萎含水率則越大，土層中能被作物吸收利用的水分則越少。本文分析認(rèn)為，由于冬小麥生長(zhǎng)期較長(zhǎng)(210 d左右)，在灌溉、降水及打藥等人類活動(dòng)影響下，土壤經(jīng)長(zhǎng)期的自然沉降后導(dǎo)致農(nóng)田表層土壤緊實(shí)度增加，凋萎含水率增大，土壤有效供水能力變差，增大了玉米的干旱風(fēng)險(xiǎn)。

3 結(jié) 語(yǔ)

(1)利用加速遺傳算法優(yōu)化率定VG模型后的擬合結(jié)果誤差小精度高，擬合效果好，用優(yōu)化后模型參數(shù)的物理含義來(lái)描述砂姜黑土實(shí)際的土壤水分特征曲線的特征是可行可靠的。

(2)經(jīng)冬小麥210 d左右的生長(zhǎng)期后，在灌溉、降水和打藥追肥等人類活動(dòng)影響下，土壤經(jīng)長(zhǎng)期的自然沉降后導(dǎo)致農(nóng)田表層土壤緊實(shí)度增加，干容重增大，飽和含水率變小，凋萎含水率增大，土壤有效供水能力變差，增大了玉米的干旱風(fēng)險(xiǎn)。農(nóng)田表層土壤緊實(shí)度的增加也會(huì)造成玉米苗期根系向深層土壤的生長(zhǎng)，增加后期的倒伏風(fēng)險(xiǎn)。因此本文認(rèn)為玉米播種前宜對(duì)農(nóng)田土壤進(jìn)行深翻處理。

(3)本文利用壓力膜儀測(cè)量土壤水分特征曲線過程中未考慮砂姜黑土的脹縮率大這一因素，這可能是導(dǎo)致0～20 cm土層在水吸力為10 200.00 cm時(shí)實(shí)測(cè)點(diǎn)與擬合曲線存在較大誤差的原因，因此在今后的土壤水分特征曲線實(shí)測(cè)

中必須要對(duì)砂姜黑土實(shí)際體積進(jìn)行測(cè)量，依此才能得到更為真實(shí)可靠的實(shí)測(cè)數(shù)據(jù)。