棚內(nèi)微咸水覆膜滴灌下帶間土壤水鹽分布特征模擬

陸培榕 邢瑋麟 楊玉杰 劉文龍 羅 紈

(揚(yáng)州大學(xué)水利科學(xué)與工程學(xué)院, 揚(yáng)州 225009)

0 引言

覆膜滴灌是將具有保水保墑作用的地表覆膜與局部控制性供水的滴灌方式相結(jié)合的高效節(jié)水技術(shù)[1]。由于覆膜限制了土壤與外部大氣的接觸,降低了灌后土壤的水分蒸發(fā)及溶質(zhì)上溯,使得即便在濕潤(rùn)淋洗程度有限的滴灌模式下,微咸水也能作為有效的灌水資源[2-3]。然而,滴灌屬于典型的點(diǎn)源入滲,灌水受土壤基質(zhì)勢(shì)的影響將以濕潤(rùn)體的形式分布于土壤內(nèi),在蒸發(fā)過(guò)程中,微咸水輸入的鹽分將滯留于濕潤(rùn)體外側(cè),并隨灌水次數(shù)增加而逐漸積聚[4],導(dǎo)致土壤中水分和鹽分的動(dòng)態(tài)分布過(guò)程存在明顯的不一致性,當(dāng)涉及地表覆膜時(shí),土表蒸發(fā)與非蒸發(fā)面交替存在,水鹽動(dòng)態(tài)運(yùn)移的復(fù)雜程度將進(jìn)一步加劇[5],尤其在缺乏雨水和漫灌淋洗的設(shè)施大棚內(nèi),鹽分的積聚將存在局部疊加,對(duì)作物生長(zhǎng)造成脅迫[6]。

現(xiàn)階段,隨著農(nóng)業(yè)種植模式的多樣化以及設(shè)施大棚內(nèi)滴灌布設(shè)模式的集約化,針對(duì)覆膜滴灌條件下的水鹽運(yùn)移研究由單點(diǎn)源入滲擴(kuò)展至多點(diǎn)源交匯。研究的對(duì)象以多點(diǎn)源入滲形成的交匯型濕潤(rùn)體幾何形態(tài)、收縮及擴(kuò)張速率、脫鹽及積鹽區(qū)的劃分為主?；⒛憽ね埋R爾白等[7]以室內(nèi)砂壤土土槽試驗(yàn)為基礎(chǔ),研究了滴頭間距為30 cm時(shí)不同滴頭流量和灌水量下濕潤(rùn)體交匯區(qū)的鹽分分布特性,發(fā)現(xiàn)滴頭下方的脫鹽深度大于交匯區(qū)的脫鹽深度,且脫鹽深度的差異隨滴頭流量的減小而增大。王衛(wèi)華等[8]以采用膜下滴灌的粉壤土棉田為研究區(qū)域,發(fā)現(xiàn)膜下濕潤(rùn)體交匯區(qū)的含水率隨滴頭流量的增大而增大,且兩滴頭間濕潤(rùn)體交匯區(qū)徑向?qū)挾扰c灌水歷時(shí)呈冪函數(shù)關(guān)系。齊智娟等[9]在輕度鹽漬化的玉米田中對(duì)比了滴頭間全膜覆蓋和半膜覆蓋的處理方式,發(fā)現(xiàn)使用微咸水滴灌時(shí)全膜覆蓋處理中兩滴頭間0～30 cm的土層內(nèi)形成了明顯的脫鹽帶,而在半膜覆蓋的條件下土壤脫鹽區(qū)僅存在于膜下局部區(qū)域。為進(jìn)一步可視化、定量化地描述土壤水鹽的分布特征,數(shù)值模擬方法伴隨計(jì)算機(jī)技術(shù)的發(fā)展得以在滴灌濕潤(rùn)體交匯下土壤水鹽動(dòng)態(tài)變化的研究中發(fā)揮作用,并以HYDRUS數(shù)值模型為主取得了一系列成果。例如,王維娟等[10]運(yùn)用HYDRUS模型,對(duì)不同間距下的雙點(diǎn)源滴灌濕潤(rùn)體交匯輪廓進(jìn)行了三維狀態(tài)下的動(dòng)態(tài)模擬,結(jié)果顯示,交匯區(qū)濕潤(rùn)體形狀從初期的1個(gè)半球體逐步演變?yōu)?個(gè)分離的近似半球體。CHEN等[11]模擬研究了間作模式下采用膜下滴灌時(shí)土壤水鹽的運(yùn)移特征,分別對(duì)不同作物的膜下根系區(qū)以及膜間裸地進(jìn)行子區(qū)劃分,分析了各子區(qū)內(nèi)的水分及溶質(zhì)的動(dòng)態(tài)形式以及各子區(qū)間的通量交互,發(fā)現(xiàn)瞬時(shí)的水平向通量交互多發(fā)生于灌水階段,而緩慢的垂向通量交互主要存在于灌后蒸發(fā)階段。GUO等[12]模擬了全膜覆蓋形式下濕潤(rùn)體存在交匯時(shí)棉田內(nèi)土壤水鹽的分布變化,結(jié)果表明,土壤鹽分的垂向分布差異隨灌水次數(shù)的增加而逐步顯現(xiàn),而改變滴頭流量造成的影響在不同生育期內(nèi)表現(xiàn)不一。上述模擬研究表明,覆膜滴灌形成的濕潤(rùn)體交匯下水鹽的分布變化具有瞬時(shí)性和局部性,在灌后蒸發(fā)過(guò)程中的變化則較為緩慢且一致。然而,個(gè)別的灌水蒸發(fā)過(guò)程對(duì)土壤水鹽分布的改變不足以對(duì)長(zhǎng)期的作物生長(zhǎng)造成顯著影響。因此,進(jìn)行針對(duì)濕潤(rùn)體交匯區(qū)域內(nèi)水分總體變化趨勢(shì)以及鹽分累積分布特征的研究將更具實(shí)用意義。

