徐存東，程 慧，劉璐瑤，謝佳琳，王國霞，田子荀

(華北水利水電大學(xué) 水資源高效利用與保障工程河南省協(xié)同創(chuàng)新中心，鄭州 450045)

灌區(qū)田間灌溉時水鹽運(yùn)移具有動態(tài)的特殊性，在漫灌壓鹽長期作用下，鹽堿成分不斷積聚于田間土壤，從而使土壤逐步演變成鹽漬土，致使干旱灌區(qū)的農(nóng)業(yè)水土環(huán)境不斷惡化[1]。鑒于此，國內(nèi)外相關(guān)學(xué)者針對此問題開展了大量研究。Jiang等研究了灌區(qū)咸水灌溉過程水鹽運(yùn)移，且發(fā)現(xiàn)在長期干旱狀況下幾年后鹽化過程便能達(dá)到平衡態(tài)[2]。Yakirevich等模擬了含水層的鹽分運(yùn)移，量化了灌區(qū)地下水受回流灌溉水鹽堿化作用程度[3]。馮紹元等針對土壤水鹽的微觀運(yùn)移進(jìn)行了轉(zhuǎn)化試驗(yàn)，發(fā)現(xiàn)了田間土壤咸水灌溉水鹽分布的特征[4]。楊金忠等構(gòu)建了區(qū)域多孔介質(zhì)飽和-非飽和溶質(zhì)運(yùn)移模型[5]。這些研究科學(xué)客觀地反映出了田間土壤水鹽的微觀運(yùn)移規(guī)律，為開發(fā)利用干旱灌區(qū)農(nóng)業(yè)水土資源提供了科學(xué)的借鑒意義。本文選取甘肅省景電灌區(qū)為試驗(yàn)研究區(qū)，基于HYDRUS-2D軟件的模擬模型研究了在田間漫灌方式下土壤的水鹽運(yùn)移，探索了干旱區(qū)在灌溉-蒸散發(fā)作用下田間土壤水鹽運(yùn)移的微觀態(tài)勢及規(guī)律。

1 計(jì)算模型

1.1 基本方程

綜合考量土壤飽和帶多相系統(tǒng)，及存在的溶質(zhì)、水分的轉(zhuǎn)化與遷移，同時存在的作物對水分、鹽分的蒸散發(fā)及吸收等過程，建立各試驗(yàn)點(diǎn)的二維垂直灌溉入滲-蒸散發(fā)橫截面模型[6]。

(1)基于水分運(yùn)移的Richards方程。

(1)

式中：S為源匯項(xiàng)；θ為土壤體積含水量;α為水流方向與垂直向夾角；z為原點(diǎn)在地表的空間坐標(biāo)，向下為正；K為水力傳導(dǎo)度;h為土體水負(fù)壓頭；t為時間變量。

(2)基于土壤鹽分運(yùn)移采用CDE對流-彌散方程。

(2)

式中：S為源匯項(xiàng)，即為解吸或吸附項(xiàng);q為水通量；ρ為土壤密度;c為溶質(zhì)濃度;D為彌散系數(shù)，D=D0t+λ|V| 。

(3)基于根系吸水項(xiàng)采用Feddes函數(shù)。

S(h)=α(h)Sp

(3)

式中：Sp為潛在吸水速率最大值；α(h)是土壤水壓力響應(yīng)值；S(h)指單位時間、體積土壤植物根系吸收的水體積。

1.2 模型條件

由于灌區(qū)高蒸發(fā)低降雨的氣候條件，加上采用輪灌灌溉方式，故模型區(qū)上邊界為開放型，可將實(shí)測的蒸散發(fā)強(qiáng)度及灌溉入滲量直接賦值，且將已知的壓力水頭邊界選定為下邊界并直接賦值實(shí)測的地下水埋深[7]。

(1)初始條件。

h(z,t)=h0(z),c(z,t)=c0(z),t=t0

(4)

(2)上邊界條件。在灌溉初期，上邊界處于壓力入滲狀態(tài)，位于土壤的表層產(chǎn)生積水且深度隨時間不斷變化，h(t)為水頭可變的第一邊界條件。

(5)

當(dāng)田間地表灌溉水完全入滲，地表處于蒸散發(fā)狀態(tài)，為蒸散發(fā)強(qiáng)度為E時的第二邊界條件。

(6)

式中：E為表土蒸散發(fā)強(qiáng)度；ti為灌水持續(xù)時間；cR為灌溉水鹽濃度；c0為初始剖面土壤溶質(zhì)的質(zhì)量濃度；ε是縱向水量交換值。

