黑河流域中游灌區含水層系統異質性對地表熱通量和溫度的影響模擬

陸 崢，孫 景 博，何 源，李 宗 超，彭 書 艷，楊 曉 帆*

(1.北京師范大學地理科學學部，地表過程與資源生態國家重點實驗室，北京 100875；2.南方海洋科學與工程廣東省實驗室(廣州)，廣州 廣東 511458；3.中國環境監測總站，北京 100012)

0 引言

地表和大氣間的動量、水熱交換主導陸面過程及其與大氣間的動力、熱力相互作用[1-4]。區域尺度上地氣間的相互作用能影響甚至改變氣候變化響應機制，引發極端干旱和洪澇事件[5-7]。理解并準確表征地氣間的動量、水熱交換規律和傳輸過程對氣候變化、農業、林業、水文和水資源管理等具有重要意義[8,9]。地下含水層系統物理參數(如滲透系數和孔隙度)的異質性對水分的交換和平衡起著重要控制作用，其中滲透系數是描述地下含水層系統空間異質性的重要參數，也是顯著影響土壤和含水層水分、溶質運移模擬的參數之一[10]。研究表明，地表感熱通量和潛熱通量是關于淺層土壤水分的函數[11-13]；此外，地下水位對地表熱通量和溫度也有控制作用[14]，如農場灌溉和地下水開采可能會分別在短期和長期改變地下水位對地表熱通量和溫度的控制模式[15]。但在地下含水層系統滲透系數空間異質性對地表熱通量的影響中，農業灌溉活動的作用尚不明確。

傳統的分布式水文模型與地下水模型主要用于刻畫水文循環現象及其物理過程，基于簡化地下水模塊的陸面過程模型則側重表征地表土壤中及地氣間的動量、水熱交換過程，三者均難以精確描述地下水位和淺層土壤水分與地表熱通量和溫度之間的定量關系。因此，本研究采用地下水—陸面過程耦合模型(Parallel Flow-Common Land Model，ParFlow-CLM)[16]定量描述大氣下邊界層—地表徑流—淺層土壤水分—深層地下水間的轉換關系，以期揭示區域三維地下滲透系數場對地表熱通量和溫度的影響規律；同時考慮黑河流域中游密集農業灌區地表水灌溉和地下水回灌的作用，定量分析灌溉效應產生的淺層土壤水分對地表熱通量的影響，以期為定量分析農業灌溉活動影響下區域滲透系數對地表熱通量的影響提供參考，并為量化不同陸面過程對土壤和含水層系統異質性的響應提供借鑒。

1 研究區與數據

1.1 研究區概況

黑河流域地處西北干旱區，核心區面積約為13萬km2，縱跨青海省、甘肅省和內蒙古自治區，是我國第二大內陸河流域，具有全球獨特的、以水為橋梁和紐帶的多層次自然景觀，對我國中西部水文循環和生態環境影響巨大[17]。本文研究區為黑河流域中游地區(100.05°～100.69°E，38.73°～39.37°N)(圖1)，該區水文地質結構較復雜，屬中溫帶干旱型氣候，下墊面類型主要為農田、草地和裸地，地帶性植被多為溫帶小灌木、半灌木荒漠植被，地勢東北、西南較高，中部為東南—西北走向、地形較平坦的狹長盆地[18-20]；區域年降雨量稀少且集中于6-10 月，平均潛在蒸(散)發量高于1 300 mm，主要水資源為南部祁連山區的冰雪融水和山區降水。特有的水土光熱條件和平坦的地形使區內形成密集農業區，約900條主要灌渠構成錯綜復雜的人工灌溉渠系[21]。自西部大開發戰略實施以來，農業灌溉和工業用水量大幅增長，地下水和地表水綜合利用已成為該區農業和生活用水的主要形式[22]，河、渠、泉等地表水源灌溉及地下水回灌使地下滲透系數空間異質性對地表熱通量和溫度的影響作用更復雜[23]。因此，厘清該區域內地下含水層系統異質性對地表熱通量和溫度的影響，對我國戰略水資源開發和生態環境保護意義重大[18,24]。