本文以設(shè)施大棚內(nèi)田間小區(qū)試驗(yàn)為基礎(chǔ),結(jié)合HYDRUS模型對(duì)微咸水膜下滴灌時(shí)土壤二維剖面內(nèi)的水鹽動(dòng)態(tài)進(jìn)行數(shù)值模擬。同時(shí),為充分考慮滴頭間濕潤(rùn)體交匯程度在灌水和蒸發(fā)過(guò)程中的空間多變性,采用剖面子區(qū)分析的方式,將存在滴頭的滴灌帶帶間區(qū)域作為主要的分析對(duì)象,探究在不同覆膜寬度及滴頭流量下帶間剖面水鹽分布的動(dòng)態(tài)變化及累積特征,以期為設(shè)施環(huán)境微咸水覆膜滴灌模式下創(chuàng)造適宜作物根系生長(zhǎng)的低鹽環(huán)境提供理論支撐與方案參考。

1 材料與方法

1.1 試驗(yàn)區(qū)概況與試驗(yàn)設(shè)計(jì)

試驗(yàn)于江蘇省揚(yáng)州市沙頭鎮(zhèn)現(xiàn)代農(nóng)業(yè)產(chǎn)業(yè)園(32°17′N(xiāo),119°29′E)內(nèi)的設(shè)施大棚中進(jìn)行。棚內(nèi)0～100 cm土層的質(zhì)地為砂壤土,平均砂粒(粒徑0.05～2 mm)占比為55.36%,粉粒(粒徑0.002～0.05 mm)占比31.78%,粘粒(粒徑0～0.002 mm)占比為12.86%;平均土壤容重為1.41 g/cm3,飽和導(dǎo)水率為47.28 cm/d,飽和含水率和田間持水率分別為0.42 cm3/cm3和0.25 cm3/cm3。試驗(yàn)期間棚內(nèi)地下水埋深的波動(dòng)范圍為1.63～1.97 m。該設(shè)施大棚長(zhǎng)期抽取地下水(平均電導(dǎo)率1.37 dS/m)進(jìn)行滴灌,且缺乏雨水或漫灌的淋洗,導(dǎo)致棚內(nèi)土壤呈輕度鹽漬化,試驗(yàn)前0～100 cm土層的平均含鹽量為1.52 g/kg。

以田間小區(qū)試驗(yàn)為研究基礎(chǔ),在棚內(nèi)劃分 2 m×2 m的地塊,采用“兩膜兩行”滴灌布置模式(圖1),膜寬20 cm,滴頭間距為35 cm,兩條滴灌帶間距 80 cm,由外接馬氏瓶對(duì)小區(qū)內(nèi)的滴灌系統(tǒng)進(jìn)行恒壓供水,并采用壓力補(bǔ)償式灌水滴頭將流量統(tǒng)一控制為1.6 L/h。同時(shí),由于取水井內(nèi)地下水含鹽量存在波動(dòng),故使用產(chǎn)業(yè)園內(nèi)提供的自來(lái)水(平均電導(dǎo)率0.22～0.30 dS/m)與氯化鈉晶體(分析純)配比質(zhì)量濃度固定為3 g/L的鹽溶液進(jìn)行灌水。各小區(qū)中均嵌入張力計(jì),用于監(jiān)測(cè)地表以下深度15 cm處的土壤水勢(shì),試驗(yàn)設(shè)計(jì)兩組灌水處理T20和T40,分別對(duì)應(yīng)上一輪灌水后土壤基質(zhì)勢(shì)降至-20 kPa和-40 kPa時(shí)進(jìn)行下一輪滴灌的灌水制度,兩組處理每輪的灌水持續(xù)時(shí)間均為2 h。此外,每組處理重復(fù)3次,在棚內(nèi)共計(jì)劃分6塊滴灌小區(qū)。試驗(yàn)于2021年7月中旬開(kāi)始至10月下旬結(jié)束,總時(shí)長(zhǎng)為90 d,為避免作物根系隨機(jī)性的擴(kuò)張對(duì)土壤水鹽動(dòng)態(tài)分布的影響,試驗(yàn)期內(nèi)未涉及作物的種植。試驗(yàn)過(guò)程中,設(shè)施大棚僅開(kāi)啟兩側(cè)部分的通風(fēng)口,無(wú)降雨或其他形式的灌水。

圖1 試驗(yàn)小區(qū)覆膜滴灌布設(shè)形式示意圖

1.2 觀測(cè)指標(biāo)與換算方法

在每輪灌水前選用內(nèi)徑為2 cm的土鉆對(duì)與滴頭水平距離10、20、30、40 cm處,深度10、20、30、40 cm處的點(diǎn)位進(jìn)行取樣,用于測(cè)定含水率及電導(dǎo)率。由于小區(qū)中各滴頭對(duì)應(yīng)的流量及灌水量相同,故每輪取樣時(shí)選取不重復(fù)的滴頭對(duì)相應(yīng)的位置進(jìn)行取樣。每次取樣完畢后,使用試驗(yàn)小區(qū)內(nèi)的土壤對(duì)取土孔進(jìn)行回填,防止大孔隙優(yōu)先流對(duì)灌水及蒸發(fā)過(guò)程中土壤水分及溶質(zhì)運(yùn)移的影響。采用干燥法測(cè)定不同取樣點(diǎn)位的土壤質(zhì)量含水率,并根據(jù)土壤容重?fù)Q算成相應(yīng)的體積含水率。部分取回的土樣經(jīng)風(fēng)干、碾磨、過(guò)篩(孔徑1 mm)后,與蒸餾水混合(土水質(zhì)量比1∶5),制備土壤浸提液。采用電導(dǎo)率儀(DDBJ-350型,上海雷磁創(chuàng)益儀有限公司)室溫(20℃)下測(cè)定浸提液的電導(dǎo)率,并與蒸發(fā)結(jié)晶法進(jìn)行比對(duì),得出土壤含鹽量與電導(dǎo)率的關(guān)系式為

