李 飛， 郭佳鍇， 張世強

（1.西北大學(xué)陜西省地表系統(tǒng)與環(huán)境承載力重點實驗室，陜西西安 710127； 2.西北大學(xué)城市與環(huán)境學(xué)院，陜西西安 710127）

0 引言

據(jù)統(tǒng)計，多年凍土約占北半球陸地面積的24%［1］，我國是世界第三凍土大國，多年凍土主要分布在青藏高原、東北北部山區(qū)和西北高山區(qū)。凍土作為一類特殊的巖土，與外界能量的交換會改變自身的熱力學(xué)和水力學(xué)性質(zhì)，因此在凍融階段的物理過程極其復(fù)雜。活動層內(nèi)部土壤的水熱遷移在地氣交互和能水平衡中起著重要作用［2］，影響冰凍圈生態(tài)環(huán)境［3］和流域水文過程［4］。

青藏高原多年凍土總面積約1.06×106km2，約占高原總面積的40%［5］。高海拔和復(fù)雜地形使青藏高原對全球氣候變化極為敏感，并通過能水交換與大氣循環(huán)作用反饋到全球氣候變化進程中［6］，是全球氣候的放大鏡和調(diào)節(jié)器［7］，而作為高原典型下墊面類型的多年凍土是影響全球氣候變化的關(guān)鍵因子。近幾十年，在全球氣候變暖的背景下，青藏高原的多年凍土呈現(xiàn)退化的趨勢，主要表現(xiàn)為多年凍土年平均地溫升高，活動層厚度加深，凍土下限海拔上升，多年凍土區(qū)面積減?。?-11］。已有研究表明，凍土退化使高原植被退化，加速了高原草地荒漠化和鹽漬化的速度［3，12］，同時流域內(nèi)夏季徑流減少、冬季徑流增加、年內(nèi)徑流過程變緩，變化強度與流域多年凍土覆蓋率有關(guān)［13-16］，部分地區(qū)地下水儲量增加［17］等。研究凍土區(qū)土壤的水熱變化可以為未來氣候模式下更準(zhǔn)確預(yù)估凍土的變化趨勢，為當(dāng)?shù)厣鷳B(tài)環(huán)境保護、水資源利用、下游社會經(jīng)濟發(fā)展提供重要決策依據(jù)。

在寒區(qū)陸面過程研究工作中，土壤內(nèi)部水熱過程的準(zhǔn)確描述一直是難點之一。國內(nèi)外科學(xué)家基于不同的算法與需求開發(fā)大量陸面過程模型，其中代表性的有SHAW（Simultaneous Heat and Water Model）［18］、CoupModel（Coupled Heat and Mass Transfer Model for Soil-Plant-Atmosphere System）［19］和VIC-3L（Variable Infiltration Capacity-3 Layer Model）［20］等，這些模型通過耦合凍土模塊，實現(xiàn)了凍土水熱過程的模擬。

陸面模型為了便于描述陸面過程，通常將復(fù)雜的陸面過程參數(shù)化，在影響凍土水熱模擬的參數(shù)中，土壤導(dǎo)熱率影響土壤和大氣的能量交換以及土壤內(nèi)部間的能量傳遞，對土壤溫度模擬的影響很大［21］。很早就有國外學(xué)者對土壤熱導(dǎo)率進行了相關(guān)研究［22-23］，Johansen［24］在Kersten［22］研究的基礎(chǔ)上于1975年首次提出了歸一化導(dǎo)熱系數(shù)的概念，建立了一個模擬非飽和土壤熱導(dǎo)率的半理論模型，適用性較廣。C?té 等［25］、陸森等［26］在Johansen［24］算法的基礎(chǔ)上對土壤熱導(dǎo)率進一步改進提高了模擬精度，擴大了適用范圍。Yang 等［27］結(jié)合全球土壤實測資料改進了Johansen［24］算法，提高了土壤在凍結(jié)狀態(tài)下的導(dǎo)熱 系數(shù)精 度。 Bao 等［28］在EBM（Eaphy-Based Model）中借鑒了Yang 等［27］的成果，并嵌入了冰的導(dǎo)熱算法［29］，在祁連山多年凍土區(qū)水熱模擬方面取得了不錯的進展。在國內(nèi)凍土導(dǎo)熱系數(shù)的研究中，李韌等［30］基于野外實測數(shù)據(jù)推導(dǎo)出適用于青藏高原北部的計算公式。而羅斯瓊等［31］、王愚等［32］在Farouki［33］、Johansen［24］、C?té等［25］參數(shù)化方案的基礎(chǔ)上發(fā)展了適用于青藏高原中部的導(dǎo)熱系數(shù)方案并取得了較好的效果。劉為民等［34］基于土體成分和特性等因素，任小鵬等［35］考慮了土壤水分的相變和巖土的特點建立了各自的導(dǎo)熱計算模型，但都存在局限性。未凍水含量對土壤熱力學(xué)參數(shù)和物理性質(zhì)有著重要意義［36-37］，影響凍融過程中水分的遷移［38］。未凍水的計算模型包括經(jīng)驗?zāi)Ｐ秃屠碚撃Ｐ?，部分模型如GIPL2（Geophysical Institute Permafrost Lab 2.0）［39］采用徐斅祖等［37］、Lovell［40］基于實驗數(shù)據(jù)提出的關(guān)于土壤溫度的未凍水經(jīng)驗公式。VIC-3L、CLM（Community Land Model）［41］、FEFLOW（Finite Element Subsur-face Flow System）［42］等模型采用了在不同的凍融方案下的未凍水計算公式［43-45］。張世強等［46-47］、Pan 等［48］和郭林茂等［49］分別利用不同的陸面模型在青藏高原展開水熱模擬研究，對各模型未凍水方案在水熱模擬方面的適用性有著很好的參考作用。

