毛俊，伍靖偉*，劉雅文，吳謀松

（1.武漢大學(xué)水資源與水電工程科學(xué)國家重點實驗室，武漢430072；2.南京大學(xué)國際地球系統(tǒng)科學(xué)研究所，南京210046）

0 引言

【研究意義】季節(jié)性凍融土壤在世界廣泛分布，是一種重要的土地利用資源。我國是世界第三大凍土國，季節(jié)性凍融土壤約占全國國土面積的53.5%[1]，主要分布在干旱、半干旱的糧食產(chǎn)區(qū)[2]，同時這些區(qū)域也是我國主要的土壤鹽漬化分布區(qū)[3]。隨著中國工業(yè)的快速發(fā)展和水資源的日益短缺，提升灌區(qū)水資源利用效率勢在必行。而凍融土壤蒸發(fā)作為季節(jié)性凍土區(qū)水文循環(huán)的重要組成部分[4]，在農(nóng)業(yè)用水的管理與制定中扮演著重要角色[5-6]。因此，研究季節(jié)性凍融土壤蒸發(fā)的機(jī)理機(jī)制，對提高水資源利用效率、促進(jìn)農(nóng)業(yè)可持續(xù)發(fā)展具有重要意義。

【研究進(jìn)展】目前已有較多針對季節(jié)性凍融土壤蒸發(fā)的野外觀測和數(shù)值模擬。Zhang 等[7]測量了不同地下水位條件下裸地與非裸地的季節(jié)性凍土蒸發(fā)量；Miao 等[8]對凍融期的潛水蒸發(fā)進(jìn)行了測量；苗春燕等[9]研究了地表覆砂對季節(jié)性凍土蒸發(fā)的影響；Wu 等[10]開展了不同含鹽量及地下水埋深下的野外土柱蒸發(fā)試驗；Zhang 等[11]開發(fā)了可計算寒區(qū)土壤蒸發(fā)的數(shù)值模擬模型；陳軍鋒等[12]使用SHAW 模型模擬了不同潛水埋深下的季節(jié)性凍融土壤蒸發(fā)；解雪等[13]利用灰色關(guān)聯(lián)分析與BP 神經(jīng)網(wǎng)絡(luò)相結(jié)合的方法，對凍融期大田土壤蒸發(fā)量進(jìn)行了模擬預(yù)報?！厩腥朦c】但這些研究均較少考慮鹽分在其中的影響。在季節(jié)性凍土區(qū)，土壤凍結(jié)與融化的交替作用是造成土壤次生鹽漬化的重要原因之一[14]，土壤中的鹽分會顯著影響并改變土壤蒸發(fā)[15-16]。因此，有必要對鹽分影響下的凍融土壤蒸發(fā)規(guī)律進(jìn)行研究。

【擬解決的關(guān)鍵問題】通過野外土柱凍融試驗與數(shù)值模擬相結(jié)合的方法，探究鹽分對季節(jié)性凍融土壤蒸發(fā)的影響規(guī)律，以期為季節(jié)性凍土區(qū)水資源計算及其高效利用提供一定的科學(xué)依據(jù)。

1 材料與方法

1.1 試驗區(qū)概況

試驗區(qū)位于內(nèi)蒙古河套灌區(qū)義長灌域永聯(lián)試驗基地，地理坐標(biāo)為東經(jīng)108°00′6.6″，北緯41°04′2.3″，海拔1020 m。地處干旱、半荒漠地帶，冬季嚴(yán)寒少雪，夏季高溫干熱。試驗區(qū)年降雨量139～222 mm，主要集中在7—8月，占全年降雨量的60%左右；年潛在蒸發(fā)量高達(dá)2200～2400 mm。通常在11月初平均氣溫降低至0℃以下，土壤進(jìn)入初始凍結(jié)階段。隨著氣溫持續(xù)降低，凍結(jié)鋒面向深層土壤推進(jìn)，并在2月下旬達(dá)到最大凍深，多年平均凍結(jié)深度約為100～110 cm。在3月上旬日平均氣溫回升至0℃以上，凍土開始從地表向下消融。與此同時，受地?zé)嶙饔?，凍土開始從凍結(jié)鋒面向上融化，至4月下旬整個土層完全融通，年內(nèi)凍土存在的時間長達(dá)約6 個月。該試驗區(qū)土壤以壤土、粉壤土、粉土和砂土為主，含黏土夾層，土壤平均體積質(zhì)量為1.5 g/cm3，凍融期地下水位埋深及電導(dǎo)率變化如圖1 所示[17]。