S=4.12EC1∶5+0.23

(1)

式中S——土壤含鹽量,g/kg

EC1∶5——土水質(zhì)量比為1∶5的土壤浸提液電導(dǎo)率,dS/m

土壤潛在蒸發(fā)量通過(guò)小型蒸發(fā)皿(DF-AM3型,內(nèi)徑20 cm,北京東方鑫鴻科技有限公司)的蒸發(fā)量進(jìn)行估測(cè)。根據(jù)已有研究理論,當(dāng)土壤含水率較高(大于70%田間持水率)時(shí),土壤處于非限制性蒸發(fā)狀態(tài),相應(yīng)的蒸發(fā)量與大氣蒸發(fā)力成正比,可用水面蒸發(fā)量乘以相應(yīng)的系數(shù)進(jìn)行估測(cè)[13-14]。然而,在實(shí)際田塊中,維持大體積土壤的高含水率,并持續(xù)地監(jiān)測(cè)蒸發(fā)造成的水分耗散量,實(shí)現(xiàn)起來(lái)較為困難。故采用100 mL圓柱形環(huán)刀進(jìn)行土壤取樣(20組),在環(huán)刀內(nèi)飽和后底部加蓋,埋置于田塊中,并保證頂部與表土齊平。每隔3 h對(duì)環(huán)刀進(jìn)行稱(chēng)量,并換算成相應(yīng)的蒸發(fā)量,直至環(huán)刀內(nèi)土壤含水率小于田間持水率時(shí)停止觀測(cè)。隨后,將觀測(cè)時(shí)段內(nèi)土樣的蒸發(fā)量與棚內(nèi)的小型蒸發(fā)皿對(duì)應(yīng)的蒸發(fā)量建立擬合關(guān)系式

Es=0.63Epan

(2)

式中Es——棚內(nèi)土壤潛在蒸發(fā)量,mm

Epan——小型蒸發(fā)皿蒸發(fā)量,mm

試驗(yàn)期內(nèi),棚內(nèi)蒸發(fā)皿的日蒸發(fā)量及T20和T40兩組處理的灌水量如圖2所示。

圖2 試驗(yàn)期棚內(nèi)蒸發(fā)皿蒸發(fā)量及T20和T40處理灌水量

1.3 HYDRUS-2D模型構(gòu)建

1.3.1模型原理及基本方程

HYDRUS-2D是以有限元計(jì)算為基礎(chǔ)的數(shù)值模型,可用于模擬二維條件下,變飽和介質(zhì)中水分、熱量及溶質(zhì)的動(dòng)態(tài)分布狀況及邊界處的瞬時(shí)及累積通量[15]。由于滴灌對(duì)應(yīng)的點(diǎn)源入滲將在土壤中形成趨于軸對(duì)稱(chēng)的半橢球狀濕潤(rùn)體[16],故滴灌對(duì)應(yīng)的含水率分布模擬可簡(jiǎn)化至二維垂向剖面內(nèi)?；诖?以土壤均勻且各項(xiàng)同性為基本假設(shè),并忽略土壤水分變化的滯后效應(yīng),模型采用修正后的Richards方程[17]描述二維模式下的土壤水分運(yùn)動(dòng),即

(3)

式中θ——土壤體積含水率,cm3/cm3

h——壓力水頭,cm

x——橫向坐標(biāo),cm

z——縱向坐標(biāo),向上為正,cm

t——時(shí)間,h

K(h)——當(dāng)土壤內(nèi)的壓力水頭為h時(shí)的土壤導(dǎo)水率,cm/h

S(h)——根系吸水項(xiàng),h-1,本研究沒(méi)有涉及作物栽培,故該項(xiàng)忽略

同時(shí),模型采用van Genuchten模型[18]來(lái)描述土壤含水率、基質(zhì)勢(shì)以及導(dǎo)水率的相互關(guān)系,即

(4)

(5)

(6)

(7)

式中θs——土壤飽和含水率,cm3/cm3

θr——土壤殘余含水率,cm3/cm3

α——進(jìn)氣吸力的倒數(shù),cm-1

Ks——土壤飽和導(dǎo)水率,cm/h

Se——土壤相對(duì)飽和度

l——孔隙關(guān)聯(lián)度系數(shù),取0.5[19]

m、n——土壤水分特征曲線(xiàn)形狀系數(shù)

本研究中,灌溉用水均為氯化鈉溶液,未涉及到具有反應(yīng)性溶質(zhì)的化肥或農(nóng)藥的輸入。因此,在模擬過(guò)程中,將土壤中的溶質(zhì)統(tǒng)一簡(jiǎn)化為保守性溶質(zhì),模型選用的保守性溶質(zhì)條件下描述溶質(zhì)運(yùn)移的對(duì)流彌散方程的二維形式為[5]

(8)

式中C——土壤溶液質(zhì)量濃度,g/L

qx、qz——橫向、縱向通量,cm/h

DT、DL——橫向、縱向彌散系數(shù),cm

1.3.2邊界及初始條件設(shè)置