圖1 黑河流域中游區域概況Fig.1 Sketch of the middle reaches of the Heihe River Basin

1.2 研究數據

研究數據包括：1)黑河流域大氣驅動數據集[25]，可提供逐小時的氣溫、地表氣壓、水汽混合比、長短波輻射、二維風場和降水數據，原始分辨率為0.05°，采用文獻[26]的方法降尺度至0.005°；2)地表覆蓋類型和高程數據，分別源于黑河流域土地利用覆被數據集[27]和航天雷達地形測繪任務(SRTM)DEM數據集(http://srtm.csi.cgiar.org/srtmdata/)，原始分辨率分別為30 m和90 m，均重采樣至0.005°；3)張掖灌溉渠系數據[28]，利用PriorityFlow[29]將DEM、河道、干渠數據處理為東西、南北向的坡度作為模型輸入參數；4)土壤質地(粒徑組成)空間分布數據集[30]；5)黑河流域30弧秒分辨率地表水及地下水灌溉量數據集[31]，利用黑河流域中游渠道流量測量數據集[32]中記錄的灌溉時間，提取小時尺度的灌溉量；6)黑河流域地表過程綜合觀測網[33]中大滿超級站(100.3722°E,38.8555°N)的渦動相關(EC)和大孔徑閃爍儀(LAS)觀測數據[34,35]。

2 模型與方法

2.1 ParFlow-CLM

ParFlow-CLM是一種開源、基于三維物理過程的地下水—陸面過程水文集成模型[36]，其將飽水帶、包氣帶和地表水視為一個完整的水文連續體，能全面刻畫不同尺度下基巖層到大氣下邊界層間的水文過程和水分、能量平衡[37]。對于地下部分，ParFlow基于三維地下水流動Richards方程求解飽水帶—中間帶—包氣帶的偏微分方程定解[38]，通過van Genuchten方程[39]描述各水動力變量的關系，并利用高效穩定的Newton-Krylov非線性求解器[40]進行矩陣求解；對于地表部分，則耦合CLM (Common Land Model,V2014修改版)[41]，利用二維地表水流動Saint-Venant方程的二階動力波近似解模擬地表徑流等過程，獲得地表下滲、蒸(散)發、植被生理及雪動力等過程中的水熱通量反饋到地下水與地表水的交界面，最終實現垂向的水分和能量平衡。此外，ParFlow-CLM采用的網格為結構化規則網格，可選正交網格或地形網格，模型詳細的求解過程可參考文獻[38,42,43]。

2.2 模型參數與情景設置

本文中模型參數設置(表1)為：水平分辨率為0.005°(約500 m)，垂直方向上網格結構分為多層，孔隙度、van Genuchten四參數(孔隙進氣值參數、孔隙分布參數、殘余含水量、飽和含水量)、曼寧粗糙度均設置為均質，其中孔隙度和van Genuchten四參數取自中國土壤數據集[44]，曼寧粗糙度修改自文獻[45]，模擬區域的底部及四周均設置為零通量邊界。使用2010年1月1日-2011年5月31日的驅動數據對模型進行起轉率定(關閉地下水側向流，計算步長為1 h)后，區域內的地下水儲量變化范圍穩定在5%之內，認為模型達到穩定[46,47]。最后，本文以起轉率定獲得的最終水位作為模型模擬的初始水位，并進行洪水測試和黑河核心干流徑流量測試，結果表明本模型模擬地表水過程的精度較高。

表1 模型參數設置Table 1 Model parameter configurations

本文采用滲透系數表征土壤和含水層系統的異質性，并設計4種情景用于量化從土壤表層到基巖的滲透系數：1)Homogeneity情景，即模型的滲透系數恒定為0.01 m/h；2)GLHYMPS情景，僅在水平方向上設置異質性，滲透系數取自文獻[48]的數據集，垂直方向上采用均質處理；3)CDHSP情景，在水平和垂直方向上均設置異質性，土壤滲透系數來源于中國土壤數據集[44,49]，包氣帶底部到基巖層滲透系數采用最底層土壤層的數值；4)MIX情景，表層土壤到包氣帶底部滲透系數的設置采用CDHSP情景，飽水帶滲透系數的設置采用GLHYMPS情景。