已有研究者在青藏高原多年凍土區(qū)探究陸面模型中參數(shù)化方案對水熱模擬的影響，Hu等［50］在唐古拉站和西大灘站對比了Johansen［24］、Farouki［33］和羅斯瓊等［31］三種熱導(dǎo)率方案，發(fā)現(xiàn)羅斯瓊等［31］的方案在未凍結(jié)土壤中表現(xiàn)最好，而Johansen［24］方案在凍結(jié)土壤中表現(xiàn)最好。Yang 等［51］在CLM 5.0［52］的基礎(chǔ)上，在6 種典型的植被下墊面下比較了9 種標(biāo)準(zhǔn)化導(dǎo)熱率方案對土壤熱狀態(tài)的模擬結(jié)果，結(jié)果表明選取的所有導(dǎo)熱率方案的計算值都偏高并低估了凍結(jié)期的土壤溫度，在這些方案中Johansen［24］方案及其派生方案在導(dǎo)熱率的數(shù)值估算以及模擬土壤溫度方面明顯更優(yōu)，其中Balland等［53］的導(dǎo)熱率方案在不同類型土壤的綜合適用性最好。劉火霖等［54］在那曲站評估了Noah-MP（The Community Noah Land Surface Model with Multiparameterization Options）陸面模型［55］下Niu 等［44］和Koren 等［56］兩種未凍水參數(shù)化方案，結(jié)果表明未凍水方案對土壤濕度模擬的影響明顯大于土壤溫度，未凍水方案的差異顯著影響凍融期土壤水分含量的變化。Hu等［57］基于土壤溫度、土壤顆粒比表面積、水的類型和土壤水分曲線建立的未凍水計算公式，歸納并總結(jié)了共34種未凍水方案，并根據(jù)實地數(shù)據(jù)進行了評估，結(jié)果表明經(jīng)驗公式可用于計算未凍水含量，而物理模型具有更高的精度，但更為復(fù)雜，難以在實際情況中應(yīng)用。在這些方案中，van Genuchten［58］、Kozlowski［59］、Zhang等［60］的方案表現(xiàn)較好。

本研究在VIC-CAS 模型的基礎(chǔ)上，綜合比較了眾多研究者在青藏高原水熱模擬的結(jié)果，分別用EBM 導(dǎo)熱率算法和CLM 5.0 未凍水算法替換原模型相應(yīng)算法，利用長江源區(qū)的沱沱河站點觀測資料進行為期一年的單點模擬試驗，分析比較了不同算法下土壤凍融過程中溫度和濕度的模擬效果，并針對模擬結(jié)果中存在的問題進行分析，為后續(xù)的模型改進提供支撐。

1 研究區(qū)概況與數(shù)據(jù)來源

沱沱河觀測站（92°37′ E，33°57′ N）位于長江源的多年凍土連續(xù)區(qū)［61］，海拔4 533.1 m，地勢平緩，氣候干燥寒冷。地屬高寒草甸地帶，植被稀疏，植被覆蓋度低。根據(jù)1997—1998年的觀測，沱沱河站點的年平均氣溫達-3.7 ℃，最低氣溫為-32.9 ℃，最高氣溫為24.3 ℃。流域內(nèi)降水較少，年降水量在300 mm 左右，降水集中在6—9 月，約占年降水量的85%，并在9月份降水量達到峰值。

本文觀測數(shù)據(jù)來自“全球能水平衡試驗-青藏高原亞洲季風(fēng)試驗”（Game-Tibet）項目1994—1999年的觀測期內(nèi)長江源區(qū)沱沱河站點的觀測資料。觀測數(shù)據(jù)來自土壤溫度和濕度觀測系統(tǒng)（SMTMS）收集的土壤溫、濕度數(shù)據(jù)。由于受限于當(dāng)?shù)氐臍夂驐l件和儀器損耗等客觀因素，觀測數(shù)據(jù)時間序列不連貫且部分時段數(shù)據(jù)誤差較大，因此選擇1997 年7月1 日至1998 年9 月29 日超過一年的觀測數(shù)據(jù)，時間分辨率為1 h。