以灌水滴頭正下方0～100 cm的土壤剖面為對(duì)象建立二維模擬域。如圖3所示,模擬域水平長(zhǎng)度與試驗(yàn)小區(qū)一致設(shè)為200 cm,垂向深度為 100 cm。模擬域由約30 000個(gè)三角形有限單元組成,對(duì)存在蒸發(fā)和灌水的上邊界加密節(jié)點(diǎn)的分布,并按照試驗(yàn)取樣點(diǎn)位在模擬域上設(shè)置觀測(cè)點(diǎn)位。模擬域的上邊界可分為覆膜區(qū)、無(wú)膜裸地以及滴頭進(jìn)水區(qū)。模擬時(shí)假設(shè)地膜完全隔絕土壤表面與外部空氣的接觸,故將覆膜區(qū)設(shè)置為零通量邊界。無(wú)膜裸地設(shè)置為大氣邊界,由于棚內(nèi)試驗(yàn)無(wú)降雨及作物蒸騰,因此,模擬域的大氣邊界上僅輸入估測(cè)所得的土壤潛在蒸發(fā)量Es。滴頭進(jìn)水區(qū)的范圍取決于滴頭的流量及土壤的水分特征,已有研究表明,滴灌時(shí)較快的出流會(huì)在土壤表面形成趨于圓形的穩(wěn)定進(jìn)水區(qū)域,采用該進(jìn)水區(qū)域作為入滲邊界,能較為精確地描述二維形式下滴灌形成的土壤水分動(dòng)態(tài)過(guò)程[20-21],相應(yīng)的進(jìn)水區(qū)域半徑可近似表示為

圖3 HYDRUS-2D構(gòu)建的試驗(yàn)小區(qū)土壤剖面模擬域示意圖

(9)

式中Rs——進(jìn)水區(qū)域半徑,cm

q——滴頭流量,cm3/h

進(jìn)水區(qū)設(shè)置為變通量邊界,長(zhǎng)度為2Rs,非灌水時(shí)通量為零,灌水時(shí)通量由滴頭流量和相應(yīng)的通量邊界尺寸確定[11]。試驗(yàn)過(guò)程中,棚內(nèi)各小區(qū)相對(duì)獨(dú)立,除滴灌外無(wú)其他形式的水分輸入,故模擬域的左右兩側(cè)設(shè)置為零通量邊界。此外,在試驗(yàn)過(guò)程中,地下水的觀測(cè)時(shí)間間隔較長(zhǎng),而模擬中以“小時(shí)”為時(shí)間步長(zhǎng),使得邊界的賦值難以動(dòng)態(tài)匹配。結(jié)合現(xiàn)有的研究結(jié)論,當(dāng)?shù)叵滤裆畛^(guò)1 m時(shí),裸地相較于作物耕地的潛水上溯蒸發(fā)量將顯著下降[22-23]。本研究未考慮作物生長(zhǎng),且觀測(cè)深度集中在0～40 cm,故模擬過(guò)程中忽略地下水波動(dòng)對(duì)近地表土壤水鹽分布的影響,將模擬域底邊設(shè)置為自由出流邊界。

模擬域剖面的初始水分及鹽分根據(jù)試驗(yàn)前的實(shí)測(cè)值分層對(duì)應(yīng)設(shè)置。在輸入含鹽量時(shí),土壤含鹽量轉(zhuǎn)換為HYDRUS-2D模型采用的溶質(zhì)液相濃度[24],關(guān)系式為

(10)

式中Sc——土壤體積含水率為θ時(shí)土壤的鹽分質(zhì)量濃度,g/L

ρ——土壤干容重,g/cm3

初始的土壤水分參數(shù)除實(shí)測(cè)所得的飽和含水率θs及飽和導(dǎo)水率Ks,其余參數(shù)均根據(jù)設(shè)施大棚內(nèi)土壤的顆粒級(jí)配及容重通過(guò)HYDRUS內(nèi)置的Rosetta軟件進(jìn)行推演。溶質(zhì)運(yùn)移相關(guān)的縱向彌散系數(shù)DL及橫向彌散系數(shù)DT與模擬域?qū)?yīng)的空間尺寸相關(guān)。參照有關(guān)HYDRUS-2D對(duì)土壤剖面水鹽運(yùn)移的模擬研究[25-26],將DL的初始值設(shè)為模擬域深度的1/10,且DT=DL/10。本研究以小時(shí)(h)為時(shí)間步長(zhǎng)模擬土壤剖面內(nèi)水分及溶質(zhì)在灌水及蒸發(fā)過(guò)程中的運(yùn)移狀況,模擬時(shí)長(zhǎng)為90 d,共計(jì)2 160 h。

1.3.3模型檢驗(yàn)

通過(guò)對(duì)比不同取樣點(diǎn)處水鹽含量的實(shí)測(cè)值與模擬值的差異來(lái)分析所建模型的模擬精度,進(jìn)而對(duì)相關(guān)土壤水分及溶質(zhì)運(yùn)移參數(shù)進(jìn)行率定和驗(yàn)證,并采用均方根誤差(Root mean square error,RMSE)、納什效率系數(shù)(Nash-sutcliffe efficiency coefficient,NSE)、平均相對(duì)誤差(Mean relative error,MRE)作為量化模擬值與實(shí)測(cè)值差異程度的評(píng)價(jià)指標(biāo)。RMSE和MRE越接近0,說(shuō)明模擬值與實(shí)測(cè)值的差異越低;NSE的波動(dòng)范圍為-∞～1,該值越接近1說(shuō)明模擬值與實(shí)測(cè)值的一致性越高,通常當(dāng)NSE大于0.5時(shí),可認(rèn)為模擬結(jié)果滿(mǎn)足精度要求[27]。

2 結(jié)果與分析

2.1 模型率定與驗(yàn)證