圖1 凍融期地下水位埋深與電導(dǎo)率變化Fig.1 Changes of water table depth and conductivity during freezing/thawing period

1.2 試驗設(shè)計與方法

1.2.1 試驗設(shè)計

此次野外土柱試驗于2018年11月16日—2019年5月1日進(jìn)行。設(shè)置2 組試驗土柱，分別記為B1和B2，每組包含4 根處理相同的土柱，雙排放置在大田PE 材質(zhì)的套管中，每根土柱之間間隔為50 cm，其野外布置如圖2 所示。土柱每根長3 m，直徑為5 cm，土柱的前2 m 為非飽和段，后1 m 為飽和段（因為在試驗區(qū)，無灌水和強(qiáng)降雨時，地下水埋深在2 m 左右，地下水的電導(dǎo)率在1500μS/cm 左右[17-18]，后1 m 的飽和段用來模擬地下水位），其土柱設(shè)計如圖3 所示。

在野外試驗大田挖取了3 種含鹽量不同的土壤，分別編號為土壤1、2 和土壤3，其土壤質(zhì)地見表1。將挖取的土壤在室內(nèi)分別進(jìn)行烘干、研磨、過2 mm篩并攪拌均勻，測定土壤1、2 和土壤3 的初始平均電導(dǎo)率分別為2438、12862 和1141 μS/cm（土壤全鹽量分別為7.80、41.16 和3.65 g/kg）。將預(yù)處理的土壤1 和土壤2 按體積含水率0.30 m3/m3進(jìn)行配水，土壤3 按飽和含水率0.43 m3/m3進(jìn)行配水，隨后將配好水的土壤在室內(nèi)用塑料薄膜密封靜置24 h。由于在配水過程中認(rèn)為土壤完全干燥，忽略了土壤吸濕水等的影響，實際配水得到的土壤1、2 和土壤3 的初始含水率分別為0.32、0.33 和0.44 m3/m3。將配好水的土壤3 按體積質(zhì)量1.5 g/cm3每10 cm 填裝在土柱的飽和段，然后將土壤1 和土壤2 按體積質(zhì)量1.5 g/cm3每10 cm分別填裝在B1 和B2 土柱的非飽和段。將填裝好的土柱放入預(yù)先埋置在大田的PE 材質(zhì)套管中進(jìn)行自然凍融試驗，并用保溫棉將套管與土柱之間的空隙進(jìn)行填充。

圖2 野外土柱布置Fig.2 Layout of soil column in field

圖3 土柱設(shè)計示意Fig.3 Diagram of soil column design

表1 土壤質(zhì)地Table 1 Soil texture

1.2.2 土壤含水率和含鹽量測定

每隔42 d 左右對B1 和B2 中的一根土柱進(jìn)行取樣，取樣時間分別為2018年12月25日（P1）、2019年2月10日（P2）、3月20日（P3）、5月1日（P4）。

采用破壞法對土柱進(jìn)行取樣，利用鋸子對土柱每10 cm 進(jìn)行切割，并將切割后土柱中的土倒入到預(yù)先準(zhǔn)備好的自封袋中。在室內(nèi)通過烘干法測定土樣的質(zhì)量含水率，取20 g 烘干后的土樣按照1∶5 的土水質(zhì)量比混合，經(jīng)振蕩、過濾后，所得清液采用上海儀電科學(xué)儀器股份有限公司生產(chǎn)的DDSJ-318 型電導(dǎo)率儀測定其電導(dǎo)率EC1:5（μS/cm），與全鹽量S（g/kg）的關(guān)系可采用在該地區(qū)標(biāo)定得到的換算[19]計算式為：

1.2.3 蒸發(fā)速率計算

在整個凍融期，試驗基地?zé)o降雨降雪，蒸發(fā)成為通過水量平衡原理計算試驗土柱蒸發(fā)速率水分消耗的唯一項，因此其平均蒸發(fā)速率E（mm/d）計算式為：