2 實驗方法

2.1 VIC-CAS模型

VIC-CAS 由趙求東等［62］在VIC-3L 的基礎(chǔ)上耦合冰川模塊構(gòu)成。VIC-CAS 模型以能量平衡和水量平衡為基礎(chǔ)，在垂直方向上形成大氣-植被-土壤的陸面模式，在單個格網(wǎng)上考慮植被、降雨、土壤、積雪等綜合因素，輸出單個格網(wǎng)的各模擬分量，并可通過匯流模型對流域內(nèi)產(chǎn)流進行匯總。

模型的輸入數(shù)據(jù)主要包括氣象驅(qū)動數(shù)據(jù)、土壤參數(shù)數(shù)據(jù)、植被參數(shù)庫和植被覆蓋數(shù)據(jù)。模型的氣象驅(qū)動數(shù)據(jù)包括氣象站觀測的逐時氣溫、風(fēng)速、氣壓等和日降水資料作為模型的氣象驅(qū)動數(shù)據(jù)輸入。土壤參數(shù)數(shù)據(jù)來自國際地圈生物圈計劃（IGBP）的全球土地覆蓋特征數(shù)據(jù)集（GLCC）并參考了美國農(nóng)業(yè)部土壤質(zhì)地分類依據(jù)（USDA）和Cosby 等［63］的成果；植被參數(shù)庫來自于GLCC 的亞洲青藏高原植被庫，部分參數(shù)參考了高原其他地區(qū)的研究資料［64-65］；植被分布來源于美國馬里蘭大學(xué)發(fā)布的全球陸面覆蓋類型數(shù)據(jù)集（UMD）。

2.2 導(dǎo)熱率算法

影響土壤導(dǎo)熱率數(shù)值的因素有很多，如土壤孔隙度、粒徑、飽和度、溫度［44，66-68］等。在對土壤導(dǎo)熱率的眾多算法中，VIC-CAS 模型中導(dǎo)熱率選擇應(yīng)用最為廣泛的Johansen［24］算法，導(dǎo)熱率K可表示為：

式中：Kdry為干土導(dǎo)熱率，只與土壤塊體密度ρd有關(guān)，可表示為：

Ksat為飽和土壤導(dǎo)熱率，主要由土壤基質(zhì)、水、冰的導(dǎo)熱率及孔隙度決定，可表示為：

式中：T和Tf分別為土壤溫度和冰點溫度；n為土壤孔隙度；Ks為土壤基質(zhì)導(dǎo)熱率；Kliq為水的導(dǎo)熱率，取值0.6 W·m-1·K-1；θliq為土壤液態(tài)水含量；Kice為冰的導(dǎo)熱率，取值2.29 W·m-1·K-1。

Ke是Kersten 系數(shù)，是關(guān)于土壤飽和度Sr的函數(shù)，可表示為：

式中：Sr與土壤中液態(tài)水含量θliq、含冰量θice以及孔隙度n有關(guān)，可表示為：

本研究選擇EBM導(dǎo)熱率算法替換原模型算法，EBM中導(dǎo)熱率K表示為：

式中：Kdry為干土導(dǎo)熱率，計算公式與式（2）相同。

Kwet［27］為濕土導(dǎo)熱率，可表示為：

式中：w為土壤有效濕度，與土壤中液態(tài)水含量θliq和土壤飽和含水量θsat有關(guān)，可表示為：

Ki為土壤中冰的總導(dǎo)熱率，與土壤含冰量有關(guān)，可表示為：

2.3 未凍水算法

VIC-CAS模型的最大未凍水含量公式［43］如下：

式中：θsat為土壤飽和含水量；g為重力加速度，取值9.81 m·s-2；φc為進氣值；T為土壤溫度；Lf為相變潛熱，取值3.36×105J·kg-1。Bp為土壤孔隙形狀參數(shù)，是關(guān)于水力傳導(dǎo)層指數(shù)函數(shù)，可表示為：

本研究選擇CLM 5.0 中的未凍水算法代替原模型算法，公式如下：

式中：θsat由土壤有機質(zhì)和礦物質(zhì)含量以及各自的孔隙度決定，可表示為：

式中：fo為土壤有機質(zhì)含量；θo為有機質(zhì)孔隙度，取值0.9［63］；θm為礦物土孔隙度，由含砂量Psand決定，可表示為：

Ψsat為飽和土壤基質(zhì)勢，由土壤有機質(zhì)和礦物質(zhì)含量以及各自的基質(zhì)勢決定，可表示為：

