李宏亮， 王璞玉， 李忠勤， 金 爽， 徐春海，梁鵬斌，岳曉英，羊 旻

（1.中國科學院西北生態(tài)環(huán)境資源研究院冰凍圈科學國家重點實驗室/天山冰川觀測試驗站，甘肅蘭州730000；2.中國科學院大學，北京100049；3.西北師范大學地理與環(huán)境科學學院，甘肅蘭州730070）

0 引言

冰川（包括冰蓋）作為冰凍圈要素之一，因其對氣候的高度敏感性和重要的反饋作用而倍受關(guān)注。隨著氣候變暖，冰川變化對海平面、淡水資源、生態(tài)、環(huán)境、氣候等產(chǎn)生了深刻的影響。根據(jù)IPCC第五次評估報告，末次盛冰期以來，全球海平面平均上升了120 m，目前全球山地冰川對海平面上升的貢獻為（0.76±0.37）mm·a-1，是僅次于海洋熱膨脹的海平面上升第二大貢獻者［1］。近50年來，我國天山山脈冰川面積縮小了18%［2］，冰川的退縮對區(qū)域生態(tài)環(huán)境、綠洲經(jīng)濟發(fā)展、工農(nóng)業(yè)生產(chǎn)以及居民用水等造成至關(guān)重要的影響［3-4］。

冰川物質(zhì)平衡是表征冰川積累和消融量值的重要冰川學參數(shù)之一，主要受控于能量收支狀況，對氣候變化有敏感的響應。物質(zhì)平衡及其動態(tài)變化是引起冰川規(guī)模和徑流變化的物質(zhì)基礎(chǔ)，是連接冰川與氣候、冰川與水資源的重要紐帶［5］。傳統(tǒng)的物質(zhì)平衡監(jiān)測主要利用花桿/雪坑法，但受限于高海拔、陡峭地勢，耗時費力。近年來利用不同時期數(shù)字高程模型（DEM）差值來計算物質(zhì)平衡，使得評估大尺度、長時間序列的冰川物質(zhì)平衡成為可能［6］，但在認識冰川變化的物理機制及對氣候變化響應過程和機理方面存在一定的限制。

目前，國內(nèi)外冰川學聚焦于冰川物質(zhì)平衡模擬。模型模擬主要有兩種方法。一種是采用半經(jīng)驗的方法，建立冰川物質(zhì)平衡和氣象指標的關(guān)系，如線性回歸模型［7］、度日因子模型［8］等。Rasmussen等［7］利用線性回歸模型在斯瓦爾巴群島的模擬結(jié)果指出，該模型對冬季平衡的模擬效果良好，同時該模型操作簡單能夠重建物質(zhì)平衡，已經(jīng)在挪威、冰島和瑞典等地得到了應用。Huintjes等［9］利用修正的度日因子模型在烏魯木齊河源1號冰川東支計算了物質(zhì)平衡，發(fā)現(xiàn)在引入潛在太陽短波輻射修正模型后，模擬的物質(zhì)平衡和消融速率有了顯著的改善，但該方法只涉及表面物質(zhì)平衡的計算，在氣候條件復雜的大區(qū)域中是不適用的［10］。另一種方法是基于能量平衡描述冰川物理過程的物質(zhì)平衡模型，如分布式能量-物質(zhì)平衡模型［11］和表層能量-物質(zhì)平衡模型（SEMB）［12］。Hock等［11］通過改進反照率參數(shù)化方案以及輻射分量計算方法，考慮地形因素如坡度、坡向等，建立的分布式能量-物質(zhì)平衡模型能夠?qū)πr尺度的冰雪消融進行準確模擬。與度日因子模型相比，該方法則更為細致地體現(xiàn)了物質(zhì)平衡的計算。此外，Yang等［13］、Li等［14］、Zhu等［15］利用表層能量-物質(zhì)平衡模型在藏東南地區(qū)進行了系統(tǒng)的模擬研究，關(guān)注云量和反照率對物質(zhì)平衡的影響，同時結(jié)合區(qū)域大氣環(huán)流如南亞季風，探討其與物質(zhì)平衡的關(guān)系，進一步在區(qū)域尺度上分析了能量-物質(zhì)平衡的差異性，從而揭示冰川對區(qū)域氣候的響應過程。這些模型僅考慮了冰川表面的能量平衡分布，而由德國柏林洪堡大學研制開發(fā)的COSIMA模型，則充分考慮了大氣、冰川表面以及表面以下10 m內(nèi)雪層的物質(zhì)和能量交換過程，通過建立物質(zhì)平衡與氣象要素（如氣溫、降水、輻射等）變化之間的關(guān)系，實現(xiàn)對整條冰川規(guī)模、小時尺度的冰川物質(zhì)平衡模擬計算。目前該模型已經(jīng)在青藏高原扎當冰川［16］和普若崗日冰帽［17］等地進行了能量-物質(zhì)平衡模擬應用，取得了很好的模擬效果。

