廖人寬，楊培嶺(.中國(guó)水利水電科學(xué)研究院，北京 00048；. 中國(guó)農(nóng)業(yè)大學(xué)，北京 00083 )

0 引 言

通過(guò)建立模型來(lái)描述根系對(duì)于土壤水分的吸收動(dòng)態(tài)以及土壤水分的運(yùn)移規(guī)律是一種較為常用的研究方法[1-3]。目前，絕大多數(shù)水分運(yùn)移模型都是基于Richards方程建立的，而考慮根系吸水的土壤水分運(yùn)移模型通常是在Richards方程中加上一個(gè)考慮根系吸水的源匯項(xiàng)(根系吸水速率)來(lái)構(gòu)建的[4-5]。由于根系吸水的不可實(shí)測(cè)性，Richards方程根系吸水源匯項(xiàng)反求估算值作為比較可靠的根系吸水“實(shí)測(cè)值”在根系吸水、吸肥規(guī)律研究方面得到了很大的認(rèn)可。前人所提出的基于兩個(gè)連續(xù)含水率剖面的根系吸水反求解法在小麥根系吸水規(guī)律研究上得到了很好的應(yīng)用[6,7]。

在土壤中應(yīng)用保水劑來(lái)保蓄土壤水分是旱地農(nóng)業(yè)節(jié)水技術(shù)中一項(xiàng)重要的非工程技術(shù)方法。目前，該項(xiàng)技術(shù)已經(jīng)在蘋(píng)果、葡萄、小麥、玉米等作物上進(jìn)行了廣泛應(yīng)用[8-10]，并取得了顯著的節(jié)水和增產(chǎn)效果。在相當(dāng)一部分旱地農(nóng)業(yè)區(qū)，由于水源條件所限，其作物生長(zhǎng)完全依賴于降雨的水分補(bǔ)給。然而，旱地農(nóng)業(yè)區(qū)降雨時(shí)空分布不均和土壤保墑能力差導(dǎo)致雨水的利用效率很低，直接影響了農(nóng)作物的生長(zhǎng)[11,12]。研究表明，在土壤中施入保水劑能夠有效提高旱地農(nóng)業(yè)區(qū)雨水利用效率，其可以反復(fù)吸持進(jìn)入土壤中的降雨，減少水分深層淋溶滲漏量，所吸持的水分隨著土壤干旱程度加劇及作物根系根壓增大會(huì)緩慢釋放以供作物生長(zhǎng)所需，從而保證了作物生育期內(nèi)的持續(xù)供水[13,14]。

然而，大量的研究還顯示[15-17]，保水劑的應(yīng)用不僅僅起到吸水保水持續(xù)供給作物生長(zhǎng)所需水分的作用，其自身獨(dú)特的反復(fù)吸水-釋水功能對(duì)土壤理化性質(zhì)、微生物群落及作物根系生長(zhǎng)等也都會(huì)產(chǎn)生一定的影響，進(jìn)而對(duì)水分在土壤中的運(yùn)移規(guī)律產(chǎn)生周期性和非穩(wěn)定性影響，這也讓保水劑應(yīng)用條件下的土壤水分運(yùn)移模型構(gòu)建變得復(fù)雜和困難。韓玉國(guó)等[18,19]將水分運(yùn)移時(shí)間t和保水劑施入土壤中的時(shí)間T視為兩個(gè)獨(dú)立的時(shí)間，建立了保水劑應(yīng)用條件下的一維土壤水分運(yùn)移模型并獲取了解析解，但該模型只針對(duì)裸土條件下土壤水分運(yùn)移的求解，并沒(méi)有考慮根系吸水源匯項(xiàng)求解的問(wèn)題，不利于深入開(kāi)展保水劑應(yīng)用條件下水分管理制度的研究。

基于此，本文以Richards水分運(yùn)移方程為基礎(chǔ)，將水分運(yùn)移時(shí)間和保水劑施入土壤中的時(shí)間視為同一時(shí)間t來(lái)建立考慮根系吸水源匯項(xiàng)的水分運(yùn)移模型，并通過(guò)反求迭代法來(lái)求解根系吸水速率的數(shù)值解，為進(jìn)一步探索和建立保水劑應(yīng)用條件下的根系吸水模型提供指導(dǎo)。