式中：ΔT為土柱取樣的時間間隔（d）；A為土柱橫截面積（m2）；i為取樣的時期；1000 為單位換算系數(shù)；V為整根土柱水的體積（m3），計算式為：

式中：θij為取樣土層體積含水率（m3/m3）；hj為取樣土層深度，0.1 m；ρw為水的密度（kg/m3）；N為土柱總的取樣土層數(shù)目；j為取樣土層編號。

1.3 考慮鹽分影響的凍融土壤蒸發(fā)模型建立

現(xiàn)有的凍融土壤模型，例如SHAW 模型，土壤蒸發(fā)E計算式[20]為：

式中：ρvs和ρva分別為土壤表面和大氣的水汽密度（kg/m3）；rv為空氣動力學(xué)阻抗（s/m）。該蒸發(fā)模型僅考慮了溶質(zhì)勢對土壤表面水汽密度的影響，并未考慮鹽分的其它作用。例如：①鹽分可影響凍融土壤中的液態(tài)含水率，含鹽量越大，液態(tài)含水率越多[21]，則越有利于土壤蒸發(fā)；②土壤含鹽量達(dá)到一定量時，易形成鹽結(jié)殼，從而阻礙土壤蒸發(fā)[22]?；谝陨峡紤]，本文提出如下考慮鹽分影響的凍融土壤蒸發(fā)模型：

式中：rs和rsc分別表示土壤表面阻抗和鹽分阻抗（s/m）。

Bittelli 等[23]對土壤表面阻抗rs進(jìn)行了研究，推薦使用Sun提出的冪函數(shù)模型計算土壤表面阻抗rs比較合理，即：

式中：θs為土壤飽和含水率（m3/m3）；θl為土壤的液態(tài)含水率（m3/m3）。其中，對于非凍土，液態(tài)含水率即為總的土壤含水率，但是對于凍結(jié)土，其液態(tài)含水率的大小與土壤負(fù)溫和含鹽量有關(guān)。經(jīng)研究發(fā)現(xiàn)，溶質(zhì)勢、基質(zhì)勢與土壤負(fù)溫之間滿足關(guān)系式[24]為：

式中：Lf為水的凍結(jié)潛熱（J/kg）；g 為重力加速度（m/s2）；T為土壤攝氏溫度（℃）；Ts為土壤絕對溫度（K），Ts=T+273.15；ψ為基質(zhì)勢（m），基質(zhì)勢與液態(tài)含水率的關(guān)系可采用VG 模型計算式為：

式中：θr為土壤殘余含水率（m3/m3）；α、m和n分別為影響土壤水分特征曲線的經(jīng)驗系數(shù)，其中m=1-1/n；溶質(zhì)勢π計算式為：

式中：R為氣體常數(shù)，為8.31 J/（mol?K）；c為鹽分濃度（mol/kg），計算式為：

式中：ρb為土壤體積質(zhì)量（kg/m3）；Mm為鹽分的摩爾質(zhì)量（g/mol）。

將式（8）、式（9）和式（10）帶入式（7）得到：

在土壤負(fù)溫T、全鹽量S和土壤水分特征曲線確定時，此時式（11）中只有θl為未知量，可通過牛頓迭代法求解出θl。

對于土壤鹽分阻抗rsc，Haruyuki F 等[25]利用土壤蒸發(fā)實驗推導(dǎo)出rsc可用對數(shù)函數(shù)表達(dá)，其計算式為：

式中：ar和br分別為鹽分控制參數(shù)；Γ為一定厚度土層的鹽分累計值（mg/cm2），其計算式為：

式中：z為一定厚度土層（cm），本文取0.25 cm。

1.4 凍融土壤水熱鹽控制方程

1.4.1 土壤水分控制方程

對于土壤水分運移方程，考慮冰水相變的影響，依據(jù)質(zhì)量守恒定律和達(dá)西定律，可得到控制方程為：

式中：C(ψ)為土壤的容水度；t為時間（s）；z為土壤層深度（m）；K為土壤水力傳導(dǎo)系數(shù)（cm/h）；θi為土壤含冰率；ρi為冰的密度（kg/m3）。

1.4.2 土壤熱控制方程

對于土壤熱控制方程，忽略氣相和液相的對流熱，考慮冰水相變的影響，依據(jù)傅里葉定律和能量守恒定律，可得到土壤熱運移方程為：