鑒于此，基于烏魯木齊河源1號冰川東支2018年消融期冰面氣象站和物質(zhì)平衡花桿觀測數(shù)據(jù)，利用COSIMA模型，開展了冰川能量-物質(zhì)平衡模擬研究，以期為模擬預估由冰川變化引發(fā)的水資源與水循環(huán)、生態(tài)環(huán)境等變化，同時揭示冰川消融機理，為研發(fā)延緩冰川消融、保護冰川的措施提供科學依據(jù)。

1 研究區(qū)概況

烏魯木齊河源1號冰川（43°06′N、86°49′E）位于我國天山中部喀拉烏成山脈主脈北坡烏魯木齊河源上游（圖1），屬冰斗山谷冰川，由于冰川消融，于1993年分離為東、西兩支。2015年1號冰川總面積為1.558 km2，其中東支面積為1.058 km2，海拔介于3 743～4 267 m［18］。受夏季氣溫升高、反照率下降、輻射增強及冰川溫度升高影響，近50年來1號冰川在面積、末端位置及冰川溫度等方面發(fā)生了顯著變化，其中末端位置和面積變化率分別為4.2 m·a-1、0.006 km2·a-1［19-20］。該區(qū)域主要受對流層上部西風急流、西伯利亞反氣旋環(huán)流和西風環(huán)流氣旋擾動三大動力要素控制［2］。大西溝氣象站（海拔3 593 m）觀測數(shù)據(jù)和物質(zhì)平衡觀測數(shù)據(jù)顯示，1959—2017年烏魯木齊河源1號冰川年均氣溫為-4.6℃，年均降水量為460 mm。5—8月是主要的消融期，降水主要來源于西風攜帶的水汽，其中5—8月降水占年降水量的78%，降水類型主要為固態(tài)降水，如雪、雹及雨夾雪［21］。

2 數(shù)據(jù)與方法

2.1 氣象數(shù)據(jù)

利用冰面自動氣象站AWS1獲取消融期（2018年4月29日14：00至9月1日10：00，北京時間）小時尺度氣象數(shù)據(jù)（圖1）。該氣象站架設(shè)于海拔4 025 m處，接近東支主流線，地勢相對平坦（坡度＜2°）。觀測的氣象要素包括氣溫、相對濕度、氣壓、風速、降水量、入射和出射的長短波輻射等，數(shù)據(jù)集存儲在Campbell CR1000中，每十分鐘記錄一個氣象要素。具體氣象站傳感器類型和指標見表1。

然而，在復雜冰面氣象環(huán)境的影響下，AWS1在觀測期間出現(xiàn)了48 h的數(shù)據(jù)缺失，因此，為了填補這一數(shù)據(jù)空白來獲得連續(xù)的氣象數(shù)據(jù)集，本文利用AWS2和AWS3相應時間段內(nèi)的數(shù)據(jù)來進行補充，具體方法見2.3節(jié)。AWS2是安裝在海拔3 882 m花桿D2處距離冰面高度2 m的便攜式氣象儀。AWS3位于海拔3 835 m冰川監(jiān)測塔處，距離1號冰川末端約100 m，該氣象站具備同時段30 min尺度的入射及出射長短波輻射、氣溫、相對濕度、風速及降水量的觀測記錄。考慮到輻射傳感器頂部會因頻繁的降雪而形成積雪和霜，從而影響短波輻射的觀測，因此，本文利用van den Broeke等［22］的方法對短波輻射進行修正，該方法則可以使凈短波輻射的計算精度提高至小時尺度?；谂R界氣溫的參數(shù)化方法來判斷降水的狀態(tài)，COSIMA中采用雙臨界氣溫指數(shù)法進行固液態(tài)降水分離：當氣溫高于5℃時為降雨，低于1℃時為降雪，氣溫介于兩臨界氣溫之間則按線性插值計算［23］。考慮到風速對固態(tài)降水的影響，采用T-200B原始降水數(shù)據(jù)的修正方法［24］，即