1 保水劑應(yīng)用條件下的土壤水分運(yùn)移模型

在保水劑應(yīng)用條件下，因?yàn)榇嬖诒Ｋ畡┍Ｋ饔眯ЧS時(shí)間變化的因素，土壤水分非飽和擴(kuò)散率D(θ)和非飽和導(dǎo)水率K(θ)不再只是一個(gè)隨含水率θ的改變而變化的函數(shù)，而是既隨含水率θ變化，又隨時(shí)間t變化，含有時(shí)間項(xiàng)t的時(shí)變關(guān)系函數(shù)式D(θ,t)和K(θ,t)。因此，在保水劑應(yīng)用條件下，考慮根系吸水及土面蒸發(fā)的一維垂直非飽和土壤水分運(yùn)動(dòng)的定解問(wèn)題應(yīng)描述為：

(1)

式中：D(θ,t)為非飽和擴(kuò)散率時(shí)變函數(shù)；K(θ,t)為非飽和導(dǎo)水率時(shí)變函數(shù)；S(z,t)為根系吸水源匯項(xiàng)；E(t)為土面蒸發(fā)強(qiáng)度。

2 保水劑應(yīng)用條件下的土壤水分運(yùn)移模型的求解

對(duì)于(式1)，由于其為非線性方程，一般在進(jìn)行解析求解時(shí)都是先用一個(gè)平均擴(kuò)散度D來(lái)近似代替Richards方程中的非飽和擴(kuò)散率D(θ,t)，再通過(guò)拉氏變換和分部積分法進(jìn)行求解。然而，由于水力參數(shù)的動(dòng)態(tài)變化規(guī)律是體現(xiàn)保水劑作用效應(yīng)的重要特征，若對(duì)其進(jìn)行平均簡(jiǎn)化處理將使得所建立的保水劑應(yīng)用條件下的水分運(yùn)移模型的意義大打折扣。因此，考慮到模型正確物理意義的體現(xiàn)和目前對(duì)于求解水力參數(shù)為時(shí)變函數(shù)的Richards方程解析解的巨大難度，本文采用可以兼顧體現(xiàn)模型正確物理意義且容易求解的隱式差分方法[20]進(jìn)行數(shù)值求解。相較于顯式差分和中心差分格式，隱式差分法在求解時(shí)是無(wú)條件收斂的。而且，本文數(shù)值試驗(yàn)中涉及對(duì)時(shí)間步長(zhǎng)和距離步長(zhǎng)取值等問(wèn)題的探討，采用顯式差分或中心差分格式的話需要時(shí)間步長(zhǎng)和距離步長(zhǎng)滿足一定的要求才能確保方程解的收斂性，會(huì)限制后文數(shù)值試驗(yàn)的探討范圍。

此外，在建立保水劑應(yīng)用條件下的一維垂向非飽和土壤水分運(yùn)移模型時(shí)，假設(shè)施入保水劑的土壤混合層為一種性質(zhì)隨時(shí)間發(fā)生動(dòng)態(tài)變化的剛性土壤，忽略體積改變帶來(lái)的影響；同時(shí)，忽略土壤水力參數(shù)測(cè)定過(guò)程的時(shí)間對(duì)參數(shù)變化的影響。

3 保水劑應(yīng)用條件下的水分運(yùn)移模型根系吸水源匯項(xiàng)求解

根系吸水源匯項(xiàng)求解借鑒Zuo Q等[6,7]提出的迭代反求法，即已實(shí)測(cè)得兩個(gè)連續(xù)的含水率剖面分別為θ(z, 0)和θ(z,T)，假定期間由于根系吸水土壤各深度含水率的減少量為Δθ(z, 0～T)，T表示連續(xù)兩次實(shí)測(cè)值的時(shí)間間隔，則0～T時(shí)段內(nèi)的平均根系吸水速率可概化表示為S(z,0～T)=Δθ(z,0～T)/T。

