基于保水劑應用的土壤水分運移模型構建方法及源匯項求解

廖人寬，楊培嶺(.中國水利水電科學研究院，北京 00048；. 中國農業大學，北京 00083 )

0 引 言

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

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

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

基于此，本文以Richards水分運移方程為基礎，將水分運移時間和保水劑施入土壤中的時間視為同一時間t來建立考慮根系吸水源匯項的水分運移模型，并通過反求迭代法來求解根系吸水速率的數值解，為進一步探索和建立保水劑應用條件下的根系吸水模型提供指導。

1 保水劑應用條件下的土壤水分運移模型

在保水劑應用條件下，因為存在保水劑保水作用效果隨時間變化的因素，土壤水分非飽和擴散率D(θ)和非飽和導水率K(θ)不再只是一個隨含水率θ的改變而變化的函數，而是既隨含水率θ變化，又隨時間t變化，含有時間項t的時變關系函數式D(θ,t)和K(θ,t)。因此，在保水劑應用條件下，考慮根系吸水及土面蒸發的一維垂直非飽和土壤水分運動的定解問題應描述為：

(1)

式中：D(θ,t)為非飽和擴散率時變函數；K(θ,t)為非飽和導水率時變函數；S(z,t)為根系吸水源匯項；E(t)為土面蒸發強度。

2 保水劑應用條件下的土壤水分運移模型的求解

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

此外，在建立保水劑應用條件下的一維垂向非飽和土壤水分運移模型時，假設施入保水劑的土壤混合層為一種性質隨時間發生動態變化的剛性土壤，忽略體積改變帶來的影響；同時，忽略土壤水力參數測定過程的時間對參數變化的影響。

3 保水劑應用條件下的水分運移模型根系吸水源匯項求解

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

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

4 反求估算方法穩定性分析

4.1 數值試驗理論根系吸水模型設定

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

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

(2)

式中：ET為騰發量，mm/d；Z為相對深度，Z=z/Lr(t)；A，B，C為經驗系數，A量綱為1/mm，B，C無量綱；Lr(t)為作物吸水層深度，cm，其表達式為：

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

(3)

式中：T為相對時間，T=t/M；M為玉米生育期天數，d；t為播種時間，d。

(2)數值試驗參數選擇。土壤物理參數及根系參數也選自邵愛軍等的文獻中。已知試驗區田間土壤分為兩層，其中0～70 cm為沙壤土，70 cm以下為粉砂土，0～120cm區域為土壤非飽和帶。沙壤土非飽和導水率為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=54.9θ3.614[圖2(a)]，水分特征曲線為θ=0.50e-0.001 46|h|(0≤|h|≤<168.64),θ=13.893|h|-0.719(|h|>168.64)。根據文獻中已知的土壤水分特征曲線和θ=θr+(θs-θr)/[1+(αh)n]m分別獲取了兩種質地土壤的非飽和導水率公式[圖1(b)和圖2(b)]。

理論根系吸水模型中的參數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本文設置為3 mm/d，迭代控制誤差為1×10-5。

圖1 沙壤土水力特征參數曲線圖

圖2 粉砂土水力特征參數曲線圖

4.2 數值試驗方案及驗證步驟

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

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

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

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

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

(4)

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

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

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

4.3 數值試驗結果

數值試驗測試結果發現，吸水間隔、水力參數、測試誤差和儀器精度對該反求方法穩定性的影響較大，時間步長、層狀土壤及邊界條件對模型穩定性的影響較小，從穩定性、計算量、精確度等方面綜合考慮，建議本文所提出的反求方法適用條件為時間間隔5～17 d，步長為1 000～10 000，誤差控制在0.9以內，儀器精度控制在0.03以內，土面蒸發強度控制在0.6 mm/d以內。

表1 各數值試驗參數選取

5 結 論

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

□

[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] 莊文化, 吳普特, 馮 浩,等. 土壤中施用聚丙烯酸鈉保水劑對冬小麥生長及產量影響[J]. 農業工程學報, 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] 廖人寬, 楊培嶺, 任樹梅. 高吸水樹脂保水劑提高肥效及減少農業面源污染[J]. 農業工程學報, 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] 韓玉國,武亨飛,楊培嶺,等. 保水劑對土壤的物理性質與水分入滲的動態影響[J]. 干旱地區農業研究, 2013,31(5):161-167.

[20] 雷志棟, 楊詩秀, 謝森傳. 土壤水動力學[M]. 北京：清華大學出版社, 1988.

[21] 邵愛軍. 土壤水分運動的數值模擬——以作物根系吸水項為例[J]. 水文地質工程地質, 1996,(2):5-8.

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