采用T20處理?xiàng)l件下與滴頭水平距離為10、20、30、40 cm處0～40 cm土層內(nèi)每輪灌水前土壤的體積含水率和電導(dǎo)率EC1∶5對(duì)模型參數(shù)進(jìn)行率定,隨后,采用T40處理下相應(yīng)的觀測(cè)值對(duì)率定后的參數(shù)進(jìn)行驗(yàn)證。如圖4所示,在率定和驗(yàn)證過(guò)程中,T20和T40處理對(duì)應(yīng)的土壤體積含水率及電導(dǎo)率的模擬值和實(shí)測(cè)值基本集中分布在1∶1線(xiàn)附近。T20和T40處理的含水率模擬值和實(shí)測(cè)值的擬合線(xiàn)斜率均大于1,而土壤電導(dǎo)率模擬值和實(shí)測(cè)值的擬合線(xiàn)斜率略小于1,說(shuō)明兩組處理中含水率實(shí)測(cè)值(T20:0.18～0.27 cm3/cm3;T40:0.19～0.23 cm3/cm3)的波動(dòng)范圍均大于模擬值(T20:0.20～0.25 cm3/cm3;T40:0.15～0.26 cm3/cm3),而土壤電導(dǎo)率實(shí)測(cè)值(T20:0.17～0.62 dS/m;T40:0.19～0.67 dS/m)的波動(dòng)范圍略高于模擬值(T20:0.16～0.56 dS/m;T40:0.21～0.62 dS/m)。表1中,率定和驗(yàn)證階段土壤含水率RMSE波動(dòng)范圍為0.01～0.02 cm3/cm3,NSE為0.52～0.63,MRE為4.39%～8.53%;土壤電導(dǎo)率對(duì)應(yīng)的RMSE、NSE和MRE在相同取樣位置處較含水率均有所提升,波動(dòng)范圍分別為0.03～0.07 dS/m,0.53～0.82以及8.36%～13.42%。同時(shí),RMSE和MRE隨與滴頭水平距離的增加而有所提升,且NSE隨之減小,說(shuō)明土壤水分和鹽分的模擬精度隨與滴頭水平距離的增加而逐漸下降?？傮w而言,RMSE接近于0,MRE較低,NSE均大于0.5,表明在構(gòu)建的二維模擬域中,所輸入的土壤水分及溶質(zhì)運(yùn)移參數(shù)能夠較為準(zhǔn)確地模擬微咸水膜下滴灌時(shí)土壤水鹽的動(dòng)態(tài)變化。經(jīng)率定及驗(yàn)證,選定的土壤殘余含水率θr及飽和含水率θs分別為0.047 cm3/cm3和0.415 cm3/cm3,土壤進(jìn)氣吸力倒數(shù)α為0.02 cm-1,土壤水分特征曲線(xiàn)形狀系數(shù)n為1.5,飽和導(dǎo)水率Ks為1.97 cm/h,孔隙關(guān)聯(lián)度系數(shù)l為0.5,縱向和橫向彌散系數(shù)DL和DT分別為 30 cm和3 cm。

表1 滴頭不同水平位置處深度0～40 cm內(nèi)土壤含水率及電導(dǎo)率模擬精度評(píng)價(jià)指標(biāo)

圖4 T20處理和T40處理下距滴頭不同水平位置處深度0～40 cm內(nèi)土壤含水率及電導(dǎo)率模擬值與實(shí)測(cè)值對(duì)比

2.2 模型應(yīng)用

2.2.1土壤剖面內(nèi)水分及鹽分的變化特征

在試驗(yàn)及模擬過(guò)程中,T20和T40處理單次灌水水量一致,但總灌水次數(shù)不同(圖2)。因此,本文以?xún)山M處理在第5次灌水后的1、3、6 d的模擬剖面含水率分布為對(duì)象,對(duì)比分析在帶間區(qū)域和滴灌帶外側(cè)區(qū)域(滴頭和試驗(yàn)小區(qū)邊界側(cè)之間,以下簡(jiǎn)稱(chēng)“帶外區(qū)”)土壤水分的變化特征。如圖5所示,由于T40處理的灌水時(shí)間間隔長(zhǎng),剖面上含水率的差異較T20處理更大。在灌水后1 d,兩組處理均在滴頭下方形成了含水率較高的圓形濕潤(rùn)體,且濕潤(rùn)體的邊緣存在交匯,使得膜間裸地下方的含水率相應(yīng)提升。在灌后第3天和第6天,受土表蒸發(fā)的影響,覆膜區(qū)和裸地的含水率差異逐漸體現(xiàn),同時(shí),土壤剖面內(nèi)的局部含水率較高的濕潤(rùn)體逐漸消失,剖面內(nèi)水平向的含水率趨于一致。

圖5 T20和T40處理灌水后1、3、6 d土壤剖面含水率分布

在土壤含鹽量變化方面,本文通過(guò)HYDRUS-2D模擬的土壤鹽分質(zhì)量濃度來(lái)反映土壤溶質(zhì)的分布及運(yùn)移狀況,并針對(duì)T20和T40處理在模擬過(guò)程中第30、60、90天的剖面鹽分分布進(jìn)行分析。由 圖6 可知,由于T20處理在相同時(shí)段內(nèi)的淋洗次數(shù)多于T40處理,使得T20處理在整個(gè)剖面的鹽分質(zhì)量濃度均低于T40處理;兩組處理受灌水淋洗作用,鹽分質(zhì)量濃度在10 g/L以下的區(qū)域自滴頭向下逐漸擴(kuò)張;在土表裸地處,土壤內(nèi)溶質(zhì)隨水分蒸發(fā)而逐漸表聚,且?guī)чg區(qū)鹽分的表聚程度和范圍均低于模擬域兩側(cè)的帶外區(qū)?？梢?jiàn),土壤剖面內(nèi)鹽分的分布狀況受滴灌后濕潤(rùn)體交匯及擴(kuò)散的影響,雖然灌水次數(shù)較多的T20處理中剖面輸入的鹽分更多,但灌水次數(shù)較低的T40處理?xiàng)l件下鹽分的表聚程度更為明顯。

圖6 T20和T40處理在第30、60、90天土壤剖面鹽分分布