式中：Cs為土壤的體積熱容量（J/（m3?℃））；λ為土壤熱傳導(dǎo)系數(shù)（W/（m?℃））。

1.4.3 土壤鹽分控制方程

凍融土壤的鹽分運移控制方程主要考慮了水動力彌散以及對流作用，控制方程為：

式中：DH和Dm分別為水動力彌散系數(shù)和離子擴(kuò)散系數(shù)（m2/s）；C和c分別為單位質(zhì)量土壤和水溶液中鹽分濃度（mol/kg）；ql為液態(tài)水通量（m/s）。

1.5 模型應(yīng)用的情景設(shè)置

為全面揭示不同鹽漬化程度下鹽分對凍融土壤蒸發(fā)的影響規(guī)律，依據(jù)表2 的土壤鹽漬化程度分級標(biāo)準(zhǔn)，對非鹽漬化、輕度和中度鹽漬化土壤分別設(shè)置2種含鹽量水平，對重度鹽漬化土壤設(shè)置3 種含鹽量水平，共9 種含鹽量水平（S1—S9），全鹽量值從小到大分別為0、2.5、5、7.5、10、15、20、35 和50 g/kg。應(yīng)用考慮鹽分影響的凍融土壤蒸發(fā)模型對這9 種含鹽量情形進(jìn)行蒸發(fā)模擬，土壤鹽漬化分級標(biāo)準(zhǔn)參考文獻(xiàn)[26]。

表2 土壤鹽漬化分級標(biāo)準(zhǔn)Table 2 Classification standard of soil salinization

2 結(jié)果與分析

2.1 鹽分對凍融土壤蒸發(fā)的影響

圖4 和圖5 分別為凍融土柱試驗的累計蒸發(fā)量及不同階段的日平均蒸發(fā)速率。由圖4 可知，含鹽量不同導(dǎo)致凍融土壤蒸發(fā)量差異顯著，含鹽量較高的B2土柱總蒸發(fā)量僅為B1 土柱的65%，且凍融土壤蒸發(fā)量主要集中在凍結(jié)初期（P1）以及融通期（P4）。從圖5 可見，在測量的P1 和P4 時期，B1 的日平均蒸發(fā)速率遠(yuǎn)高于B2 土柱，分別約是B2 土柱的1.5 倍和1.8 倍。而在測量的P2 和P3 時期，B1 和B2 的日平均蒸發(fā)速率相當(dāng)，這表明含鹽量差異在該2 個階段的作用不顯著。這可能是因為在該2 個時期，土壤溫度和含水率分別成為了影響蒸發(fā)的主要因素，從而影響到了鹽分的作用。

圖4 凍融土壤累計蒸發(fā)量Fig.4 Accumulated evaporation in freezing/thawing soil

圖5 凍融土壤日平均蒸發(fā)速率Fig.5 Daily average evaporation rate in freezing/thawing soil

2.2 考慮鹽分影響的凍融土壤蒸發(fā)模型驗證與評價

將考慮鹽分影響的凍融土壤蒸發(fā)模型與凍融土壤的水熱鹽運移模型相結(jié)合，以B1 土柱的試驗結(jié)果進(jìn)行參數(shù)率定，以B2 土柱的試驗結(jié)果進(jìn)行模型驗證。率定得到的土壤水分特征曲線如表3 所示，率定得到的鹽分阻抗控制參數(shù)ar=1.05，br=-4.225。模型模擬設(shè)置的上邊界為大氣邊界，土壤溫度的下邊界為給定溫度邊界，水分和鹽分運移的下邊界為零通量邊界。使用的氣象數(shù)據(jù)資料來自內(nèi)蒙古河套灌區(qū)義長灌域試驗站在試驗區(qū)布設(shè)的自動氣象站。

表3 土壤水分特征曲線參數(shù)Table 3 Parameters of soil water characteristic curve

圖6 為率定期和驗證期凍融土壤水鹽運移的模擬值與實測值相關(guān)圖。引入均方根誤差（RMSE）、決定系數(shù)（R2）以及相對百分誤差絕對值的平均值（MAPE）3 種指標(biāo)作為水鹽運移模型率定和驗證的合理性評價指標(biāo)。