圖1 烏魯木齊河源1號冰川觀測示意圖Fig.1 The map showing the location of the Urumqi Glacier No.1（a）and the observation system on the east branch（b），together with the photos of automatic weather station（AWS）（c），weighing precipitation/snow gauge（d）and ablation stake（e）

表1 自動氣象站傳感器信息Table 1 Sensor information of the AWSs used in this study

式中：CE為捕捉率（%）；Ua為風速（m·s-1）。一般認為液態(tài)降水的捕捉率為90%，進行液態(tài)降水的修正。

2.2 冰川物質(zhì)平衡數(shù)據(jù)

為了進行冰川物質(zhì)平衡的驗證，在AWS1同海拔附近布設(shè)了5個花桿點（圖1）。2018年消融期，每隔半個月進行了詳細的觀測，選取消融期的7次，分別是6月20日、7月20日、8月1日、8月10日、8月14日、8月18日以及9月1日。具體每根花桿物質(zhì)平衡觀測內(nèi)容包括：花桿至冰川表面的垂直高度、附加冰的厚度、粒雪層的厚度和密度及雪坑剖面結(jié)構(gòu)。此外，密度測量必須要求雪坑垂直高度大于5 cm，否則冰雪密度采用野外觀測經(jīng)驗值，冰和雪的密度分別取值為900 kg·m-3和300 kg·m-3。單點物質(zhì)平衡計算公式為

式中：bs、bice、bsi分別為雪、冰川冰以及附加冰的物質(zhì)平衡，具體參見文獻［25］。由此獲得5根花桿的單點物質(zhì)平衡值，取其平均值作為冰川物質(zhì)平衡的驗證數(shù)據(jù)。

然而，在野外觀測過程中，由于氣候和人為等因素會造成物質(zhì)平衡的觀測存在誤差，如觀測人員花桿讀數(shù)、雪/粒雪密度的測量及假設(shè)、花桿的傾倒及下沉、雪/粒雪層的誤判、附加冰難以辨別以及重復測量等［26-27］。綜上所述，據(jù)Thibert等［26］提出的單點物質(zhì)平衡誤差計算方法，該方法綜合考慮了單點物質(zhì)平衡測量時因雪/粒雪層誤判以及重復測量造成的誤差，最終確定冰川物質(zhì)平衡計算結(jié)果誤差為±0.16 m w.e.。

2.3 COSIMA模型

基于COSIMA（COupled Snowpack and Ice surface energy and MAss balance model）開展1號冰川的單點物質(zhì)平衡模擬，該模型是由表面能量平衡模型和表層下雪冰模型耦合形成的。其中，表面能量平衡模型主要基于如下能量平衡等式。

式中：SWin、LWin、LWout、Qsens、Qlat、QG、α分別為入射短波輻射、入射與出射長波輻射、感熱與潛熱通量、地熱通量以及反照率。需要說明的是，只有當冰川表面溫度達到消融點時消融才會發(fā)生，此時F的數(shù)值就是各能量項的代數(shù)和，規(guī)定方向指向冰川表面的能量通量為正值。

COSIMA中利用塊體空氣動力學方法計算湍流熱通量［28］，則冰川表面感熱通量和潛熱通量可表示為

式中：ρair為空氣密度（kg·m-3）；K為卡爾曼常數(shù)（0.4）；hz為儀器架設(shè)高度，本文為2 m；zo為表面粗糙長度，新雪、粒雪和冰的粗糙長度分別為0.24 mm、4 mm及1.7 mm［29-31］；u為風速（m·s-1）；Tair、Ts分別為2 m處和冰川表面氣溫（℃）；LE為蒸發(fā)潛熱（2.514×106J·kg-1）；LS為升華潛熱（2.849×106J·kg-1）；qair與qs分別為2 m處和冰川表面比濕（kg·kg-1）；RHair/s為2 m處和冰川表面相對濕度（%），其中冰川表面的相對濕度假設(shè)為100%；p為氣壓（hPa）；E和Es分別為2 m處和冰川表面飽和水汽壓（hPa）。

針對湍流通量計算的不確定性，COSIMA利用Braithwaite［32］描述方法對湍流通量進行穩(wěn)定性修正，如下所示。