首先令S=0，根據(jù)初始的含水率實(shí)測(cè)值(t=0時(shí)刻)，采用隱式差分法求解(式1)，算出t=T時(shí)刻第1次迭代的θ(1)(zi,T)值(第1次迭代時(shí)未考慮S的影響，即設(shè)S=0)，把此時(shí)計(jì)算出的Δθ(1)(zi,0～T)=θ(1)(zi,T)-θ(T)(zi,T)除以時(shí)間T后近似為平均根系吸水量S(1)=Δθ(1)/T并帶入原方程繼續(xù)進(jìn)行求解，求解出第2次迭代出的θ(2)(zi,T)值(此后的迭代均考慮S的影響，即設(shè)S≠0)，此時(shí)進(jìn)行判斷Δθ(2)(zi,0～T)=θ(2)(zi,T)-θ(T)(zi,T)是否滿足<誤差控制標(biāo)準(zhǔn)(相對(duì)誤差平方和的平均值)，若滿足則結(jié)束迭代，若不滿足則有S(2)=S(1)+Δθ(2)/T，將其繼續(xù)帶入原方程進(jìn)行迭代求解直至滿足迭代控制標(biāo)準(zhǔn)為止。

4 反求估算方法穩(wěn)定性分析

4.1 數(shù)值試驗(yàn)理論根系吸水模型設(shè)定

(1)理論根系吸水模型。理論根系吸水模型采用邵愛(ài)軍等[21,22]提出的模型，其模型主要為基于蒸騰作用與作物根系吸水關(guān)系而建立。擬采用的根系吸水模型表達(dá)式為：

S=ETA[e-B(lnZ-C)2]/Z

(2)

式中：ET為騰發(fā)量，mm/d；Z為相對(duì)深度，Z=z/Lr(t)；A，B，C為經(jīng)驗(yàn)系數(shù)，A量綱為1/mm，B，C無(wú)量綱；Lr(t)為作物吸水層深度，cm，其表達(dá)式為：

Lr(t)=150(1-e-1.380 813 T)

(3)

式中：T為相對(duì)時(shí)間，T=t/M；M為玉米生育期天數(shù)，d；t為播種時(shí)間，d。

(2)數(shù)值試驗(yàn)參數(shù)選擇。土壤物理參數(shù)及根系參數(shù)也選自邵愛(ài)軍等的文獻(xiàn)中。已知試驗(yàn)區(qū)田間土壤分為兩層，其中0～70 cm為沙壤土，70 cm以下為粉砂土，0～120cm區(qū)域?yàn)橥寥婪秋柡蛶?。沙壤土非飽和?dǎo)水率為D=0.001 48 e15.99θ[圖1(a)]，水分特征曲線為θ=0.48 e-0.001 43|h|(0≤|h|≤129.2)，θ=3.915|h|-0.475(|h|>129.2)；粉砂土非飽和導(dǎo)水率為D=54.9θ3.614[圖2(a)]，水分特征曲線為θ=0.50e-0.001 46|h|(0≤|h|≤<168.64),θ=13.893|h|-0.719(|h|>168.64)。根據(jù)文獻(xiàn)中已知的土壤水分特征曲線和θ=θr+(θs-θr)/[1+(αh)n]m分別獲取了兩種質(zhì)地土壤的非飽和導(dǎo)水率公式[圖1(b)和圖2(b)]。

理論根系吸水模型中的參數(shù)A=8.600 827×10-4-1.189 26×10-3(T-0.591)2，B=1.662 603，C=-1.308 06；M=88 d，ET本文設(shè)置為3 mm/d，迭代控制誤差為1×10-5。

圖1 沙壤土水力特征參數(shù)曲線圖

圖2 粉砂土水力特征參數(shù)曲線圖

4.2 數(shù)值試驗(yàn)方案及驗(yàn)證步驟

本文根系吸水反求估算方法穩(wěn)定性驗(yàn)證主要考慮的影響因素有：①時(shí)間間隔和時(shí)間步長(zhǎng)；②土壤水力特征參數(shù)；③土壤含水率測(cè)量誤差和儀器精度；④土壤層狀性及邊界條件。數(shù)值試驗(yàn)參數(shù)設(shè)置見(jiàn)表1。數(shù)值試驗(yàn)流程如下，其驗(yàn)證的步驟為：

(1)首先輸入水分運(yùn)移參數(shù)及理論根系吸水模型參數(shù)，并確定時(shí)間及空間步長(zhǎng)、儀器精度、邊界條件等反求方法的控制條件；