式中：Ψo為飽和有機質(zhì)基質(zhì)勢，取值-10.3 mm［69］；Ψm為飽和礦物土基質(zhì)勢，由土壤中砂土含量Psand決定，可表示為：

B為Clapp-Hornberger 經(jīng)驗系數(shù)［63，70］，可 表示為：

式中：Bo取2.7［69］；Bm由土壤中黏土含量Pclay決定，可表示為：

2.4 模型參數(shù)率定

本研究中VIC-CAS 模型模擬的初始態(tài)文件根據(jù)站點的觀測數(shù)據(jù)得到。模型的率定期選擇1997年7月31日—1997年9月31日，驗證期選擇1997年10 月1 日—1998 年9 月29 日。率定后的主要參數(shù)如表1所示。

表1 VIC-CAS模型長江源區(qū)沱沱河站點率定后主要參數(shù)Table 1 Main parameters of VIC-CAS model after calibration at Tuotuohe station in the source area of Yangtze River

2.5 模型評估指標(biāo)

利用率定后的模型進行土壤水熱數(shù)值模擬，將1997年10月1日—1998 年9 月29 日作為模型的驗證期分析模型的模擬效果和可靠性。根據(jù)丁永建等［71］對凍土區(qū)土壤水分含量變化的描述將模擬期土壤的溫度和濕度變化分為四個階段，分別是秋季降溫階段，1997 年9 月21 日到10 月21 日；冬季凍結(jié)階段，1997年10月21日到1998年3月21日；春季升溫階段，1998 年3 月21 日到4 月21 日；夏季融化階段，1997年4月21日到9月21日。

為了更直觀地比較土壤溫、濕度的模擬值和觀測值，將土壤溫度模擬值線性插值到觀測值所在的深度。同時在分析土壤濕度時，將每層包含的各個深度土壤濕度觀測的平均值與模擬值進行比較。

本文采用納什效率系數(shù)（NSE）和均方根誤差（RMSE）來作為評價指標(biāo)。

其中，NSE 是評價水文模型模擬結(jié)果好壞的常用指標(biāo)，計算公式如下：

NSE 的取值范圍從負(fù)無窮到1。NSE 越接近1，模擬結(jié)果越好；接近0 表示模擬結(jié)果接近觀測值的平均值水平，即總體結(jié)果可信，但模擬過程誤差大；遠小于0，則模型不可信。

RMSE 作為評估預(yù)測模型精度的指標(biāo)，通常用來衡量真實值與模擬值之間的偏差，公式如下：

相關(guān)系數(shù)的取值范圍從0 到正無窮。RMSE 越接近0，真實值與模擬值的偏差越小，模擬精度越高。

3 結(jié)果分析

3.1 導(dǎo)熱率算法結(jié)果對比

圖1為兩種導(dǎo)熱率算法下不同深度模擬的土壤逐時體積熱通量對比。EBM 算法的熱通量模擬值在7 個深度均低于原算法，熱通量的減少表明EBM算法相較原算法阻礙了土壤凍融過程，表現(xiàn)為土壤溫度在降溫凍結(jié)階段會偏高，在升溫融化階段會偏低。兩者在20～60 cm 深度凍結(jié)期的熱通量差值低于融化期，表明對凍結(jié)期土壤溫度變化的影響小于融化期溫度變化，而在60 cm 深度以下，凍結(jié)期和融化期的差值接近，對土壤溫度變化影響差別較小。

圖1 沱沱河站VIC-CAS算法和EBM算法模擬的7個深度土壤逐時熱通量對比Fig.1 Comparisons in hourly soil heat fluxes between simulated by VIC-CAS algorithm and EBM algorithm at seven depths in Tuotuohe station

圖2為兩種導(dǎo)熱率算法下不同深度模擬與觀測的土壤日均溫度對比。在80 cm 深度以上，VICCAS 算法溫度模擬值在夏季融化期前40 天左右接近觀測值，其余時段略低于觀測值。在80 cm 及80 cm 深度以下，從秋季降溫至冬季土壤溫度下降至0 ℃，模擬值接近觀測值；冬季土壤開始凍結(jié)至凍結(jié)期末，模擬值高于觀測值；春季升溫至夏季融化期末，模擬值低于觀測值。7 個深度模擬結(jié)果總體較好，20～120 cm 深度模擬的土壤溫度的NSE 大于0.9，130～160 cm 深度模擬的土壤溫度的NSE 大于0.8。

圖2 沱沱河站觀測與VIC-CAS算法和EBM算法模擬的7個深度土壤日均溫度對比Fig.2 Comparisons in daily soil average temperature between observed and simulated by VIC-CAS algorithm and EBM algorithm at seven depths in Tuotuohe station