此外,為進(jìn)一步比對(duì)帶間區(qū)與帶外區(qū)內(nèi)沿水平方向的鹽分動(dòng)態(tài)變化特征,對(duì)T20和T40處理土壤剖面分別距滴頭兩側(cè)10、20、30、40 cm水平位置處0～40 cm的常規(guī)耕作土層的水鹽含量變化進(jìn)行模擬分析。如圖7所示,在土壤含水率變化方面,與滴頭的水平距離越小,灌水前后土壤含水率的波動(dòng)范圍越大,但帶間區(qū)和帶外區(qū)相同水平位置處含水率的差異越小;同時(shí),在不同水平位置處,各處理中每輪灌水前后土壤含水率的變化幅度趨于一致,沒(méi)有明顯的總體上升或下降趨勢(shì)。此外,距滴頭10、20、30 cm水平位置處帶外區(qū)與帶間區(qū)在灌水后含水率的峰值基本相同,而在40 cm位置處T20和T40處理帶間區(qū)的含水率峰值分別比帶外區(qū)高11.55%～17.78%和13.34%～16.52%,說(shuō)明在40 cm處兩滴頭形成的濕潤(rùn)體發(fā)生交匯重合,造成含水率上升。

圖7 與滴頭不同水平位置處0～40 cm土層平均土壤含水率與鹽分質(zhì)量濃度變化曲線(xiàn)

在土壤含鹽量變化方面,HYDRUS-2D模型通過(guò)土壤鹽分質(zhì)量濃度波動(dòng)來(lái)反映觀測(cè)點(diǎn)位含鹽量的變化。T20的帶間區(qū)與帶外區(qū)的含鹽量均低于T40處理,且差異隨水平距離的增加而逐漸擴(kuò)大;兩組處理內(nèi)各水平位置處外側(cè)區(qū)域的含鹽量均高于帶間區(qū),且差異隨灌水次數(shù)的增加而擴(kuò)大。距滴頭10、20、30、40 cm水平位置處T20處理帶間區(qū)的土壤含鹽量總體呈下降趨勢(shì),模擬時(shí)段末的土壤含鹽量均低于初始值。相比之下,在T40處理下帶間區(qū)內(nèi)距滴頭水平位置30、40 cm處土壤鹽分呈總體上升趨勢(shì),尤其在距離40 cm處灌水淋洗后土壤含鹽量的最低值與初始值幾乎一致。上述現(xiàn)象表明,相比于帶外區(qū),濕潤(rùn)體的交匯能夠抑制土壤鹽分的集聚程度,但當(dāng)灌水頻率降低時(shí)帶間區(qū)內(nèi)膜間裸地的積鹽程度隨之上升。

2.2.2情景模擬

為進(jìn)一步探究帶間區(qū)內(nèi)土壤水鹽變化受滴灌強(qiáng)度和覆膜寬度的影響,本研究以率定及驗(yàn)證后的模型參數(shù)和T40處理的灌水制度為基礎(chǔ),利用HYDRUS-2D模型模擬了不同設(shè)計(jì)情景下帶間區(qū)深度0～40 cm范圍內(nèi)土壤水鹽的變化過(guò)程。共設(shè)計(jì)了6組滴頭流量(0.5、1.0、1.5、2.0、2.5、3.0 L/h)遞增的情景,且各流量情景對(duì)應(yīng)的單次灌水量(3.2 L)、模擬期內(nèi)總灌水量(28.8 L)以及灌水含鹽量(3 g/L)均與T40處理保持一致。覆膜寬度的情景同樣為6組,在帶間區(qū)一側(cè)的膜寬分別為15、20、25、30、35、40 cm,相應(yīng)的帶間區(qū)內(nèi)膜間裸地間距分別為50、40、30、20、10、0 cm(圖8)。

圖8 不同覆膜寬度的設(shè)計(jì)情景

(1)帶間區(qū)土壤水鹽總體變化趨勢(shì)

圖9為設(shè)計(jì)的36組模擬情景下帶間區(qū)內(nèi)深度0～40 cm范圍內(nèi)土壤含水率及鹽分質(zhì)量濃度在 90 d 內(nèi)的變化過(guò)程。圖9中,滴頭流量對(duì)帶間區(qū)內(nèi)的含水率影響程度較低,不同覆膜寬度下帶間區(qū)內(nèi)的含水率在每輪灌水前后的變化量基本一致,僅有0.02～0.04 cm3/cm3的差異。帶間區(qū)內(nèi)土壤含水率隨膜間裸地間距的縮小而略有上升,平均含水率最高值(28.76 cm3/cm3)在裸地間距為0 cm的處理中,最低值(25.12 cm3/cm3)在膜間裸地間距50 cm的處理中。相比之下,帶間區(qū)的含鹽量隨灌水次數(shù)的增加呈逐步下降趨勢(shì),且下降程度隨膜間裸地間距的減小而增大,當(dāng)膜間裸地間距由50 cm減至0 cm時(shí),平均土壤鹽分質(zhì)量濃度由9.53 g/L下降至6.25 g/L。此外,不同的滴頭流量?jī)H在灌水間歇期內(nèi)土壤含水率和鹽分質(zhì)量濃度呈輕微差異,所有情景中平均土壤含水率和鹽分質(zhì)量濃度的最大差異僅為0.14 cm3/cm3和0.22 g/L,出現(xiàn)在膜間裸地間距50 cm時(shí)流量為0.5 L/h和3.0 L/h的處理之間。上述模擬結(jié)果表明,帶間區(qū)在多輪滴灌后均能得到鹽分淋洗,且淋洗程度隨地膜間距的增大而增加,但滴頭流量對(duì)水鹽分布的影響相對(duì)較弱。

圖9 不同設(shè)計(jì)情景下帶間區(qū)深度0～40 cm內(nèi)土壤含水率及鹽分質(zhì)量濃度變化曲線(xiàn)

(2)膜間裸地及滴頭流量對(duì)帶間區(qū)土壤鹽分分布的影響