式中：Qsens/lat和Q'sens/lat分別為穩(wěn)定性修正后和修正前的感熱和潛熱通量；Ri為總體理查遜數(shù)；g為重力加速度（9.81 m·s-2）；Tair（c）、Tair分別為攝氏度氣溫（℃）和開爾文氣溫（K）。

COSIMA模型表層下雪冰模型部分是將冰面以下10 m范圍內(nèi)的距離以0.2 m等分，由于每一層的水熱組成存在差異，從而使得每一層的物理性質(zhì)（冰溫、密度以及液體水含量）不盡相同［16］?；谝巴鈱崪y數(shù)據(jù)，初始雪深設(shè)置為0.32 m。初始溫度剖面利用1號冰川實測的表面溫度和冰川底部溫度（-5.2℃）通過線性插值獲得。初始的雪坑密度剖面也是基于1號冰川密度經(jīng)驗值以插值的方式得到，其中表層雪和粒雪密度分別采用經(jīng)驗值280 kg·m-3和265 kg·m-3；雪坑以下冰川冰的密度采用經(jīng)驗值900 kg·m-3。

采用COSIMA模型進行1號冰川物質(zhì)平衡模擬研究具有如下兩方面優(yōu)勢：第一，在各能量通量作用下，物質(zhì)平衡以物質(zhì)流的形式發(fā)生變化，COSIMA模型可以將各物質(zhì)流作用下導致的逐小時物質(zhì)平衡變化進行區(qū)分，實現(xiàn)物質(zhì)平衡精細化研究；第二，與以往只考慮表層物質(zhì)平衡相比，COSIMA模型可以較為系統(tǒng)的呈現(xiàn)物質(zhì)平衡在表層下的變化過程，比如融水滲浸、再凍結(jié)以及密實化過程［16］，使得物質(zhì)平衡的計算更為全面、合理。

下面詳細介紹模型中需要的數(shù)據(jù)，輸入數(shù)據(jù)包括入射短波輻射、云因子、氣溫（2 m）、降水量、風速（2 m）、相對濕度（2 m）及氣壓。并且在復雜冰面氣象環(huán)境的影響下，AWS1在8月20日12：00至8月22日12：00內(nèi)的數(shù)據(jù)缺失，為了獲得連續(xù)的氣象數(shù)據(jù)集，本文利用AWS2及AWS3相應時間段內(nèi)的各氣象要素來進行補充。

（1）入射短波輻射

COSIMA模型中入射短波輻射由CNR4傳感器直接測得。雖然受大氣透明度及復雜地形的影響，入射短波輻射在空間上呈現(xiàn)差異性［11，33］，但是，入射短波輻射在較小海拔間仍具有強相關(guān)性，因此，嘗試建立了AWS1與AWS3入射短波輻射的擬合關(guān)系，兩者擬合相關(guān)系數(shù)達0.95，使用AWS3的入射短波值來代替缺失時段。

（2）云因子

云對輻射的吸收和反射程度不同，其中10%～20%輻射量被云吸收，同時36%～80%的輻射量被云反射［34］，所以考慮云對短波輻射的影響是十分必要的。由于缺乏云量數(shù)據(jù)的直接觀測，選擇利用式（12）計算得到云因子，其中SWTOA為大氣頂層太陽輻照度，取其常數(shù)為1 368 W·m-2［35］。

（3）氣溫和降水量

氣溫和降水量作為氣候變化最為直接的參數(shù)，其變化程度深刻影響著大氣與冰川區(qū)之間的能量交換，進而影響冰川區(qū)物質(zhì)的積累和消融。由于存在缺失數(shù)據(jù)，因此，利用AWS1、AWS2和AWS3的氣溫和降水數(shù)據(jù)，考慮海拔梯度-1.1℃·（100m）-1、1.6 mm·（100m）-1，填補AWS1氣溫和降水量空白值。

（4）風速、相對濕度及氣壓