率定期凍融土壤水鹽模擬的RMSE分別為0.021 m3/m3和1.89 g/kg，R2為0.87 及0.52，MAPE為5.01%和24.04%；驗證期凍融土壤水鹽模擬的RMSE分別為0.023 m3/m3和10.28 g/kg，R2為0.87 及0.72，MAPE為6.52%和34.76%。模擬數(shù)據(jù)顯示：率定期鹽分的RMSE和MAPE均優(yōu)于驗證期，但是R2卻劣于驗證期，其原因在于R2表示的是整體的數(shù)據(jù)點在回歸線的集中程度，受極端值的影響較小。從圖6 可以看出，驗證期存在幾個較大的極端值，但是其數(shù)據(jù)更加集中在45°線附近，這使得驗證期的鹽分R2高于率定期?？傮w來看，凍融土壤水鹽的模擬值與實測值吻合得較好，這為驗證鹽分影響下的凍融土壤蒸發(fā)模型奠定了基礎(chǔ)。

圖6 土柱水鹽動態(tài)模擬值與實測值Fig.6 Correlation between simulated and measured values of water and salt

圖7 為率定期和驗證期凍融土壤蒸發(fā)模擬值與實測值。從圖7 可看出，SHAW 模型無法對鹽分影響下的凍融土壤蒸發(fā)過程進(jìn)行較好的模擬，在較長時間段低估了B1 土柱蒸發(fā)，而始終高估了B2 土柱蒸發(fā)，模擬結(jié)果的MAPE值分別為26.67%和41.27%，這充分表明SHAW 模型并不適用于鹽漬化季節(jié)性凍土區(qū)的土壤蒸發(fā)模擬。與SHAW 模型的模擬結(jié)果相比，考慮鹽分影響的凍融土壤蒸發(fā)模型明顯改善了對鹽漬化凍融土壤蒸發(fā)的模擬效果，對B1 和B2 土柱均進(jìn)行了較好模擬，模擬結(jié)果的MAPE分別為9.74%和6.35%，這表明建立的考慮鹽分影響的凍融土壤蒸發(fā)模型更加適用于鹽漬化季節(jié)性凍土區(qū)。

圖7 凍融土壤蒸發(fā)模擬值與實測值Fig.7 Simulated and measured values of evaporationfrom freezing / thawing soil

2.3 考慮鹽分影響的凍融土壤蒸發(fā)模型的應(yīng)用

圖8 為9 種含鹽量下模擬得到的凍融土壤在整個凍融期的總蒸發(fā)量。從圖8 可看出，凍融土壤總蒸發(fā)量隨含鹽量的增大表現(xiàn)出先增加后減少的變化趨勢，在土壤含鹽量為S6 時，土壤蒸發(fā)量最大，為136.3 mm。與無鹽土S1 相比，S2—S6 的凍融土壤蒸發(fā)量分別增加了1.09%、13.68%、56.22%、73.73%以及86.46%。而與含鹽土S6 相比，S7—S9 的凍融土壤蒸發(fā)量分別減少了10.20%、26.34%和42.55%。由此可見，鹽分對凍融土壤蒸發(fā)的作用是先促進(jìn)后抑制。

圖8 不同含鹽量下的土壤總蒸發(fā)量Fig.8 Soil total evaporation under different salt contents

3 討論

土壤蒸發(fā)是地表能量和水量平衡重要的組成部分，其準(zhǔn)確計算對于干旱、半干旱季節(jié)性凍土區(qū)土壤鹽漬化防治具有重要的現(xiàn)實意義。而土壤中的鹽分對土壤蒸發(fā)的影響是十分顯著的，這與他人研究結(jié)果一致[16,22]。凍融土壤蒸發(fā)量的大小受到凍融期時段的影響，在凍結(jié)初期和融通期土壤蒸發(fā)量較大，在其余階段土壤蒸發(fā)量較小，這與Kaneko 等[27]研究結(jié)果一致。此外，土壤含鹽量不同導(dǎo)致的土壤蒸發(fā)強(qiáng)度差異也主要表現(xiàn)在凍結(jié)初期和融通期，這與Wu 等[10]研究結(jié)果一致。