為契合地表觀測數據時段，模擬時段為2012年6月1日-10月31日，步長為1 h。ParFlow-CLM輸出參數有地表感熱通量(H)、潛熱通量(LE)、上行長波輻射、下滲量、總蒸(散)發、土壤蒸發、植被蒸騰、不考慮升華的地表蒸發、雪水當量、地表溫度(LST)、各土壤層溫度(St)等。利用大滿超級站渦動相關(EC)和大孔徑閃爍儀(LAS)觀測數據對感熱通量、渦動相關觀測數據對潛熱通量、地表溫度觀測數據(采用文獻[47]方法由大滿超級站四分量儀器觀測數據計算得到)對LST、地下土壤溫(濕)度儀觀測數據對St進行檢驗和情景對比分析。評價指標為均方根誤差(RMSE)和Pearson相關系數(R)。

3 結果分析

3.1 地表熱通量的驗證與情景對比

如圖2(彩圖見封3)所示，地表感熱通量和潛熱通量模擬結果時序變化幅度均與觀測數據相似，說明4種情景的模擬精度較高[15]，能在不同程度上捕捉到儀器觀測數據的日內及季節變化趨勢，較好地反映黑河流域中游水動力過程與能量傳輸過程。需要說明的是，相比潛熱通量，感熱通量與觀測數據的一致性更高，R值達0.85以上(表2)，原因是潛熱通量僅為點尺度觀測，表征空間尺度的地表參量具有一定局限性。通過與觀測數據對比可知，4種情景的模擬值均高估了感熱通量、低估了潛熱通量，這是由于模型中植被參數為固定值。從圖2c 10月10-20日的數據可以看出，潛熱通量模擬值在午后驟降，導致GLHYMPS情景的模擬結果與觀測數據誤差較小，其他3種情景誤差較大，原因在于長期、持續的灌溉導致淺層土壤處于穩定的高飽和狀態，打破了潛水埋深對潛熱通量的控制，因而在飽和滲透系數情景(GLHYMPS)下模擬結果的誤差較小。同時注意到，在4種情景的感熱通量和潛熱通量模擬中，Homogeneity情景的模擬精度最低，這充分說明在大尺度區域模擬中，考慮表層土壤—包氣帶—飽水帶不同的水力性質并設置垂向和水平向非均質滲透系數是必要的，其能有效降低地表熱通量的模擬誤差。總體而言，Homogeneity情景的峰值較觀測值有一定滯后性，GLHYMPS情景到達峰值后下降時間最早，而CDHSP和MIX情景的模擬值與觀測值一致性很高，說明包氣帶的滲透系數對地表熱通量具有決定性控制作用。

圖2 ParFlow-CLM感熱通量與潛熱通量模擬值與大滿超級站EC和LAS觀測數據的驗證與對比Fig.2 Validation and comparison of simulation values of ParFlow-CLM sensible heat flux and latent heat flux with EC and LAS observation data of Daman station

表2 感熱通量(H)、潛熱通量(LE)、地表溫度(LST)和分層土壤溫度模擬值與觀測數據的驗證及情景比較Table 2 Validation and scenario comparison of simulated and observed data for sensible heat flux,latent heat flux,land surface temperature and soil layer temperature

3.2 地表溫度和分層土壤溫度的驗證與情景對比