(3)下邊界條件。

(7)

2 應(yīng)用實(shí)例

2.1 研究區(qū)域概況

甘肅省景電灌區(qū)為典型的高蒸發(fā)低降雨大陸干旱性氣候，其地下水深均超過4 m。在種植期內(nèi)灌溉均采用周期為12～15 d的輪灌方式，且采取地面漫灌的方式，因此選取入滲-蒸散發(fā)的模擬模型符合客觀實(shí)際條件。灌區(qū)的作物處于生育期時，土壤的水鹽運(yùn)移最具代表性且最為頻繁，故試驗(yàn)時間可取春小麥拔節(jié)期到孕穗期，即2013年5月10日到2013年6月5日,共26 d。依據(jù)灌區(qū)的農(nóng)業(yè)用水灌溉制度，春小麥從拔節(jié)到成熟一般情況下灌溉定額為3 000 m3/hm2，為達(dá)到灌溉的壓鹽改善效果，本試驗(yàn)灌溉定額選為4 000 m3/hm2，采取灌溉間隔12 d，兩次灌水且每次灌水20 cm深的方式[8]。

選取灌區(qū)內(nèi)地形較平坦、田塊較規(guī)則的典型地塊作為現(xiàn)場試驗(yàn)區(qū)，在該地塊上取邊長1.5 km的等邊三角形作為研究對象，且分別將三頂點(diǎn)處的格田作為灌溉試驗(yàn)點(diǎn)1、2、3，格田面積為0.1～0.15 hm2，同時，格田應(yīng)起埂50 cm以防側(cè)漏。試驗(yàn)點(diǎn)的土壤基本為粉壤土及沙壤土，物理特性各參數(shù)值見表1。此外，土壤內(nèi)鹽分以硫酸鹽-氯化物型為次要構(gòu)成，以硫酸鹽為主要構(gòu)成[9]。選取的試驗(yàn)地塊能較好的代表景電灌區(qū)的水土環(huán)境特征。

表1 土壤水分運(yùn)動特征參數(shù)

2.2 模型及網(wǎng)格劃分

選定試驗(yàn)點(diǎn)田間土體，忽略側(cè)向水分交換，計(jì)算垂直最大剖深Z=350 cm的模擬區(qū)域；利用三角形等間距網(wǎng)格離散化模型區(qū)域，水平方向X軸節(jié)點(diǎn)為51個，垂直方向Z軸節(jié)點(diǎn)為101個，即生成1萬個網(wǎng)格及5 151個節(jié)點(diǎn)；在進(jìn)行迭代計(jì)算時，離散時間單位取小時(h)，最短時間步長設(shè)定為2.4×10-5h；在允許偏差下，土壤含水量迭代精度為1×10-4，壓力水頭迭代精度為1×10-2；迭代次數(shù)最大取20，若迭代次數(shù)超過了20則將新時間步長取原時間步長的1/3重新開始迭代計(jì)算[10]；此外，增大迭代比例取1.3，增大步長迭代信號取3，縮小迭代比例取0.7，縮小步長迭代信號取7。

2.3 確定參數(shù)

分析試驗(yàn)點(diǎn)的土壤物理、化學(xué)指標(biāo)的模擬值和實(shí)測值不同之處，采用Rosetta網(wǎng)絡(luò)神經(jīng)模型，不斷試算調(diào)整參數(shù)值，得到作物的根系吸水特性、溶質(zhì)運(yùn)移特性及土壤的水動力特性3個參數(shù)值，分別見表2～表4。

表2 根系吸水特性參數(shù)

表3 溶質(zhì)運(yùn)移特性參數(shù)

表4 土壤的水動力特性參數(shù)

3 結(jié)果討論

3.1 含鹽量變化討論

以試驗(yàn)點(diǎn)1為分析實(shí)例，研究輪灌方式下土壤鹽分的運(yùn)移動態(tài)。參數(shù)率定模擬田間水鹽的非遲滯，得到土壤在蒸散發(fā)及灌水入滲過程下含鹽量分布，見圖1～圖4。

圖1 第一輪灌水入滲過程土壤含鹽量分布圖

圖3 第二輪灌水入滲過程土壤含鹽量分布圖

圖4 第二輪蒸散發(fā)過程土壤含鹽量分布圖