EBM 算法模擬值在20 cm 和60 cm 深度整體接近觀測值，在40 cm 深度，夏季融化期前三個月左右高于觀測值，其余時段略低于觀測值。在60 cm 深度以下，從秋季降溫至春季升溫，模擬值高于觀測值；從春季升溫至夏季融化期末，模擬值低于觀測值。20～80 cm，NSE 大于0.9，模擬效果較好；在100～160 cm，NSE 低于0.8，且隨深度增加，NSE 下降，模擬效果較差。

兩種算法模擬的土壤溫度趨勢接近，從差值（EBM-VIC-CAS）上來看，在80 cm 深度以上，兩種算法的差值在冬季凍結(jié)期和夏季融化期均為正值，原因是EBM算法對導(dǎo)熱率的計算低于原算法，減緩了土壤的凍結(jié)過程，使冬季模擬的土壤溫度偏高，這種影響持續(xù)到了夏季，盡管導(dǎo)熱率的降低減緩了融化過程，但夏季模擬的溫度仍高于原算法。冬季凍結(jié)期前兩個月和夏季融化期前兩個月差值較大，且冬季差值高于夏季，主要原因是土壤水分相變的初始階段能量變化較大，對土壤的溫度差異影響明顯。同時在春末夏初土壤溫度升至0值前后小段時間，差值表現(xiàn)為負(fù)值，在這段時間內(nèi)土壤水分相變吸收大量熱量，短時間內(nèi)使土壤迅速降溫。 在80 cm 深度及以下，隨深度增加，冬季模擬的土壤溫度對夏季土壤溫度的影響在減小，夏季融化期兩者差值表現(xiàn)為負(fù)值。同時由于深層土壤能量傳遞的滯后性，冬季和夏季差值的峰值時間在推遲，發(fā)生在凍結(jié)期和融化期的中段。

在80 cm 深度以上，EBM 算法的NSE 大于VICCAS 方案，同時RMSE 小于后者，EBM 算法的效果較優(yōu)。在80 cm 深度，EBM 算法的模擬效果下降，NSE 小于VIC-CAS 算法，但仍大于0.9。在80 cm深度以下，兩種算法的NSE 隨深度增加而下降，EBM 算法在100 cm 深度最高為0.7858，小于0.8；而VIC-CAS 算法在3 個深度的NSE 均超過0.8。同時，3 個深度EBM 算法的RMSE 也大于后者，EBM算法在深層土壤模擬效果較差。

圖3為兩種導(dǎo)熱率算法下三層土壤模擬與觀測各層的日平均土壤濕度對比。表層（0～0.1 m）和中層土壤（0.1～0.6 m）兩種算法的模擬值接近，而在下層土壤（0.6～1.6 m）在冬季凍結(jié)期有明顯差異。導(dǎo)熱率差異引起的土壤水分含量變化在表層不明顯。兩種算法在中層土壤的模擬值在凍結(jié)期整體接近觀測值，模擬效果較好，而在其余時段低于觀測值。兩者的差值在冬季凍結(jié)期開始一個月內(nèi)較明顯，主要原因可能是導(dǎo)熱率差異在凍結(jié)初期對土壤溫度影響較明顯，減緩了土壤水分的相變過程，使土壤中滯留的液態(tài)水含量較高。而在下層土壤，兩種算法模擬結(jié)果皆高于觀測值，VIC-CAS 算法模擬值在凍結(jié)期后半段至夏初土壤開始融化這段時間低于EBM 算法，更接近觀測值。由于土壤能量傳遞的滯后性因素，下層土壤水分變化影響要延遲一段時間顯現(xiàn)出來，土壤在冬季增加的水分是由于兩種算法在凍結(jié)期土壤溫度差值高于中層土壤，表明對土壤的凍結(jié)過程影響更大。在這種情況下，土壤水分的相變速率遠小于中層土壤，滯留的液態(tài)水含量更高。

圖3 沱沱河站觀測與VIC-CAS算法和EBM算法模擬的三層土壤日均濕度對比Fig.3 Comparisons in daily soil average moisture between observed and simulated by VIC-CAS algorithm and EBM algorithm at three layers in Tuotuohe station

圖4為三層土壤兩種導(dǎo)熱率算法下模擬的日均含冰量對比。兩種算法在表層土壤的含冰量方面差別不大。EBM 算法在冬季前兩個月的中層土壤含冰量明顯低于原算法，原因是EBM 算法模擬的土壤溫度在凍結(jié)初期略高于后者，表現(xiàn)出土壤水分相變過程減緩，使土壤含冰量計算偏低。EBM 算法在下層土壤的含冰量低于原算法，時間正是下層土壤水分含量變化的時段，冬季含冰量降低會減小凍結(jié)土壤的整體導(dǎo)熱率，進一步減緩凍結(jié)過程。

圖4 沱沱河站VIC-CAS算法和EBM算法模擬的三層土壤日均含冰量對比Fig.4 Comparisons in daily soil average ice content between observed and simulated by VIC-CAS algorithm and EBM algorithm at three layers in Tuotuohe station