近地面湍流交換受大氣和冰川表面濕度差及風速等因素的影響［36］，利用Hock等［11］、Braithwaite等［37］的方法插補AWS1風速，其理論依據(jù)是假設(shè)整個冰川面上風速相同的情況下建立氣象站點之間的關(guān)系，其中，Hock等［11］通過比較Storglacia¨ren冰川上分布在不同海拔的風速，發(fā)現(xiàn)各氣象站間散點分布明顯，風速的高低與各站點所處位置有關(guān)，因此，無法得出一個一般性的方案來定量分析風速在冰川面上的空間變化。所以，在考慮地形對風速的影響下，建立布設(shè)在冰舌和平衡線處風速的正相關(guān)關(guān)系，獲得了良好的分布式模擬結(jié)果，因此，利用6月24日至8月20日AWS1和AWS3風速，其相關(guān)系數(shù)為0.41，雖然兩站點風速相關(guān)性較低，但也可用于插補AWS1站點缺失值。相對濕度的插補同樣依據(jù)Hock等［11］描述的方法，研究指出相對濕度在氣象站點間的變化相對較小，所以，可以假設(shè)兩站點間相對濕度是相同的，基于這一假設(shè)可利用氣象站點相對濕度拓展為整個冰川面上的相對濕度。并且實際情況中AWS1與AWS2站點之間距離較近，因此，可利用AWS2同時段相對濕度獲得AWS1連續(xù)數(shù)據(jù)序列。為了精確地獲得氣壓與高度的對應關(guān)系，利用等溫大氣壓高方程插補缺失的氣壓值，式中是把大氣當成干空氣處理的，并且實際大氣并非等溫大氣，重力加速度也隨緯度和高度而有變化，但在大氣底層（＜100 km）重力加速度隨高度的變化不大，通過計算標準大氣中氣壓與高度的對應值，發(fā)現(xiàn)氣壓仍是滿足隨高度按指數(shù)遞減的規(guī)律，且等溫大氣與實際大氣非常接近［38］。

式中：PAWS1、PAWS3分別為AWS1和AWS3處氣壓值；t為氣層間平均氣溫；HAWS1、HAWS3為AWS1和AWS3處海拔，分別為4 025 m和3 835 m；α為常數(shù)，約為

0.036。

2.4 誤差評估

根據(jù)Huintjes等［16］研究，COSIMA模擬物質(zhì)平衡的日誤差（Uday）評估采用如下方法。

式中：Xj，max、Xj，min分別為利用6次測量間隔獲得的降水梯度運行COSIMA得到的物質(zhì)平衡最大值和最小值；nj為測量間隔j的天數(shù)。利用降水梯度［（1.6±0.3）mm·（100m）-1］分 別 運 行COSIMA，其 中，1.6 mm·（100m）-1作為參照降水梯度。綜合評估表明，消融期冰川物質(zhì)平衡的日誤差為±0.03 m w.e.。

3 結(jié)果與分析

3.1 消融期冰面氣象狀況

自動氣象站觀測結(jié)果顯示，2018年消融期，1號冰川冰面氣溫波動存在顯著的差異性。消融期初（4月29日至6月8日），氣溫波動幅度劇烈，日均值為-3.2℃；6月9日至8月28日，氣溫大于0℃；消融期末（8月19日至9月1日）氣溫持續(xù)遞減［圖2（a）］。整個消融期相對濕度均值約為69%，其變化態(tài)勢與氣溫較為類似［圖2（b）］，其中，5月相對濕度低，且低于均值的比重為91%。整個消融期，風速在0.6～5.3 m·s-1內(nèi)波動變化，均值為2.3 m·s-1［圖2（c）］，以6月27日為分界，4月29日至6月27日風速相對較高，之后變化小。氣壓在消融期內(nèi)的變化相對穩(wěn)定，均值為625 hPa［圖2（d）］。消融期降水量為692.2 mm，主要降水集中在6—8月［圖2（e）］。

圖2 2018年消融期氣象狀況日均值Fig.2 Daily mean values of meteorological condition during the ablation season in 2018：air temperature（a），relative humidity（b），wind speed（c），air pressure（d），precipitation（e）

3.2 COSIMA模型驗證

利用消融期物質(zhì)平衡模擬結(jié)果對比實測的物質(zhì)平衡，從散點圖結(jié)果可以看出［圖3（a）］，模擬的物質(zhì)平衡值與實測的物質(zhì)平衡值有非常好的一致性，相關(guān)系數(shù)達0.96，離散程度（RMSE）為±0.22 m w.e.。4月29日至8月17日物質(zhì)平衡先緩慢減小，然后受強烈消融影響快速減小，8月17日后物質(zhì)平衡變化逐漸減小。與實測值相比，7月20日與8月1日兩次模擬物質(zhì)平衡被低估，小于實測物質(zhì)平衡，可能源于模型計算的湍流通量被低估，導致冰川消融減少［圖3（b）］。盡管如此，整體來看模擬結(jié)果能很好的呈現(xiàn)消融期1號冰川的物質(zhì)平衡變化，為進一步模擬其能量通量及物質(zhì)平衡的變化特征提供了有效支撐。