分析圖 1可知，受蒸散發(fā)影響田間水鹽隨水分上移，致使地表產(chǎn)生積鹽，即灌水前土壤表層產(chǎn)生鹽分峰。灌溉后土壤表層水力梯度及含水率迅速增加，鹽分不斷被水流溶解且隨水流下滲，土層0～10 cm處含鹽量由初始的9.38%降到1 h后7.02%，6 h后降為5.49%，12 h后降為4.36%。經(jīng)過24 h灌水后，積于地面的水落干，土壤表層的含鹽量降到2.89%，其脫鹽率為69.19%，即表明表層土壤的鹽分受灌水淋洗作用明顯。在灌水作用下，土壤鹽分整體呈現(xiàn)明顯的下移趨勢，其鹽分峰值由表層下移為1 h后的9 cm，6 h后15 cm，12 h后20 cm，24 h后23 cm。土壤的含鹽量分布在灌水24 h后實(shí)測值與模擬值變化趨勢基本相同，可設(shè)為蒸散發(fā)階段模擬的第一輪初始條件。

分析圖2可知，在地面水落干后，蒸散發(fā)作用使溶鹽水不斷上移，土壤表層在水分蒸干后稀釋出鹽分，土壤表層含鹽量從灌水2 d后3.67%回升至6 d后4.49%，12 d后6.42%，返鹽效果明顯，其返鹽率為122.15%。經(jīng)過12 d的灌水，土層1～100 cm的不同深度均出現(xiàn)一定程度返鹽現(xiàn)象，但低于各層土壤含鹽量，即表明土壤耕作層受蒸散發(fā)作用返鹽。第一輪灌溉后，試驗(yàn)點(diǎn)1的耕作層出現(xiàn)脫鹽現(xiàn)象，平均脫鹽率20.2%。土壤的含鹽量分布在灌水12 d后實(shí)測值與模擬值變化趨勢基本相同，可設(shè)為灌水階段模擬的第二輪的初始條件。

分析圖3可知，經(jīng)過第二輪灌溉，灌水24 h后，表層土壤初始含鹽量由6.42%降為1.18%，其脫鹽率為81.62%，土壤鹽分峰由表層下移到100 cm處，各層土壤的鹽分均出現(xiàn)下移。土壤含鹽量在灌水24 h后的分布實(shí)測值與模擬值變化趨勢基本相同，可設(shè)為蒸散發(fā)階段模擬的第二輪初始條件。

分析圖4可知，蒸散發(fā)作用使田間土壤0～50 cm處含鹽量從灌水2 d后0.88%升為12 d后1.36%，其返鹽率為54.55%，田間土壤出現(xiàn)積鹽現(xiàn)象。田間土壤的鹽分峰從灌水2 d后110 cm降為12 d 160 cm處，其峰值從灌水2 d后2.62%降為12 d后1.86%，田間土壤鹽分緩慢降低。在第二輪灌溉后，試驗(yàn)點(diǎn)1處土壤出現(xiàn)脫鹽現(xiàn)象，其平均脫鹽率為52.6%。田間土壤含鹽量在灌水12 d后的分布實(shí)測值與模擬值變化趨勢基本相同，即表明實(shí)際現(xiàn)象與模擬結(jié)果基本一致。

綜上所述，試驗(yàn)點(diǎn)1在灌水期土壤鹽分受到淋洗，在輪灌的間歇期又產(chǎn)生返鹽，在兩輪灌溉后，試驗(yàn)點(diǎn)1的田間土壤平均脫鹽率為52.6%，表明土壤耕作層的含鹽量變化受輪灌的淋洗作用明顯。此外，灌區(qū)地下水水位在夏灌期相對較低，未頂托地下水的水鹽運(yùn)移。

3.2 模型討論

將3個試驗(yàn)點(diǎn)得到的平均土壤含水量的實(shí)測值與模擬值進(jìn)行對比，驗(yàn)證所構(gòu)建的數(shù)值模型，驗(yàn)證結(jié)果見圖5。實(shí)測值與模擬值的貼合程度可采用標(biāo)準(zhǔn)差(SD)來衡量，計(jì)算式如下所示[11]。

(8)

式中：i為樣本號；n為觀測樣本數(shù)；X為樣本實(shí)測值；Xi為樣本模擬值。

圖5 土壤平均含鹽量模擬值與實(shí)測值對比

計(jì)算可得土壤含鹽量實(shí)測值與模擬值SD值及各個灌溉階段的脫鹽率模擬值，見表5。分析表5可得，第二輪灌水SD值均小于第一輪灌水，表明第二輪灌水的模擬準(zhǔn)確度更好。究其原因，可能是第一輪灌水時，土壤鹽分尚未受到水的作用鹽分積聚于土壤上層，外界環(huán)境及根系的吸水作用也會造成影響，使第一輪灌水模擬準(zhǔn)確度不高。同時，表5中，土壤含鹽量SD值均處于0.026 5～0.075 1之間，即模擬值與實(shí)測值偏差小，表明該模型擬合效果較好，結(jié)果可靠合理。