進(jìn)一步研究表明，SHAW 模型不能對受鹽漬化影響的凍融土壤蒸發(fā)進(jìn)行較好模擬，而考慮鹽分影響的凍融土壤蒸發(fā)模型則明顯改善了對鹽漬化凍融土壤蒸發(fā)的模擬效果。由此可見，在對含鹽土壤進(jìn)行蒸發(fā)模擬時，考慮鹽分在其中的影響有利于提升模擬計算的精度。張少文等[15]和彭振陽等[16]在對含鹽土進(jìn)行模擬時，通過考慮鹽分對土壤蒸發(fā)的影響提升了模擬精度。并且應(yīng)用該模型發(fā)現(xiàn)，隨著土壤含鹽量的增大，凍融土壤蒸發(fā)呈現(xiàn)先增大后減小的變化趨勢，在土壤含鹽量不高于中度鹽漬化程度時，土壤中的鹽分可促進(jìn)土壤蒸發(fā)，這與Wu 等[10]研究結(jié)果一致。而在重度鹽漬化下，土壤易形成鹽結(jié)殼，對土壤蒸發(fā)具有抑制作用，這與張建國等[22]研究結(jié)果一致。

為應(yīng)對秋澆或冬灌后產(chǎn)生的土壤次生鹽漬化問題，對于排水能力不足的地區(qū)，利用干排水控鹽技術(shù)是治理該問題的重要舉措[28]。該技術(shù)通過地下水的天然運動將耕地多余水分和鹽分排至相鄰鹽荒地，然后通過鹽荒地的蒸發(fā)作用消耗掉，而將鹽分儲存在土壤和地下含水層中。由此可見，鹽荒地的蒸發(fā)能力是干排水控鹽技術(shù)的核心?？墒躯}分對土壤蒸發(fā)能力具有重要影響，隨著干排水控鹽技術(shù)的持續(xù)使用，鹽荒地的土壤含鹽量將持續(xù)累積并達(dá)到重度鹽漬化，此時鹽荒地的土壤蒸發(fā)能力將受到嚴(yán)重限制，通過干排水控鹽防止土壤次生鹽漬化的效果也將大打折扣。因此，對于利用干排水控鹽技術(shù)治理土壤鹽漬化的地區(qū)，鹽荒地的鹽漬化程度應(yīng)控制在中度，可使鹽荒地土壤蒸發(fā)能力達(dá)到最大，也可增強(qiáng)干排水控鹽效果。

值得注意的是，在試驗研究過程中，土柱蒸發(fā)的測量是破壞式的，這可能增加測量的不確定性。為減小蒸發(fā)測量的不確定性，建議后期的土柱試驗可采用原位監(jiān)測的方法。與此同時，野外土柱試驗剖面的水鹽含量是一個動態(tài)變化過程，但是在此次研究過程中采用的是靜態(tài)分析方法，并未突出考慮鹽分的動態(tài)變化過程對土壤蒸發(fā)的影響，這可在后期利用凍融水熱鹽運移模型探究該影響。此外，對于考慮鹽分影響的凍融土壤蒸發(fā)模型中的鹽分阻抗參數(shù)，其值的影響因素也值得進(jìn)一步研究。

4 結(jié)論

1）鹽分對季節(jié)性凍融土壤蒸發(fā)影響顯著，在整個凍融期，含鹽量較高的B2 土柱總蒸發(fā)量僅為B1土柱的65%。在受鹽漬化影響的季節(jié)性凍土區(qū)，考慮鹽分對土壤蒸發(fā)的影響是很有必要的。

2）SHAW 模型對凍融土壤蒸發(fā)試驗?zāi)M結(jié)果較差，率定期和驗證期的MAPE值分別為26.67%和41.27%。而建立的考慮鹽分影響凍融土壤蒸發(fā)模型率定期和驗證期的MAPE值分別為9.74%和6.35%，提升了鹽漬農(nóng)田蒸發(fā)的模擬精度。

3）當(dāng)土壤含鹽量S＜15 g/kg 時，凍融土壤蒸發(fā)隨含鹽量增大而增大，而當(dāng)S＞15 g/kg 時，凍融土壤蒸發(fā)隨含鹽量增大而減小。這表明鹽分對凍融土壤蒸發(fā)隨含鹽量增加呈現(xiàn)出先促進(jìn)后抑制的作用。