3.3 能量通量特征

圖3 模擬與實測累積物質(zhì)平衡Fig.3 Modelled and observed accumulative mass balances：scatter plot（a），change trend（b）

表2 消融期烏魯木齊河源1號冰川東支AWS1處月尺度能量組成變化及占比Table 2 Monthly variation and contributions of energy components at the AWS1 site on the east branch of the Urumqi Glacier No.1 during the ablation season

消融期各能量項月尺度能量組成如表2所示。消融期冰川表面能量收入項中凈短波輻射、感熱通量所占比重的均值為84%、16%，所以，凈短波輻射是冰川表面主要的能量來源，其次是感熱通量。凈短波輻射（SWnet）在月尺度上表現(xiàn)出先減小后增加的趨勢，6月和8月分別出現(xiàn)最小值（+30.48 W·m-2）和最大值（+108.43 W·m-2），影響因素除太陽傾斜角外，實測入射短波輻射在8月顯著增加是其主要原因，而6月出現(xiàn)最小值主要是源于頻繁的降雪事件使得冰川表面的反照率增加，冰川表面接收的太陽輻射減小。在扎當冰川上的研究也得到了類似的結(jié)果［39］，說明消融期冰川表面接收的太陽輻射受表面反照率反饋機制的調(diào)控。感熱通量在消融期內(nèi)始終是正值（+10.89 W·m-2），表明冰川上空2 m處氣溫高于冰川表面溫度，大氣向冰川表面?zhèn)鬏斈芰浚?9］。而凈長波輻射、潛熱通量、用于冰雪消融的能量及地熱通量作為能量支出項所占比重依次為55%、7%、32%及6%，可以看出，能量的支出項主要是凈長波輻射和用于冰雪消融的能量。消融期，凈長波輻射月均值為-39.3 W·m-2，6—7月相對濕度較高（RH=69%），并且氣溫和模擬冰川表面溫度等于0℃，導致凈長波輻射（LWnet）在6—7月相對較高且存在少數(shù)正值（圖4）。觀測期間，氣溫較高，空氣水汽壓小于冰川表面水汽壓，使得多數(shù)潛熱通量在消融期內(nèi)是負值，表明冰川表面以升華和蒸發(fā)為主，同時向大氣輸出熱量［23］，潛熱通量從負值轉(zhuǎn)變?yōu)檎档奶鞌?shù)中，地氣溫差和風速較?。▓D2），不利于湍流的發(fā)展，此時大氣向冰川表面?zhèn)鬏斈芰浚銜l(fā)生凝結(jié)或凝華［40］。地熱通量（QG）是透射短波輻射（Qps）和熱傳導通量（QC）的和，均值為-4.85 W·m-2，此外，受模型中冰川最底部冰溫為常數(shù)的影響，極有可能高估了地熱通量。當冰溫達到0℃時，透射短波輻射將會用于冰下消融。此外，消融期由升華或蒸發(fā)消耗的能量比例在46%（5月）至5%（8月）之間，均值為30%。

圖4 消融期烏魯木齊河源1號冰川東支日均能量組成狀況Fig.4 Daily energy components on the east branch of the Urumqi Glacier No.1 during the ablation season