由圖3(彩圖見封2)可知，4種情景模擬的地表溫度總體上具有一定的契合度和相關性，R值均高于0.91，RMSE均小于3 K(表2)。其中，Homogeneity情景誤差最大，MIX情景模擬值與觀測數據誤差最小、相關性最高，表明同時設置水平和垂直向異質性并考慮包氣帶—飽水帶三維水力特征情景下地表溫度的模擬精度最高。

圖3 ParFlow-CLM地表溫度模擬值與大滿超級站觀測數據的驗證與對比Fig.3 Validation and comparison of simulation values of ParFlow-CLM land surface temperature with observation data of Daman station

利用大滿超級站2、4、10、20、40、80、120、160(cm)多層土壤溫度觀測數據對4種情景模擬的不同土壤層溫度進行日尺度比較(圖4，彩圖見封2)，可以看出，4種情景模擬的土壤溫度廓線與觀測值整體趨勢一致，尤其是CDHSP和MIX情景下二者溫度廓線形狀最接近，但Homogeneity情景下土壤溫度模擬值的垂向差異較大，尤其是在0.2～0.5 m土層深度范圍內，GLHYMPS情景溫度廓線的峰值出現更早且持續時間更久，與觀測值差異最大。

圖4 大滿超級站觀測數據與ParFlow-CLM土壤層溫度模擬值對比Fig.4 Comparison of simulation values of ParFlow-CLM soil layer temperature with observation data of Daman station

進一步選取代表性較強的第一層、第四層和第七層土壤溫度，與0.02 m、0.2 m和1.2 m土層深度觀測值進行定量比較(圖5)，可以看出，Homogeneity情景模擬值存在低估趨勢，MIX情景模擬值與觀測值最接近；結合表2中分層土壤溫度模擬值與觀測值的R和RMSE，可知MIX情景模擬的分層土壤溫度精度最高，GLHYMPS情景精度最低，原因在于其使用的滲透系數難以反映真實的土壤水力傳導情況。

3.3 空間分布特征

為探究土壤和含水層系統水力特征對地表熱通量和溫度空間分布的影響，以Homogeneity情景模擬結果為基準，進一步繪制GLHYMPS、CDHSP、MIX情景與Homogeneity情景模擬的地下水位埋深(WTD)、感熱通量、潛熱通量和地表溫度空間差值圖(圖6，彩圖見附錄2)。1)從地下水位埋深的空間差值分布看，CDHSP情景模擬結果與Homogeneity情景較接近，MIX情景模擬結果在地形起伏較大的山區與Homogeneity情景差異較大，原因在于山區垂向滲透系數異質性較高，且存在地層分層現象。4種情景在平原區的整體差異較小。2)從地表熱通量和溫度的空間差值分布看，與Homogeneity情景相比，GLHYMPS情景的感熱通量模擬值整體偏高，尤其是東北裸地和西南林草區，而潛熱通量在上述兩個區域基本為負值；CDHSP和GLHYMPS情景模擬結果與Homogeneity情景相反，在西部裸地及稀疏植被覆蓋地區尤為明顯，原因在于二者不同的土壤滲透系數參數化設置方案；相比Homogeneity情景，MIX情景稍顯復雜，感熱通量在東北裸地和西南林草區較高，在西北和西南裸地區偏低，地表溫度和潛熱通量值整體呈現偏高趨勢，僅在西南黑河干、支流沿岸出現低值。綜上，模型模擬的地下水埋深、地表熱通量和地表溫度的空間差值分布具有一定的敏感性，受下墊面、地形等環境因素影響而產生不同的陸面過程響應。整體上，4種情景在中部平坦的農田區差異較小，而在東北海拔較高的荒漠區差異明顯。

圖6 研究區平均地下水位埋深、感熱通量、潛熱通量和地表溫度的空間差值對比Fig.6 Spatial difference comparison of averaged groundwater table depth,sensible heat flux,latent heat flux and land surface temperature

