蒸發條件下寧夏銀北地區鹽漬土水分運移數值模擬

趙文娟+李春光

摘要：對鹽漬土中的水鹽動態進行定量分析是提出適宜于寧夏銀北地區鹽漬土調控措施的基礎工作。蒸發條件下土壤水分運動模擬初步研究表明，將有限體積法應用到求解蒸發條件下的土壤水分運動模擬計算中，模擬值與實測值吻合較好；春季土壤含水率模擬曲線在縱斷面較為一致，地表蒸發量對0.2～0.6 m的土壤含水率有較大的影響。

關鍵詞：土壤水分；蒸發；數值模擬；鹽漬土

中圖分類號：S152.7 文獻標識碼：A 文章編號：0439-8114（2015）15-3638-04

DOI：10.14088/j.cnki.issn0439-8114.2015.15.013

Abstract： Its the basis of the work to use the numerical simulation method for finding appropriate control measures in Yinbei region of Ningxia province. The results showed that there was a good agreement between simulated and measured values when the finite volume method was selected to simulate the movement of soil moisture evaporation. The simulated curve of soil moisture in the horizontal section was consistent with that of the real values. Surface evaporation had a greater impact on the soil moisture of 0.2～0.6 m.

Key words： soil moisture； evaporation； numerical simulation； saline soil

寧夏回族自治區屬于黃河中上游半干旱-半漠境鹽漬區，鹽漬化土地面積為20 718 km2，占自治區總面積的40%。其中尤以銀北地區分布最廣，鹽漬化面積為17 070 km2，約占全區鹽漬化土地總面積的81%左右，該地區含鹽量0.6%～1.0%的重度鹽漬土就占到本區鹽漬化土壤面積25.66%[1]。由于銀北地區地形低洼，排水不暢，地下水位高致使水中的鹽分在強烈的蒸發作用下堆積于土壤內，這是形成土壤鹽漬化的主要原因。通過對該地區的鹽化土地實施水利工程改造和農業生物改良調控措施，鹽漬化得到了整體控制[2]。但近幾年來該地區局部的鹽漬化現象明顯加重，主要是由于打淺井使用高礦化度地下水進行灌溉及淋洗，加劇了土壤的鹽害[3]。這就需要提出更為適宜銀北地區土壤鹽漬化調控措施，作為該地區農業可持續發展和生態環境保護的基礎。

本研究是通過建立能夠描述寧夏銀北地區鹽漬土中水鹽運移的數學模型，采用數值計算方法進行模擬計算，對鹽漬土中的水鹽動態進行定量分析。分析模擬結果，掌握水分、鹽分在鹽漬土體中的基本運動規律，為確定合理的土壤鹽漬化調控模式提供科學依據。

1 材料與方法

1.1 土壤水分運動數學模型

土壤水分運移模型是基于對流擴散方程和初邊值條件共同構建的模型，模型中的各類土壤水動力參數均需要通過試驗來確定。該模型屬于非線性偏微分方程，當前最有效的求解方法是數值計算方法。有限差體積法[4]是當前求解流動問題數值計算方法中最成功的方法，該方法是將方程中的每一項在所劃分的控制體積進行積分，各項積分均有明確的物理意義，構成線性方程組進行求解計算。

1.2.1 有限體積法離散過程 采用有限體積方法對二維非飽和土壤水分運動方程進行線性離散。

其中，？準0為初值含水率（cm3/cm3），D（？準）為非飽和土壤擴散率（cm2/min），為非飽和土壤導水率（cm/min），？著為水量交換強度（cm/d）。

對方程在圖1所示的控制體內積分可得：

由于全隱式格式可以在較大時間步長保持計算結果的穩定性[5]，所以采用全隱格式對擴散項進行時間上的積分，設？駐x=？駐z=？駐d則？駐V=A·？駐x=A·？駐z，得：

1.2.2 邊界條件處理 上邊界條件為地表，與當地的氣象因素相關聯，即與垂直向的水量交換強度？著（cm/d）有關。當？著為負數時為該時期的蒸發強度。根據監測蒸發量數據即可得到整個監測區域內的水量交換強度值。