圖5為兩種導(dǎo)熱率算法下模擬與觀測的逐時融化深度對比。觀測的最大融化深度5 m 左右，VICCAS 算法的模擬最大融化深度8 m 左右，約高于觀測值3 m；而EBM 算法模擬最大融化深度2.2 m 左右，約低于觀測值2.8 m。EBM 算法模擬的融化深度變化趨勢大于VIC-CAS 算法，且在數(shù)值上更接近觀測值。VIC-CAS 算法對導(dǎo)熱率的計算值較高，使融化期土壤的暖舌較深，EBM 算法在融化期對土壤的水熱變化速率弱于VIC-CAS 算法，暖舌較淺，但更接近觀測值。

圖5 沱河站觀測與VIC-CAS算法和EBM算法模擬的土壤逐時融化深度對比Fig.5 Comparisons in hourly soil melting depth between observed and simulated by VIC-CAS algorithm and EBM algorithm in Tuotuohe station

3.2 未凍水算法結(jié)果對比

圖6為兩種未凍水算法下模擬的不同深度土壤逐時土壤熱通量對比?？傮w上7個深度兩種算法的體積熱通量差異并不明顯，模擬趨勢相近，結(jié)果表明未凍水含量差異對土壤能量傳遞數(shù)值的影響很小，對土壤溫度的變化影響有限。

圖6 沱沱河站VIC-CAS算法和CLM 5.0算法模擬的7個深度土壤逐時熱通量對比Fig.6 Comparisons in hourly soil heat fluxes between simulated by VIC-CAS algorithm and CLM 5.0 algorithm at seven depths in Tuotuohe station

圖7為兩種未凍水算法下模擬與觀測的不同深度土壤日均溫度對比。7 個深度下兩種算法模擬結(jié)果整體趨勢接近，差異總體較小。兩種算法的NSE和RMSE 差異較小，CLM 5.0 未凍水算法的NSE 和RMSE 在20 cm 和40 cm 深度稍優(yōu)于VIC-CAS 算法，NSE 在0.9 以上。在60 cm 深度以下，稍差于原算法，但NSE仍大于0.8。

圖7 沱沱河站觀測與VIC-CAS算法和CLM 5.0算法模擬的7個深度土壤日均溫度對比Fig.7 Comparisons in daily soil average wtemperature between observed and simulated by VIC-CAS algorithm and CLM 5.0 algorithm at seven depths in Tuotuohe station

圖8為兩種未凍水算法下三層土壤模擬與觀測的各層平均土壤濕度對比。CLM 5.0 未凍水算法模擬的三層土壤的濕度均低于VIC-CAS 算法。CLM 5.0算法模擬的表層和中層土壤濕度的RMSE和NSE低于VIC-CAS算法，而下層土壤模擬效果優(yōu)于后者。從差值（CLM 5.0-VIC-CAS）上來看，在表層和中層土壤，差異主要集中在冬季凍結(jié)期，原因是未凍水算法影響冬季土壤水分含量，且計算值較低。下層土壤差值低于前兩層，主要原因是前兩層土壤中的含冰量增加減少了下層土壤水分在凍結(jié)期的遷移，相比于原算法在下層土壤的模擬濕度偏高，在冬季更接近觀測值。

圖8 沱沱河站觀測與VIC-CAS算法和CLM 5.0算法模擬的三層土壤日均濕度對比Fig.8 Comparisons in daily soil average moisture between observed and simulated by VIC-CAS algorithm and CLM 5.0 algorithm at three layers in Tuotuohe station

圖9為兩種未凍水算法模擬的三層土壤日均含冰量對比。三層土壤模擬的含冰量變化趨勢接近，且表現(xiàn)為隨深度增加，土壤中冰的出現(xiàn)時間延遲。CLM 5.0 算法在表層和中層土壤模擬的含冰量在凍結(jié)期高于VIC-CAS 算法。而在下層土壤的含冰量相差不大。土壤中冰含量差異主要受冬季土壤中未凍水含量的影響，模擬的土壤中含冰量變化引起的相變能量傳遞差異也是間接因素。

圖9 沱沱河站VIC-CAS算法和CLM 5.0算法模擬的三層土壤日均含冰量對比Fig.9 Comparisons in daily soil average ice content between observed and simulated by VIC-CAS algorithm and CLM 5.0 algorithm at three layers in Tuotuohe station

圖10 為兩種未凍水算法模擬與觀測的逐時融化深度對比。兩種算法融化深度模擬的趨勢接近，CLM 5.0 算法融化深度略低于VIC-CAS 算法。兩種算法融化深度在融化期中期高于觀測值，在末期低于觀測值，數(shù)值變化較大。上述結(jié)果表明，CLM 5.0 算法相較原算法對土壤中水熱過程在土壤層中的能量傳遞速率和影響深度上變化不明顯。