(2)輸入初始含水率剖面θ(z, 0)和理論根系吸水值S(z, 0)，采用隱式差分法求解水分運(yùn)移方程，得到時(shí)刻t土壤剖面的水分分布理論值θ(z,t)，由于土壤含水率實(shí)際測(cè)試過(guò)程中往往存在誤差，因此本文采用一個(gè)服從正態(tài)分布的誤差向量(vector error-VE)來(lái)對(duì)其進(jìn)行誤差擾動(dòng)(程序中通過(guò)per值的大小來(lái)體現(xiàn))，使其更加接近實(shí)際的情況，即θ*(z,t)=θ(z,t)+誤差向量；

(3)采用代數(shù)多項(xiàng)式方程(式4)對(duì)第(2)步得到的θ*(z,t)進(jìn)行擬合，從而獲得連續(xù)而平滑的θ**(z, t)分布曲線;

θ**(z,t)=r1+r2(z-Z)+r3(z-Z)2+…+rn(z-Z)n-1

(4)

式中：r1，r2，…，rn為擬合參數(shù)；Z為z1，z2，…，zn累加值的均值，擬合從n=3開(kāi)始，直至計(jì)算值和理論值之間的最大絕對(duì)誤差小于誤差控制要求w(w為儀器測(cè)試精度)為止；

(4)采用提出的反求解法，結(jié)合第(3)步得到的θ*(z,t)分布曲線計(jì)算不同時(shí)刻t1(初時(shí))和t2(末時(shí))間的平均根系吸水模擬值；

(5)分析平均根系吸水理論值同模擬值之間的誤差。

4.3 數(shù)值試驗(yàn)結(jié)果

數(shù)值試驗(yàn)測(cè)試結(jié)果發(fā)現(xiàn)，吸水間隔、水力參數(shù)、測(cè)試誤差和儀器精度對(duì)該反求方法穩(wěn)定性的影響較大，時(shí)間步長(zhǎng)、層狀土壤及邊界條件對(duì)模型穩(wěn)定性的影響較小，從穩(wěn)定性、計(jì)算量、精確度等方面綜合考慮，建議本文所提出的反求方法適用條件為時(shí)間間隔5～17 d，步長(zhǎng)為1 000～10 000，誤差控制在0.9以內(nèi)，儀器精度控制在0.03以內(nèi)，土面蒸發(fā)強(qiáng)度控制在0.6 mm/d以內(nèi)。

表1 各數(shù)值試驗(yàn)參數(shù)選取

5 結(jié) 論

本文將保水劑應(yīng)用條件下的土壤水分非飽和擴(kuò)散率及非飽和導(dǎo)水率分別表示為D(θ,t)和K(θ,t)，且本文將保水劑施入土壤中的時(shí)間和試驗(yàn)進(jìn)行時(shí)間視為同一時(shí)間t，解決了之前提出的模型中將保水劑施入土壤中的時(shí)間和試驗(yàn)進(jìn)行時(shí)間視為相對(duì)獨(dú)立而阻礙模型進(jìn)一步發(fā)展的問(wèn)題，本文提出的模型構(gòu)建思路和數(shù)值求解方法為建立保水劑應(yīng)用條件下考慮作物根系吸水的水分運(yùn)移模型提供了基礎(chǔ)條件；提出了保水劑應(yīng)用條件下土壤水分運(yùn)移模型根系吸水源匯項(xiàng)的反求估算方法，并設(shè)置了數(shù)值試驗(yàn)對(duì)該反求方法進(jìn)行了檢驗(yàn)，結(jié)果發(fā)現(xiàn)吸水間隔、水力參數(shù)、測(cè)試誤差和儀器精度對(duì)該反求方法穩(wěn)定性的影響較大，時(shí)間步長(zhǎng)、層狀土壤及邊界條件對(duì)模型穩(wěn)定性的影響較小。

□

[1] Lekakis E H, Antonopoulos V Z. Modeling the effects of different irrigation water salinity on soil water movement, uptake and multicomponent solute transport[J]. Journal of Hydrology, 2016,530:431-446.

[2] Williamson T N, Lee B D, Schoeneberger P J, et al. Simulating soil-water movement through loess-veneered landscapes using nonconsilient saturated hydraulic conductivity measurements[J]. Soil Science Society of America Journal, 2014,78(4):1 320-1 331.