以試驗(yàn)結(jié)束的時(shí)間(第90天)為分析對(duì)象,圖10(圖中不同顏色的實(shí)線(xiàn)為不同滴頭流量下相同膜間裸地間距對(duì)應(yīng)的平均值,實(shí)線(xiàn)周邊的陰影為相應(yīng)的標(biāo)準(zhǔn)偏差)為各模擬情景下土壤鹽分質(zhì)量濃度在深度0～40 cm處帶間區(qū)鹽分的垂向分布形式。由圖10可知,在相同的膜間裸地間距內(nèi)滴頭流量的改變對(duì)鹽分垂向分布造成的差異較小,所有模擬情景中最大偏差出現(xiàn)在深度5 cm處裸地間距為50 cm的處理中,僅為0.046 g/L。此外,除膜間裸地間距為0 cm的處理外,其余情景中土壤的含鹽量均隨深度增加而降低,隨膜間裸地間距的擴(kuò)大而增大。各類(lèi)覆膜寬度條件下土壤鹽分質(zhì)量濃度的差異隨深度的加深而逐漸縮小,在深度 5 cm 處每增加10 cm的膜間裸地,土壤鹽分質(zhì)量濃度將提升約1.73 g/L,而在深度40 cm時(shí),提升幅度降為1.34 g/L。

圖10 不同膜間裸地間距下帶間區(qū)內(nèi)土壤深度0～40 cm土壤鹽分質(zhì)量濃度變化曲線(xiàn)

以帶間區(qū)水平向的中點(diǎn)為原點(diǎn),分析從中點(diǎn)至兩側(cè)滴頭處深度0～40 cm土壤鹽分質(zhì)量濃度在各類(lèi)情景下的分布狀況(圖11)。由圖11可知,不同水平位置處的土壤鹽分質(zhì)量濃度隨膜間裸地間距的增加而降低;除膜間裸地間距為0 cm處理外,土壤鹽分質(zhì)量濃度在不同的覆膜寬度處理下的最大值均出現(xiàn)在中點(diǎn)位置,并向兩側(cè)滴頭處逐漸降低。此外,不同滴頭流量下土壤鹽分的分布趨勢(shì)基本一致,相應(yīng)的差異主要存在于中點(diǎn)處,且鹽分濃度隨滴頭流量的增加而增加。例如,在膜間裸地間距為50 cm處理中,中點(diǎn)處的鹽分質(zhì)量濃度從流量3.0 h/L對(duì)應(yīng)的12.82 g/L逐漸降至流量 0.5 h/L 對(duì)應(yīng)的12.15 g/L。隨著膜間裸地間距的降低,含鹽量將在靠近兩側(cè)滴頭處呈現(xiàn)抬升,并在距兩側(cè)滴頭 10～20 cm 處回落。該現(xiàn)象主要是由帶外區(qū)的覆膜寬度較窄(圖8),集聚的鹽分逐步擴(kuò)散,并推進(jìn)至帶間區(qū)所致。

圖11 滴頭流量和膜間裸地間距對(duì)帶間區(qū)內(nèi)不同水平位置處深度0～40 cm土壤鹽分質(zhì)量濃度的影響

3 討論

采用數(shù)值模型的方式對(duì)土壤水鹽動(dòng)態(tài)進(jìn)行模擬是當(dāng)前分析及評(píng)估滴灌實(shí)施效益的重要方式[28],而所建模型的模擬精度保障是決定模擬研究充分可行的前提條件。在本研究的模型率定與驗(yàn)證過(guò)程中,含水率模擬值的變化量小于實(shí)測(cè)值的變化量,原因是模型假設(shè)土壤的剖面為均質(zhì)且各向同性,而實(shí)地條件下,由于取樣或土壤干縮裂縫等因素容易在土體局部形成優(yōu)先流[29],導(dǎo)致土壤中模擬與實(shí)測(cè)的含水率差異擴(kuò)大。同時(shí),淺層地下水波動(dòng)會(huì)引起上層土壤水分分布的變化,而在模擬域構(gòu)建時(shí)將下邊界簡(jiǎn)化為自由出流邊界,這將導(dǎo)致模擬時(shí)水分的下滲速率過(guò)高,影響了對(duì)水分分布的模擬效果。模擬所得的土壤含鹽量略高于實(shí)際觀測(cè)值,原因可能是在模擬域構(gòu)建過(guò)程中,含鹽量假設(shè)為自上而下均勻線(xiàn)性分布的狀態(tài),與實(shí)際取樣所得的觀測(cè)值之間必然存在差異,進(jìn)而影響后續(xù)的模擬精度。另外,在膜下滴灌時(shí),土壤的鹽分遷移與水分運(yùn)動(dòng)并不同步[8],灌水和蒸發(fā)造成的濕潤(rùn)體交匯與收縮,加劇了土壤水鹽運(yùn)移及分布的復(fù)雜性,使得越接近帶間區(qū)中間部位的模擬結(jié)果誤差越大(表1)?？紤]到相應(yīng)的誤差分析指標(biāo)(RMSE、NSE、MRE)均處于精度可接受的范圍內(nèi),表明驗(yàn)證后的模型能夠有效描述微咸水膜下滴灌過(guò)程中的水鹽分布狀況。