圖10 沱沱河站觀測與VIC-CAS算法和CLM 5.0算法模擬的土壤逐時融化深度對比Fig.10 Comparisons in hourly soil melting depth between observed and simulated by VIC-CAS algorithm and CLM 5.0 algorithm in Tuotuohe station

4 討論

4.1 模型對土壤溫度的模擬差異

多位研究者采用不同的陸面模型模擬了青藏高原多年凍土區(qū)的土壤溫度，在唐古拉站高寒草甸區(qū)，劉洋等［72］基于SHAW 模型2007 年的結(jié)果顯示，50～150 cm 深度土壤溫度誤差較大，而0～50 cm 和150 cm 深度以下模擬值和觀測值擬合較好。胡國杰等［73］基于CoupModel模型在2008年的結(jié)果顯示，0～105 cm 深度土壤模擬值和觀測值吻合得較好，在105 cm 以下深度土壤的誤差增大。馬啟民等［74］基于一維水熱耦合模型2008 年的結(jié)果表明，0～50 cm深度土壤溫度模擬效果較好，但90 cm 以下偏差增大。秦慧艷等［75］基于GIPL2 模型在可可西里站高寒草甸區(qū)2011—2012年的結(jié)果顯示，0～100 cm模擬值與觀測值相關(guān)性較高，120 cm 深度以下兩者偏差增大，且模擬值稍高于觀測值。郭林茂等［49］基于FEFLOW 模型在風(fēng)火山高寒草甸區(qū)2006—2008 年的土壤溫度模擬結(jié)果顯示，0～20 cm 深度土壤的模擬精度低于20 cm深度以下土壤。

本文替換EBM 導(dǎo)熱率算法后，在20～60 cm 深度，模擬值與觀測值間的誤差減?。辉?0～160 cm 深度，模擬值與觀測值的誤差增大。替換CLM 5.0 未凍水算法后與原算法的結(jié)果差異較小。與上述結(jié)果對比可以看出，GIPL2、SHAW 等模型與本文的模擬結(jié)果較為接近，即整體趨勢模擬上相關(guān)性較一致，淺層土壤的精度高于深層土壤，且深層土壤的精度隨深度增加而較差，其中，劉洋等［72］和馬啟民等［74］的結(jié)果在深度上更接近本文。FEFLOW 模型的結(jié)果與本文結(jié)果相反，在深層土壤溫度模擬精度優(yōu)于淺層。

造成這種結(jié)果可能的原因是：（1）各層輸入的土壤水熱參數(shù)存在誤差。多年凍土區(qū)土壤存在異質(zhì)性，不同深度土壤特性差異較大。本文模擬時將土壤柱劃分為10 cm、50 cm、100 cm 三層土壤；秦艷慧等［75］根據(jù)文獻［76］按5 cm、20 cm、大于100 cm 劃分土壤層；胡國杰等［73］、馬啟民等［74］按不同深度劃分土壤厚度，在淺層土壤劃分厚度較淺，而深層土壤劃分較厚，輸入的土壤參數(shù)不能代表各層的真實情況，這導(dǎo)致了淺層土壤的模擬精度高于深層土壤。（2）模型中凍土算法熱量平衡方程大都來自于Stefan 的熱傳導(dǎo)方程，在計算中各層土壤參數(shù)除因深度劃分不均勻引起的誤差外，自身與實際值也存在偏差，熱傳導(dǎo)方程再計算土壤中熱量由上而下或由下而上傳遞過程時，這種誤差會不斷累積。（3）研究區(qū)植被為高寒草甸，根系在0～30 cm 左右，影響淺層土壤的熱性質(zhì)。相比較其他模型，郭林茂等［49］認(rèn)為FEFLOW 模型缺少植被輸入是影響淺層土壤的模擬結(jié)果的重要原因。

4.2 模型對土壤濕度的模擬差異

對青藏高原多年凍土濕度的模擬結(jié)果中，劉洋等［72］的結(jié)果表明SHAW 模型在5～20 cm 土壤的模擬值低于觀測值，而70～210 cm 土壤凍融期模擬值高于觀測值。胡國杰等［73］的結(jié)果表明CoupModel模型在5～20 cm 的模擬值低于觀測值，而在70 cm以下接近觀測值。馬啟民等［74］的結(jié)果顯示0～50 cm的模擬值在凍結(jié)期接近觀測值，而50 cm 以下的模擬值在凍結(jié)期高于觀測值。