表5 實(shí)測、模擬SD值及脫鹽率

4 模型應(yīng)用

應(yīng)用驗(yàn)證后的模型，計(jì)算停灌期典型地塊的土壤平均含鹽量，見圖6。由圖6可知，在停灌后一個月，土壤仍有脫鹽趨勢，但趨勢不明顯，可能是由于土壤中殘余水流的下滲水力梯度和蒸散發(fā)互相作用，造成土壤鹽分下降趨勢放緩。鹽分峰穩(wěn)定在165 cm深處附近，峰值緩慢降低，由停灌1個月后的2.36%降至3個月后的2.03%，土壤含鹽量分布曲線坦化變形，鹽分分布逐漸均勻。持續(xù)的地表蒸散發(fā)及潛水蒸發(fā)導(dǎo)致土壤表層含鹽量持續(xù)上升，由停灌1個月后的1.42%上升至3個月后的1.75%，耕作層積鹽趨勢明顯，即停灌期內(nèi)土壤耕作層出現(xiàn)了積鹽過程。

圖6 典型地塊土壤平均含鹽量分布預(yù)測

5 結(jié) 語

(1)兩次輪灌后，田間地塊中選取的典型試驗(yàn)點(diǎn)土壤的耕作層相繼出現(xiàn)了兩次脫鹽與返鹽，試驗(yàn)點(diǎn)土壤的平均脫鹽率都高于50%，表明在種植期內(nèi)，采用的輪灌方式淋洗脫鹽效果明顯，停灌期內(nèi)土壤耕作層則呈現(xiàn)返鹽積鹽的態(tài)勢。

(2)HYDRUS-2D模型能較好的模擬出土壤微觀水鹽運(yùn)移的變化動態(tài)。夏灌期對小麥進(jìn)行兩次的定額間歇性灌水能較徹底的淋洗積聚于土壤的鹽分，排水積鹽，在一定程度上緩解了土壤的鹽漬化問題。

(3)探索灌區(qū)不同灌溉模式下的土壤水鹽運(yùn)移的特征，尋求農(nóng)田閑置期土壤鹽分控制的灌溉方式，制定合理的灌溉洗鹽制度將會是今后研究的重點(diǎn)。

□

[1] 徐存東. 景電灌區(qū)水鹽運(yùn)移對局域水土資源影響研究[D]. 蘭州：蘭州大學(xué),2010.

[2] Jing Jiang,Shaoyuan Feng,Zailin Huo, et al. Application of the SWAP model to simulate water-salt transport under deficit irrigation with saline water[J]. Mathematical and Computer Modelling,2010,54:902-911.

[3] A Yakirevich, N Weisbrod,M Kuznetsov,et al. Modeling the impact of solute recycling on groundwater salinization under irrigated lands: a study of the Alto Piura aquifer, Peru[J]. Journal of Hydrology,2013,482:25-39.

[4] 林 琳,楊金忠,史良勝,等. 區(qū)域飽和-非飽和多孔介質(zhì)的溶質(zhì)運(yùn)移簡化模型[J]. 水利學(xué)報(bào),2007,38(3):342-348.

[5] 趙志才,馮紹元,霍再林,等. 咸水灌溉條件下土壤水鹽分布特征[J]. 應(yīng)用生態(tài)學(xué)報(bào),2010,21(4):945-951.

[6] 彭建平,邵愛軍. 土壤水鹽運(yùn)移數(shù)值模擬----三峽工程對長江河口地區(qū)土壤水鹽動態(tài)影響評價[J]. 中國農(nóng)村水利水電,2006,(2):50-53.

[7] 劉婷姍,高艷明,李建設(shè). 不同礦化度微咸水入滲下土壤水鹽運(yùn)移特征研究[J]. 灌溉排水學(xué)報(bào),2014,33(3):68-72.

[8] 劉 瑛. 甘肅景電灌區(qū)水資源利用現(xiàn)狀分析[J]. 地下水,2011,33(1):127-128,137.

[9] 楊思存,車宗賢,王成寶,等. 甘肅沿黃灌區(qū)土壤鹽漬化特征及其成因[J]. 干旱區(qū)研究,2014,31(1):57-64.

[11] 李 亮,史海濱,賈錦鳳,等. 內(nèi)蒙古河套灌區(qū)荒地水鹽運(yùn)移規(guī)律模擬[J]. 農(nóng)業(yè)工程學(xué)報(bào),2010,26(1):31-35.