現(xiàn)有的試驗(yàn)及數(shù)值模擬研究均證實(shí),滴灌過(guò)程中濕潤(rùn)體交匯程度越高,對(duì)應(yīng)的帶間區(qū)內(nèi)土壤含水率的提升程度越大[10,30-31]。本研究中,不同覆膜寬度下,每輪灌水過(guò)程中帶間區(qū)內(nèi)含水率的抬升程度趨于一致,說(shuō)明覆膜寬度的變化對(duì)濕潤(rùn)體交匯程度的影響較小。原因可能是灌水歷時(shí)較短,期間的水分蒸發(fā)量不足以體現(xiàn)地膜抑制蒸發(fā)的效果。在鹽分分布方面,SELIM等[4]以?xún)傻晤^間距為40 cm的壤砂土剖面為對(duì)象,同樣利用HYDRUS模型模擬了微咸水(0～2 dS/m)滴灌下水鹽的分布狀況,發(fā)現(xiàn)以滴頭流量1.0 L/h灌水后,將在濕潤(rùn)體交匯區(qū)以下10～40 cm范圍內(nèi)形成了超過(guò)初始含鹽量2～4倍的錐形積鹽區(qū)。然而,在本研究中,滴灌形成的濕潤(rùn)體交匯區(qū)下方并未出現(xiàn)明顯的積鹽區(qū)(圖6),主要原因可能是設(shè)置的土壤初始含水率較低,滴灌使用的微咸水濃度低于土壤的鹽分質(zhì)量濃度;另一方面,在灌水后的蒸發(fā)階段,地表覆膜抑制了土壤水分的散失[32],限制了溶質(zhì)的上溯,減緩了濕潤(rùn)體外圍積鹽區(qū)的擴(kuò)張。在模擬末期,由于帶外區(qū)覆膜寬度較窄,外側(cè)裸地集聚的鹽分向膜內(nèi)擴(kuò)散,使得部分模擬情景中,含鹽量的最低值出現(xiàn)在距滴頭10～20 cm的范圍內(nèi)。該現(xiàn)象說(shuō)明覆膜雖能夠抑制蒸發(fā)導(dǎo)致的垂向鹽分上溯,但并不能限制鹽分在土壤剖面內(nèi)的側(cè)向擴(kuò)散,尤其在淋洗強(qiáng)度和范圍有限的滴灌模式下,鹽分的空間分布由早期的“膜外表聚型”逐漸向后期的“膜內(nèi)底聚型”轉(zhuǎn)變[33]。

灌水量一致的情況下,不同滴頭流量下帶間區(qū)內(nèi)的總體含水率變化幅度的差異較低,且水分的差異將隨膜間裸地面積的降低而進(jìn)一步縮小;土壤鹽分的差異則主要集中在膜間裸地內(nèi),且積聚程度隨滴頭流量的增大而增加(圖11),原因是滴頭流量越大形成的濕潤(rùn)體深度范圍越小,但水平向范圍越大[34],使得滴頭流量較大的處理中土壤鹽分的淋洗深度有限。同時(shí),在點(diǎn)源入滲過(guò)程中,土壤中濕潤(rùn)體的擴(kuò)張速率和含水率隨與滴頭距離的增加而降低[35],即便當(dāng)兩滴頭形成的濕潤(rùn)體存在交匯時(shí)土壤水平向的淋洗程度也相對(duì)有限,導(dǎo)致灌水間歇期內(nèi)滴頭流量較高處理的蒸發(fā)返鹽程度高于滴頭流量較低的處理。另外,本研究的灌水制度是根據(jù)設(shè)計(jì)的土壤含水率閾值而制定的,由于研究中沒(méi)有考慮根系吸水的情況,土壤達(dá)到含水率閾值主要依靠灌后的表層蒸發(fā),這將導(dǎo)致灌水間歇期相對(duì)延長(zhǎng),加劇了蒸發(fā)造成的土壤鹽分集聚,且模擬過(guò)程中土壤的潛在蒸發(fā)是根據(jù)小型蒸發(fā)皿估算所得,造成與實(shí)際耕種時(shí)土壤水鹽環(huán)境的偏差。同時(shí),帶間區(qū)內(nèi)垂向及水平向的水分及溶質(zhì)通量交換有待進(jìn)一步細(xì)化,并加強(qiáng)對(duì)淺層地下水波動(dòng)形式的高頻動(dòng)態(tài)監(jiān)測(cè),從而提升模型的實(shí)用性以及棚內(nèi)微咸水覆膜滴灌方案設(shè)計(jì)的合理性。

4 結(jié)論

(1)以HYDRUS-2D軟件為基礎(chǔ)構(gòu)建了“兩膜兩行”的微咸水覆膜滴灌二維土壤剖面有限元模擬域,并引用設(shè)施大棚內(nèi)的蒸發(fā)狀況作為大氣邊界條件模擬棚內(nèi)環(huán)境。率定及驗(yàn)證的結(jié)果表明,構(gòu)建的土壤水鹽模型精度可靠,在滴灌帶帶間區(qū)內(nèi)土壤含水率及含鹽量的均方根誤差、納什效率系數(shù)以及平均相對(duì)誤差均在合理的范圍內(nèi)波動(dòng),模擬精度隨與滴頭水平距離的減小而增加,且含水率的模擬精度較含鹽量高。

(2)針對(duì)不同膜間裸地間距及滴頭流量展開(kāi)的情景模擬結(jié)果表明,隨灌水次數(shù)的增加,帶間區(qū)域的土壤含鹽量總體呈下降趨勢(shì),裸地間距越小下降程度越大;當(dāng)灌水量不變時(shí),不同滴頭流量下帶間區(qū)內(nèi)土壤水鹽含量的差異較低,且膜間裸地間距越小差異將進(jìn)一步減小。

(3)多輪微咸水滴灌后,膜間裸地間距越大土壤的近地表積鹽現(xiàn)象越突出,不同滴頭流量造成的土壤鹽分含量差異在膜間裸地內(nèi)較為明顯,且提升滴頭流量將加劇積鹽程度;在膜間裸地距離為0～30 cm處理中,帶外區(qū)積聚程度較高的土壤鹽分將向含鹽量較低的帶間區(qū)擴(kuò)散,使得水平向的含鹽量最低值并非出現(xiàn)在淋洗程度最大的滴頭處,而存在于與滴頭水平距離10～20 cm范圍內(nèi)。