3.4 不同下墊面的影響

為深入探究下墊面因素通過土壤和含水層系統異質性對地表熱通量的影響程度，本文繪制了模擬區域內裸地(荒漠及高山稀疏植被)和農田(農作物及農用地混合像元)兩種下墊面對應的平均感熱通量、潛熱通量和地表溫度箱型圖(圖7)。可以看出，整體上Homogeneity、GLHYMPS和MIX 3種情景結果較接近；CDHSP情景與其他情景模擬結果差異(標準差和中位數)較大，感熱通量模擬值較高，潛熱通量和地表溫度模擬值偏低，這說明CDHSP情景下土壤和含水層系統的水力特性設置使得地表溫度和熱通量的模擬值域更廣泛、分布更離散。橫向比較，對于農田下墊面，4種情景下的平均熱通量和地表溫度模擬結果差異較小，這與前面分析結果一致；而對于裸地下墊面，Homogeneity和GLHYMPS情景的模擬結果中位數相近，MIX情景模擬的感熱通量中位數比前兩種情景高、潛熱通量的中位數則較低。從模擬結果的值域范圍看，MIX情景與Homogeneity情景最為接近。綜上可知，裸地下墊面對地表熱通量和溫度的影響較大，MIX情景模擬結果的穩定性最好。

圖7 裸地和農田下墊面的平均感熱通量、潛熱通量和地表溫度箱型圖Fig.7 Box plots of average sensible heat flux,latent heat flux and land surface temperature of bare land and cropland

4 結論

理解滲透系數的異質性對區域尺度水文集成模擬尤為重要。本文使用基于三維物理過程的地下水—陸面過程耦合模型ParFlow-CLM，依據不同土壤和含水層滲透系數情景設置，定量描述土壤和含水層系統異質性對地表熱通量的影響。經過對地下水位埋深、地表感(潛)熱通量、地表溫度和分層土壤溫度與觀測數據的比較以及典型下墊面影響下時空分布特征的統計分析，主要結論如下：1)在黑河流域中游密集農業區，忽略水平向異質性會導致地表熱通量模擬精度降低和不確定性增加；忽略垂直向異質性會導致地表溫度模擬精度降低。同時，長期、持續的灌溉打破了潛水埋深對潛熱通量的絕對控制，淺層水分成為主導土壤異質性影響地表能量分布及變化的因素之一。2)3種異質性情景與均值情景的地表熱通量和溫度差值空間分布對下墊面、地形等環境因素具有一定敏感性，下墊面因素對地表能量受土壤和含水層系統異質性的影響程度具有控制作用，尤其是對感熱通量和地表溫度的模擬影響較大。3)滲透系數的空間異質性在地形起伏較大的山區對地表熱通量和溫度的模擬影響遠大于平原區。

本研究存在以下不足：1)模擬結果受輸入數據(滲透系數數據集、驅動數據等)質量限制，需利用長期監測數據，進行不同尺度的對比和驗證；同時，滲透系數數據集需利用野外勘察及輔助信息進行校核。2)模擬時長仍需繼續擴展，最好能涵蓋一個旱澇年周期，可深入探究滲透系數對地表熱通量的年際及年內影響規律。3)黑河流域中游模擬區域滲透系數的異質性變化范圍有限，仍需擴展到整個黑河流域甚至西北地區，以更全面地模擬土壤和含水層系統異質性對地下水—地表水交互的影響。

國家青藏高原科學數據中心提供黑河流域大氣驅動數據集(2000-2018)、土地利用覆被數據集、30弧秒分辨率地表水和地下水灌溉量數據集、流域中游渠道流量測量數據集以及張掖灌溉渠系數據集、中國土壤數據集和祁連山綜合觀測網“黑河流域地表過程綜合觀測網數據集”，CGIAR-CSI提供SRTM DEM數據，此致謝忱！