1號冰川是典型的大陸型冰川，為了進一步研究該冰川能量通量的變化特征，將其與我國其他大陸型冰川已有研究開展對比分析（表3）。雖然不同研究所涉及的年份的氣象狀況各異、研究時段和模擬方法存在差異，但仍然能夠從整體上揭示大陸型冰川能量通量的變化特征。與1號冰川早期能量平衡研究相似［41-42］，造成冰川消融的能量主要來源于凈短波輻射，消融冰面由于其逆溫層結(jié)，感熱通量占比次之；冰川能量支出主要為凈長波輻射、冰川消融耗熱。由于消融冰面的蒸發(fā)受制于表面水汽壓，使得冰面蒸發(fā)受到一定的抑制作用，潛熱通量負值較大。與1號冰川早期能量平衡參數(shù)化方案相比，COSIMA更具系統(tǒng)性，且考慮冰面以下10 m范圍內(nèi)的滲浸短波輻射和熱傳導通量，更為細致地反映能量物質(zhì)平衡過程機理。同一冰川不同海拔位置，海拔越高，凈輻射越小［41-46］，凈輻射值也存在明顯的年際變化，同時，不同時期各能量通量的差異涉及冰面氣象狀況、參數(shù)化方案及研究時段等因素。此外，李晶等［43］指出表磧表面凈輻射小于雪面凈輻射，說明冰川下墊面差異對凈輻射的影響也是不容忽視的。利用扎當冰川與其他區(qū)域冰川做比較，發(fā)現(xiàn)凈輻射是扎當冰川最大的能量收入項（+39 W·m-2），占總的能量收入項的82%，同時考慮到扎當冰川的研究點海拔顯著高于其他冰川，所以，不能將扎當冰川與其余冰川進行對比，但是凈輻射在不同類型冰川上的比重差異是明顯的［47］。相似海拔下冰川所處氣象狀況可直接影響湍流通量的大小，其中七一冰川（+14.2 W·m-2）及科契卡爾巴西冰川（+14.4 W·m-2）的感熱通量大于烏魯木齊河源1號冰川（+10.89 W·m-2）及老虎溝12號冰川（+5.7 W·m-2），表明氣溫越高，風速越大，感熱通量越大［43，46，48］，而潛熱通量的低值說明潛熱交換劇烈消耗的能量多，與其余冰川相比，消融期科契卡爾巴西冰川潛熱交換劇烈（-23 W·m-2），說明科契卡爾巴西冰川表面飽和水汽壓和空氣水汽壓差異明顯，水汽壓梯度最大，蒸發(fā)最為強烈［43］。消融能量也具有隨海拔增加而減小的特征，但只涉及烏魯木齊河源1號冰川以及老虎溝12號冰川。此外，基于扎當冰川和帕隆4號冰川及其他區(qū)域典型監(jiān)測冰川能量通量的對比結(jié)果表明單條冰川所處氣候環(huán)境對消融能量的影響的確存在，但難以定量化，而要準確衡量對消融能量的影響程度還有賴于從控制研究時段和采用相同模型等兩方面進行著手研究［47］。

表3 烏魯木齊河源1號冰川與其他大陸型冰川能量組成對比Table 3 Comparison of energy components of the Urumqi Glacier No.1 with other continental glaciers

3.4 物質(zhì)平衡特征

冰川物質(zhì)平衡是冰川積累和消融的代數(shù)和，取決于降水、消融、升華、蒸發(fā)以及再凍結(jié)過程。1號冰川累積物質(zhì)平衡及其組成顯示（表4），模擬的消融期累積物質(zhì)平衡為-0.671 m w.e.，其中，表面消融是最主要的物質(zhì)損失項（-0.742 m w.e.），其次是冰下消融（-0.114 m w.e.），升華/蒸發(fā)也會引起物質(zhì)損失，但值相對較?。?0.022 m w.e.）。固態(tài)降水是主要的物質(zhì)收入項（0.196 m w.e.），再凍結(jié)帶來的物質(zhì)收入僅為0.007 m w.e.?？偟膩碚f，幾乎全部表層和表層下的消融量都隨著冰川徑流而損失掉。由于地表水汽壓小于空氣水汽壓，加之風速較小，產(chǎn)生了0.005 m w.e.水汽凝結(jié)量。Weidemann等［49］指出積雪場的能量收支會受到風吹雪的影響，COSIMA模型中缺少風吹雪的參數(shù)化方案，因此會導致上述升華作用被高估，但消融期1號冰川消融強烈，冰川表面積雪性狀會因氣溫和融水等因素發(fā)現(xiàn)改變，因此，只有當風速較大且頻繁有降雪發(fā)生時升華作用則被認為高估。此外，Oerlemans等［50］在西格林蘭冰蓋及瑞士的冰川上的研究指出冰川重力流在消融冰川中非常常見，能夠增強湍流交換，而在應用COSIMA模型時忽略了冰川重力流，因此可能會造成表面消融量的低估。雖然會存在上述問題，但總體來看，消融期1號冰川物質(zhì)平衡的變化主要取決于表面消融和固態(tài)降水。

表4 烏魯木齊河源1號冰川累積物質(zhì)平衡組成及與其他大陸型冰川對比Table 4 Accumulative mass balance components of the Urumqi Glacier No.1 and comparison with other continental glaciers