采用上述相同方法給出邊界條件在計算區域內離散得到各個端點的離散系數，見表1。以上系數通過編寫MATLAB程序進行迭代運算，即可得到蒸發條件下的土壤水分運動曲線。運行前，需給出初始值。假設土壤的含水率隨監測深度的不同而變化，呈現為隨深度變化的曲面，即有：

？準0（SD，0）=HSL （6）

其中，SD為土壤不同深度（m），以向量表示；HSL為不同深度的土壤含水率（cm3/cm3），以向量表示。

1.2 試驗區情況

試驗區位于寧夏銀川市賀蘭縣南梁臺子鐵西村。該村地理坐標為東經105°58，北緯38°41，位于賀蘭縣的西北方向。試驗區為一面積約459 m2的條形區域。距離地表20 cm試驗土壤的平均初始含水率為15.843%，土壤容重為1.6 g/cm3。

2 結果與分析

利用時間延遲反射儀對試驗區的15個以三行五列形式布置的測點進行土壤含水率和電導率分時段、分層的測定工作。其中相關的水動力學參數有：

選擇2012年5月16日10點到5月17日10點的時段對自西向東的由2號、5號、8號、11號和14號、17號、20號及23號組成的縱斷面數值模擬土壤水的運動過程，其他測點作為邊界條件進行模擬計算。當日蒸發所引起的土壤水通量為6.4 cm/d，蒸發強烈。模擬結果表明，同一深度下的土壤含水率自西向東運動趨勢基本保持一致，但模擬數值自西向東逐漸增加，土壤含水率最低點均出現在0.2 m附近，最大值出現在1.0 m附近。模擬值與實測值吻合較好，模擬曲線可以基本表達土壤水分在鹽漬土內的運動趨勢，但個別測點的實測值與模擬值的偏移量較大，可能的原因是假設模擬區域的土壤質地為均質土壤，即土壤容重是相同，在模擬過程出現了偏移現象，另一方面也可能是測量誤差造成實測值與模擬曲線的偏移。

將5月17日的模擬值作為初始值，模擬土壤水分在5 d后的變化過程，擬合曲線見圖3。由于氣象數據的變化，土壤含水率也發生變化。含水率的最小值出現在了0.4 m附近，土壤監測斷面含水率與5月17日的相比是增加的，且測定點的模擬數值自西向東的增量更為明顯。在0.8 m附近土壤含水率逐漸進入穩定狀態，且隨著迭代時間的增長，數值解與實測值吻合較好。說明運用有限體積法計算土壤縱斷面含水率是有效的，模擬曲線也顯示出0.2～0.6 m處土壤含水率較低。

3 結論

采用數值模擬的方法研究蒸發條件下的寧夏銀北鹽漬土區土壤水分運動，在整個模擬過程中，模擬值與實測值吻合較好，將有限體積法應用到求解土壤水分運動方程是可行的。模擬結果表明，土壤的含水率在春季的變化較為頻繁，在同一深度的土壤含水率從數值上呈現西低東高的現象；縱斷面上同一時刻的土壤含水率模擬曲線形態是基本一致的，土壤蒸發量使0.2～0.6 m的土壤含水率較小。

參考文獻：

[1] 何文壽，劉陽春，何進宇.寧夏不同類型鹽漬化土壤水溶鹽含量與其電導率的關系[J].干旱地區農業研究，2010，28（1）：111-113.

[2] 俞仁培.對鹽漬土資源開發利用的思考[J].土壤通報，2001，32S1：139-140.

[3] 史曉杰，萬 力，張永庭，等.銀北地區土壤鹽漬化形成機理與模擬研究[J].水文地質工程地質，2006，34（6）：116-117

[4] VERSTEEG H K，MALALASEKERA W.An introduction to computational fluid dynamics：The Finite V Method[M].London：Longman Group Ltd，1995.

[5] 姜啟源，謝金星，葉 俊.數學模型[M].北京：高等教育出版社，2011.