本文替換未凍水算法對土壤濕度的差異主要表現(xiàn)在土壤凍結(jié)期到融化期前，融化期濕度變化較小，與劉火霖等［54］未凍水參數(shù)化方案替換的結(jié)果基本一致。對比他人的模擬結(jié)果可以看出，替換EBM導(dǎo)熱率算法后，0～60 cm 土壤濕度在凍結(jié)期結(jié)果上優(yōu)于上述三種模型；60～160 cm 凍結(jié)期的模擬值高于觀測值，與劉洋等［72］、馬啟民等［74］的結(jié)果接近。替換CLM 5.0 未凍水算法后，0～60 cm 土壤模擬值低于觀測值，與劉洋等［72］、胡國杰等［73］、馬啟民等［74］的結(jié)果相同；60～160 cm 土壤模擬值接近觀測值，與胡國杰等［73］的結(jié)果接近。

出現(xiàn)上述差異的原因可能是：（1）VIC-CAS 模型輸出的是三層土壤的平均濕度，反映的是三層土壤中各層整體的水分變化情況，不能表現(xiàn)出土壤水分隨深度的變化及各深度間水分的差異性。而劉洋等［72］、胡國杰等［73］、馬啟民等［74］利用的模型能夠模擬各深度的土壤水分，與觀測值的比較在結(jié)果上可能更好。（2）在替換兩種算法前，VIC-CAS 模型默認(rèn)的土壤濕度模擬值在凍結(jié)期接近觀測值，表現(xiàn)了VIC-CAS 模型在凍結(jié)期能夠較好模擬土壤水分的優(yōu)點。替換EBM 導(dǎo)熱率算法主要對60～160 cm 土壤有影響，使水分模擬偏高，表現(xiàn)出SHAW 等模型在深層土壤模擬較差的缺點，原因可能是替換的導(dǎo)熱率方案計算值相比原方案與實際值偏差增大，并在深層表現(xiàn)得更為明顯。替換CLM 5.0 未凍水算法，降低了0～60 cm 土壤濕度的模擬精度，但提高了60～160 cm 的模擬精度，與胡國杰等［73］的結(jié)果接近?？赡茉蚴荂LM 5.0 未凍水算法水分模擬自身的缺點，熊健勝等［77］認(rèn)為CLM 模型模擬多年凍土?xí)r，凍結(jié)期含水量被低估，而含冰量被高估，但偏低的模擬值彌補了VIC-CAS 模型在深層土壤模擬值較高的缺點，使下層土壤的模擬結(jié)果反而變好。胡國杰等［73］利用的CoupModel 模型的土壤水分模擬值整體偏低，且深層的土壤模擬精度高于淺層土壤，因此與本文的結(jié)果接近。

4.3 模型中冰含量對土壤熱特性的影響

冰作為不透水面，顯著減低土壤水分的遷移，同時冰的導(dǎo)熱率是水的4 倍，影響凍結(jié)土壤的總體導(dǎo)熱性。本文在替換未凍水算法后，模擬結(jié)果圖9中0～60 cm 土壤含冰量有所增加，但圖6 和圖7 中，20～160 cm 各深度替換未凍水方案前后模擬的熱通量和溫度的差異較小。

李振坤［78］在改則地區(qū)的研究表明CLM 3.0 參數(shù)化方案改進前后0～80 cm 土壤溫度模擬值間的相關(guān)系數(shù)大于0.9，含冰量變化對淺層土壤溫度變化較小。本文的結(jié)果對比上述結(jié)果的差異更小，兩種未凍水方案下觀測值與模擬值間相關(guān)系數(shù)差值最大在0.03 左右，由CLM 未凍水方案的缺點來看，實際的含冰量比模擬值更低，含冰量對土壤導(dǎo)熱性能的影響可能更小。

5 結(jié)論

本研究分別用EBM 導(dǎo)熱率算法和CLM 5.0 的未凍水算法替換以VIC-CAS 模型為基礎(chǔ)的導(dǎo)熱率和未凍水算法，對長江源區(qū)沱沱河站進行了為期一年的單點數(shù)值模擬試驗，主要結(jié)論如下：

（1）EBM 算法導(dǎo)熱率低于原算法，會減緩?fù)寥赖膬鋈谶^程，對不同深度土壤溫度影響較明顯，在土壤濕度模擬上，對深層土壤計算數(shù)值較高，而淺層變化不大。算法提升了淺層土壤溫度的模擬效果，但使深層土壤溫度和濕度模擬效果變差。

（2）CLM 5.0 算法未凍水?dāng)?shù)值低于原算法，主要減少了凍結(jié)期土壤的水分，但土壤水分的變化對土壤能量傳遞影響有限，與原算法土壤溫度較接近。算法對深層土壤濕度的模擬效果變好，但淺層土壤濕度模擬效果變差。

（3）本文雖然替換了兩種土壤參數(shù)算法，但在土壤濕度的模擬效果上并沒有顯著的改進。通過對比可以看出兩種算法與原算法各有優(yōu)缺點，EBM導(dǎo)熱率算法需要改進誤差產(chǎn)生的累積效應(yīng)對深層土壤的影響，而CLM 5.0 未凍水算法在對淺層土壤的濕度模擬上過低仍是個需要進一步考慮的問題。