為進一步研究烏魯木齊河源1號冰川物質(zhì)平衡組成，將其與其他區(qū)域監(jiān)測的大陸型冰川進行比較。結(jié)果表明，與羌塘1號冰川、扎當冰川相比，烏魯木齊河源1號冰川物質(zhì)損失明顯（表4）。夏季積累型冰川深受降水的季節(jié)變化及固態(tài)降水比重的影響［51］，而烏魯木齊河源1號冰川固態(tài)降水顯著小于其余兩者，主要原因是羌塘1號冰川和扎當冰川分別在消融期深受西風和印度季風的影響，使得固態(tài)降水量增加。另一方面，盡管表面消融和升華/蒸發(fā)作用差異很小，但烏魯木齊河源1號冰川再凍結(jié)作用顯著小于其余兩者。已知再凍結(jié)作用主要發(fā)生在消融期［40］，取決于冰/積雪溫度及其物理性質(zhì)［52］，并且冰/積雪溫度越高，再凍結(jié)過程越弱［16］，所以，從定性的角度上可以得出烏魯木齊河源1號冰川的冰/積雪溫度高于其余兩者，而再凍結(jié)過程的物理機制還有待進一步實測驗證。

4 結(jié)論與展望

針對烏魯木齊河源1號冰川東支，以氣象數(shù)據(jù)驅(qū)動COSIMA模型，結(jié)合冰川物質(zhì)平衡觀測數(shù)據(jù)，實現(xiàn)了該冰川2018年消融期單點能量-物質(zhì)平衡模擬。得到如下結(jié)論：

（1）物質(zhì)平衡模擬值為（-0.67±0.03）m w.e.，模擬與實測物質(zhì)平衡較為一致，相關(guān)系數(shù)達0.96，因此，COSIMA應用在烏魯木齊河源1號冰川上能夠較為準確地計算單點能量-物質(zhì)變化。

（2）凈短波輻射是主要的能量來源（84%），其次為感熱通量（16%）。凈長波輻射是主要的能量支出項（55%），其次分別是冰川消融耗熱（32%）、潛熱消耗（7%）及地熱通量（6%）。受能量收支的影響，模擬的冰川物質(zhì)平衡主要取決于表面消融和固態(tài)降水量，同時在模型自身局限性的影響下，模擬升華和表面消融量存在不可避免的誤差。

（3）與其他大陸型冰川已有研究結(jié)果對比發(fā)現(xiàn)，不同海拔以及氣候條件下能量通量的差異明顯，同時，物質(zhì)平衡對比結(jié)果表明烏魯木齊河源1號冰川固態(tài)降水量和再凍結(jié)顯著小于羌塘1號冰川和扎當冰川，推測與單條冰川所處大氣環(huán)流密切有關(guān)。

鑒于此，為了更為深入地研究1號冰川能量-物質(zhì)平衡特征及其對氣候變化的響應，深刻認識氣候機制作用下的冰川變化過程及機理，下一步將開展單點到分布式的能量-物質(zhì)平衡模擬研究。第一，完善觀測數(shù)據(jù)。需要在1號冰川消融區(qū)補充氣象和四分量輻射數(shù)據(jù)的觀測儀器，以便獲得更高精度的分布式輸入數(shù)據(jù)。第二，完善觀測技術(shù)。除常規(guī)觀測數(shù)據(jù)外，國際上對于冰川能量-物質(zhì)平衡的模擬驗證通常將遙感影像和延遲攝影技術(shù)結(jié)合，提取不同時間同一時刻雪線的瞬時變化，來獲得研究期觀測及模擬雪線的時空變化，借以驗證模型模擬效果，從而有助于提高模型模擬精度。第三，修正輸入數(shù)據(jù)以驅(qū)動模型，獲得模擬結(jié)果。基于1號冰川DEM數(shù)據(jù)提取掩膜，利用WRF模式修正氣象要素以獲得不同海拔上氣象要素梯度變化，通過雙曲正切函數(shù)計算固態(tài)降水量比重，最后結(jié)合入射短波輻射數(shù)據(jù)驅(qū)動二維模型，獲得1號冰川分布式能量-物質(zhì)平衡模擬結(jié)果。

致謝：感謝西北師范大學沈思民碩士研究生在野外觀測和論文撰寫過程中的大力幫助。