[3] Draye X, Yangmin K, Lobet G, et al. Model-assisted integration of physiological and environmental constraints affecting the dynamic and spatial patterns of root water uptake from soils[J]. Journal of Experimental Botany, 2010,61(8):2 145-2 155.

[4] Javaux M, Schr?der T, Vanderborght J, et al. Use of a three-dimensional detailed modeling approach for predicting root water uptake[J]. Vadose Zone Journal, 2008,7(3):1 079-1 088.

[5] Li H, Yi J, Zhang J, et al. Modeling of soil Water and salt dynamics and its effects on root water uptake in Heihe arid wetland, Gansu, China[J]. Water, 2015,7(5):2 382-2 401.

[6] Zuo Q, Zhang R. Estimating root-water-uptake using an inverse method[J]. Soil Science, 2002,167(9):561-571.

[7] Zuo Q, Meng L, Zhang R. Simulating soil water flow with root-water-uptake applying an inverse method[J]. Soil Science, 2004,169(1):13-24.

[8] Eneji A E, Islam R, An P, et al. Nitrate retention and physiological adjustment of maize to soil amendment with superabsorbent polymers[J]. Journal of Cleaner Production, 2013,52(4):474-480.

[9] Namvar H, Asgharzade A, Babaeian M, et al. The effects of superabsorbent polymer on yield and yield component of two grape varieties.[J]. International Journal of Biosciences, 2014.

[10] 莊文化, 吳普特, 馮 浩,等. 土壤中施用聚丙烯酸鈉保水劑對(duì)冬小麥生長(zhǎng)及產(chǎn)量影響[J]. 農(nóng)業(yè)工程學(xué)報(bào), 2008,24(5):37-41.

[11] Turner N C. Agronomic options for improving rainfall-use efficiency of crops in dryland farming systems.[J]. Journal of Experimental Botany, 2004,55(407):2 413-2 425.

[12] McCutcheon, M. C, Farahani, et al. Effect of soil water on apparent soil electrical conductivity and texture relationships in a dryland field[J]. Biosystems Engineering, 2006,94(1):19-32.

[13] 廖人寬, 楊培嶺, 任樹(shù)梅. 高吸水樹(shù)脂保水劑提高肥效及減少農(nóng)業(yè)面源污染[J]. 農(nóng)業(yè)工程學(xué)報(bào), 2012,28(17):1-10.

[14] Liu F, Ma H, Xing S, et al. Effects of super-absorbent polymer on dry matter accumulation and nutrient uptake of Pinuspinaster container seedlings[J]. Journal of Forest Research, 2013,18(3):220-227.

[15] Han Y G, Yang P L, Luo Y P, et al. Porosity change model for watered super absorbent polymer-treated soil[J]. Environmental Earth Sciences, 2010,61(6):1 197-1 205.

[16] Guilherme M R, Aouada F A, Fajardo A R, et al. Superabsorbent hydrogels based on polysaccharides for application in agriculture as soil conditioner and nutrient carrier: A review[J]. European Polymer Journal, 2015.

[17] Bai W, Zhang H, Liu B, et al. Effects of super-absorbent polymers on the physical and chemical properties of soil following different wetting and drying cycles[J]. Soil Use & Management, 2010,26(3):253-260.

[18] Han Y G, Yu XX, Yang P L, et al. Dynamic study on water diffusivity of soil with super-absorbent polymer application[J]. Environmental Earth Sciences, 2012, 69(1): 289-296.

[19] 韓玉國(guó),武亨飛,楊培嶺,等. 保水劑對(duì)土壤的物理性質(zhì)與水分入滲的動(dòng)態(tài)影響[J]. 干旱地區(qū)農(nóng)業(yè)研究, 2013,31(5):161-167.

[20] 雷志棟, 楊詩(shī)秀, 謝森傳. 土壤水動(dòng)力學(xué)[M]. 北京：清華大學(xué)出版社, 1988.

[21] 邵愛(ài)軍. 土壤水分運(yùn)動(dòng)的數(shù)值模擬——以作物根系吸水項(xiàng)為例[J]. 水文地質(zhì)工程地質(zhì), 1996,(2):5-8.

[22] 邵愛(ài)軍, 李會(huì)昌. 野外條件下作物根系吸水模型的建立[J]. 水利學(xué)報(bào), 1997,(2):68-72.