水穩(wěn)定同位素示蹤的冰川流域水文模擬及不確定性研究
——以烏魯木齊河源1號(hào)冰川為例

何天豪，高紅凱，2，李向應(yīng)，韓添丁，賀志華，張志才，段崢，劉敏，2，丁永建，9

（1.華東師范大學(xué)地理科學(xué)學(xué)院，上海200241；2.華東師范大學(xué)地理信息科學(xué)教育部重點(diǎn)實(shí)驗(yàn)室，上海200241；3.西北大學(xué)陜西省地表系統(tǒng)與環(huán)境承載力重點(diǎn)實(shí)驗(yàn)室，陜西西安710127；4.西北大學(xué)城市與環(huán)境學(xué)院，陜西西安710127；5.中國科學(xué)院西北生態(tài)環(huán)境資源研究院冰凍圈科學(xué)國家重點(diǎn)實(shí)驗(yàn)室，甘肅蘭州730000；6.Centre for Hydrology，University of Saskatchewan，Saskatoon SK S7N 5C9，Saskatchewan，Canada；7.河海大學(xué)水文水資源學(xué)院，江蘇南京210098；8.Department of Physical Geography and Ecosystem Science，Lund University，Lund 7222，Sweden；9.China-Pakistan Joint Research Center on Earth Sciences，CAS-HEC，Islamabad 45320，Pakistan）

0 引言

冰川作為固體水庫，顯著減少了徑流變差系數(shù)，對(duì)穩(wěn)定寒旱區(qū)水資源起到了重要作用［1］。冰川也具有遠(yuǎn)高于非冰川區(qū)的徑流系數(shù)，為我國西北干旱區(qū)工農(nóng)業(yè)生產(chǎn)和生活用水，以及生態(tài)環(huán)境保護(hù)提供了重要水源［2-3］。因此，冰川水文效應(yīng)的研究不僅是重要的科學(xué)問題［4-7］，還關(guān)系到干旱區(qū)水資源管理相關(guān)政策的制定［4，8］。模型研究是科學(xué)認(rèn)識(shí)寒區(qū)水文過程機(jī)理的重要手段［9-10］。水文模擬也是量化各個(gè)徑流組分的有效方法，還是科學(xué)預(yù)估未來變化的重要工具［11-12］。目前，水文模型已成功應(yīng)用于研究不同地區(qū)高寒流域各水源的空間變化和時(shí)間演變，以及融雪融冰和降雨徑流對(duì)總徑流量貢獻(xiàn)的模擬等研究［13-14］。然而，冰川地區(qū)由于特殊的地理?xiàng)l件和復(fù)雜的水文過程，傳統(tǒng)水文學(xué)方法在該區(qū)域并不適用［15-17］。另外，眾多模型研究采用單一指標(biāo)（流量）進(jìn)行模型參數(shù)率定，容易造成異參同效［18］，導(dǎo)致徑流模擬和各徑流成分的量化仍然存在很大不確定性［19］。

水穩(wěn)定同位素作為徑流數(shù)據(jù)之外的重要信息，為限制模型不確定性提供了新的約束。穩(wěn)定氫氧同位素（2H和18O）具有“保守”特征，并且穩(wěn)定氫氧同位素比率在流域尺度各水體中差異明顯［20-22］，具有進(jìn)一步限制徑流組分互相妥協(xié)的潛力［15，23］。也有研究證實(shí)了同位素示蹤數(shù)據(jù)具有減少參數(shù)不確定性的能力，從而有助于改善模型結(jié)構(gòu)［24-25］。然而，現(xiàn)有的研究主要是定性研究和以端元法為主的半定量研究［26-28］，缺乏對(duì)小流域穩(wěn)定氫氧同位素時(shí)空尺度效應(yīng)的整體性和系統(tǒng)性認(rèn)識(shí)，限制了其作為重要信息源的作用。穩(wěn)定氫氧同位素示蹤水文模型則能提供更有效的研究手段。穩(wěn)定氫氧同位素示蹤水文模型在非寒區(qū)降雨徑流模擬中已有較多研究［23，29-30］，但是寒區(qū)同位素示蹤水文模型卻不多見［25，31］。冰川流域的同位素水文模擬研究在國際上也處于剛起步階段［15，18］。因此，本文將分布式冰川流域水文模型（FLEXG）和同位素分餾混合過程進(jìn)行了耦合，構(gòu)建穩(wěn)定氫氧同位素示蹤的冰川流域水文模型（FLEXG-iso），于烏魯木齊河源1號(hào)冰川流域進(jìn)行了模擬檢驗(yàn)和徑流分割，以期加深對(duì)寒區(qū)水文過程更全面的認(rèn)識(shí)。

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

1.1 研究區(qū)概況

研究區(qū)位于中國西北新疆維吾爾自治區(qū)的烏魯木齊河源頭（43°50′N、86°49′E），處于中國天山山脈天格爾山北坡，深居歐亞內(nèi)陸，具有大陸性氣候和中國西部和亞洲中部干旱區(qū)高寒環(huán)境，對(duì)氣候變化反應(yīng)敏感。該地區(qū)山勢險(xiǎn)峻，海拔范圍為3 712～4 457 m。烏魯木齊河源區(qū)年平均氣溫-5.2℃，最熱月份一般為7月，平均氣溫約5℃，最冷月份一般為1月，平均氣溫約-15℃。一年中溫度低于0℃的時(shí)段長達(dá)7～8個(gè)月。年均降水量約441 mm，一般集中在消融季（5月至10月），消融季降水量占年降水量的90%［32］。自觀測以來，年均溫和年降水呈持續(xù)升高趨勢，河源區(qū)大西溝氣象站數(shù)據(jù)顯示，1959—1996年的年平均溫度-5.3℃，較1997—2010年的-4.3℃增加了1℃，并且各個(gè)季節(jié)溫度都比以往顯著升高。1986—2010年平均降水量約為493 mm，較1959—1985年的426 mm增加了約15.7%［33］。烏魯木齊河源1號(hào)冰川（圖1）是朝東北的山谷冰川，為該流域面積最大的冰川。在氣候變暖背景下，1號(hào)冰川自20世紀(jì)90年代以來表現(xiàn)為加速退縮趨勢，冰川面積、長度、物質(zhì)平衡和平衡線高度都發(fā)生了顯著變化［34］。冰川的總面積正逐年減小，冰川面積由1962年的1.95 km2退縮為2017年的1.54 km2［35］。平衡線高度由1959—1996年的年均4 042 m上升至1997—2008年的年均4 106 m［36］。由于冰川末端不斷消融退縮，于1993年分離為2條獨(dú)立的東西方分支山谷冰川，海拔范圍東支在3 752～4 225 m之間，西支在3 848～4 445m之間，總面積為1.59 km2，平均厚度為44.5 m，最大長度為2.1 km，總儲(chǔ)量為8 334.1×104m3［37］。氣溫升高也導(dǎo)致了1號(hào)冰川表面雪層特征及成冰帶的變化，大量滲入的冰川融水導(dǎo)致粒雪再凍結(jié)，導(dǎo)致冷滲浸帶逐漸被凍結(jié)粒雪和滲浸冰取代［38-39］。變得更為單一的冰內(nèi)、冰下排水道減弱了對(duì)融水的阻滯和貯存作用，加快了匯流過程的速度［40］。在1號(hào)冰川東支頂部源頭區(qū)發(fā)現(xiàn)小型冰面湖的存在［34］。烏源1號(hào)冰川流域冰川融水徑流比例高，已有研究表明，水文站觀測的總徑流中有70%是來自于冰川區(qū)徑流，其中44%來自冰川區(qū)由降水產(chǎn)生的徑流，26%來自冰川自身消融產(chǎn)生的徑流［14］。非冰川地區(qū)的土地覆蓋主要是裸露的土壤/巖石和稀疏的草甸，季節(jié)性積雪覆蓋［17］。

圖1 烏魯木齊河源1號(hào)冰川的位置（左上），流域詳細(xì)海拔、冰川覆蓋區(qū)域等高線以及徑流觀測站（右）Fig.1 Location of Glacier No.1 in headwater of Urumqi River（top left），detailed altitude of the Glacier No.1 basin，contour lines in the ice-covered area，and runoff gauging station in the basin（right）

1.2 數(shù)據(jù)來源

本次研究所用的冰川物質(zhì)平衡和水文氣象數(shù)據(jù)來自2013—2016年的天山冰川觀測試驗(yàn)站年報(bào)數(shù)據(jù)［41-42］。冰川物質(zhì)平衡觀測從1號(hào)冰川冰舌開始自下而上，西支布設(shè)A～J橫剖面，海拔約3 850～4 150 m，東支布設(shè)有A～K橫剖面，海拔約3 800～4 150 m。并于每個(gè)剖面等距離布設(shè)測桿，一般為3根。對(duì)于積累區(qū)的觀測，采用挖取雪坑的方法。觀測始于每年的4月，每月底或月初觀測一次，至8月底為完整的年度觀測［43］。2013—2016年逐日的水文氣象數(shù)據(jù)（徑流、降水和溫度）來自于烏魯木齊河源區(qū)1號(hào)冰川水文站。1號(hào)冰川水文站設(shè)在離1號(hào)冰川末端300 m的河道上，斷面海拔3 695 m。穩(wěn)定同位素?cái)?shù)據(jù)（δ18O）由中國科學(xué)院冰凍圈科學(xué)國家重點(diǎn)實(shí)驗(yàn)室提供［44-46］。2013—2016年的消融季，在1號(hào)冰川斷面采集了融水徑流，并于氣象觀測站采集了降水。冰面樣品于2016年消融季采集。地下水、快響應(yīng)徑流（主要來自冰川融水和非冰川區(qū)的地表徑流）和慢響應(yīng)徑流（主要來自非冰川區(qū)的地下水補(bǔ)給）的δ18O來自于Sun等［9］研究的實(shí)測值。地形數(shù)據(jù)和冰川輪廓［47］來源于全球數(shù)字高程模型（STRMDEM，http：//srtm.csi.cgiar.org，V4，分辨率90 m）和第二次冰川編目數(shù)據(jù)（http：//westdc.westgis.ac.cn）。

1.3 研究方法

1.3.1 基于穩(wěn)定同位素的分布式寒區(qū)水文模型

（FLEXG-iso）的構(gòu)建

（1）水文模塊

模型的寒區(qū)水文部分模擬基于Gao等［17］提出的分布式的FLEXG模型，該模型優(yōu)勢在于充分考慮了流域降水、產(chǎn)匯流等異質(zhì)性，使模型更符合實(shí)際。FLEXG將流域劃分為不同高程帶和坡向，以及冰川區(qū)和非冰川區(qū)兩種水文響應(yīng)單元。

根據(jù)FLEXG模型獲得的流域降水、徑流和蒸發(fā)數(shù)據(jù)，根據(jù)水量平衡原理計(jì)算冰川物質(zhì)平衡。

式中：GMB為冰川物質(zhì)平衡；P為冰川區(qū)降水；Q為冰川區(qū)徑流；E為冰川區(qū)融雪徑流和降雨徑流的蒸發(fā)。本研究未考慮升華、雪崩和風(fēng)吹雪的影響。

（2）同位素模塊

流域水文和同位素過程極為復(fù)雜，模型是對(duì)流域異質(zhì)性過程的概化描述，因此模型假設(shè)和簡化在所難免。同位素水文主要伴隨著分餾和混合這兩個(gè)主要過程。本研究假定降雪和降雨的δ18O相同，且不考慮降水δ18O隨高程的變化。將新降雪和現(xiàn)有積雪的δ18O進(jìn)行混合來更新現(xiàn)有積雪的δ18O（δ18OSW）。通過式（2）中給出的瑞利分餾方法模擬融雪的δ18O（δ18OSM）［48］。然后通過融雪融化的δ18O更新剩余的積雪［δ18OSR，式（3）］。采用分餾的校正因子（CFs）來改善融雪δ18O的模擬［18］。瑞利分餾因子（α）被假定為溫度的函數(shù)［式（4）］［49］。采用另一個(gè)校正因子（CFe），考慮了由地表蒸發(fā)引起的δ18O瑞利分餾過程（公式和融雪一樣）。

式中：SM為融雪量；SW為雪水當(dāng)量；SR為融雪后剩余的雪水當(dāng)量；f為積雪中剩余的雪水當(dāng)量與原雪水當(dāng)量的比值［式（5）］；T為當(dāng)?shù)氐拇髿鉁囟龋↘）。

FLEXG為分布式的冰川水文過程模型，將冰川、積雪、降雨徑流等水文過程基于概化的物理過程概念進(jìn)行描述，比如線性水庫和不同景觀的響應(yīng)單元等［17，50］。流入各水庫的δ18O與現(xiàn)有水庫的δ18O混合，更新各個(gè)水庫的δ18O。直接徑流的δ18O來源于冰川區(qū)和非冰川區(qū)徑流的δ18O。假設(shè)進(jìn)入各水庫的水量與水庫原本水量快速完全混合?？紤]了Ss和Sf的額外儲(chǔ)量，該額外蓄量可以理解為非飽和區(qū)中的蓄水容量或參與降雨混合過程中僅在重力作用下無法自由排出的其他水量，由參數(shù)Us（mm）和Uf（mm）表示［30］。假設(shè)儲(chǔ)量固定，僅參與混合過程。由于野外實(shí)地觀測發(fā)現(xiàn)了冰川區(qū)表面存在小型冰面湖［34］，并且冰川區(qū)至觀測站之間有部分非冰川區(qū)域。由此將這樣的情況概化為在冰川區(qū)Sf，g額外儲(chǔ)量（Uf，g），儲(chǔ)量相比土壤較小。Uf，g的先驗(yàn)范圍假定為0到Uf最小值的1/4。最終得到的模型FLEXGiso。本研究未考慮δ18O隨高程變化。

1.3.2 徑流成分的定量化

Weiler等［51］列出了三種徑流成分來源的定義：①來自源區(qū)的貢獻(xiàn)，即來自每類土地覆蓋的貢獻(xiàn)；②來自產(chǎn)流過程的貢獻(xiàn)（如地表徑流，壤中流和地下水徑流）；③輸入的貢獻(xiàn)（融冰，融雪和降雨）。根據(jù)定義2進(jìn)行基于產(chǎn)流過程［式（6）～（14），F(xiàn)LEXGiso水文模塊）對(duì)總徑流的貢獻(xiàn)估算。總徑流包括地表徑流、壤中流和地下水徑流。本研究中，Qf為融雪和降雨混合后的貢獻(xiàn)量，是非冰川區(qū)的地表徑流和壤中流的總和；Qs為地下水的貢獻(xiàn)量；Qg為融冰和降雨混合后的貢獻(xiàn)量，是冰川區(qū)的地表徑流。因此，本研究將徑流成分定義為融雪、融冰、降雨和地下水，以更加詳細(xì)地定量評(píng)估水文斷面年季尺度徑流成分。假設(shè)水進(jìn)入水箱快速完全混合，使原來水箱的組分比例改變，從而可以得到從水箱流出水的組分比例（等于水箱的組分比例）。各貢獻(xiàn)量均為當(dāng)日產(chǎn)流的部分，有一部分輸入水量會(huì)留在各水箱中，參與下一時(shí)間的混合以及產(chǎn)流、蒸發(fā)等水文過程。于是，得以追蹤整個(gè)過程的各組分變化。

式中：Qin為總輸入；Ms為融雪；Mg為融冰；Pl為降雨；Ps為降雪；Rrf為從積雪內(nèi)部液態(tài)水（Swl）到積雪（Sw）的再凍結(jié)水；Pe為土壤或冰表面產(chǎn)生的徑流；Su為非冰川區(qū)土壤非線性儲(chǔ)層；Ea為實(shí)際蒸發(fā)；Ru為Su的過剩水量；Sf和Ss為快響應(yīng)水庫和慢響應(yīng)水庫；Rf和Rs為Sf和Ss的補(bǔ)給；Qf為非冰川區(qū)快速響應(yīng)徑流；Qs為非冰川區(qū)慢響應(yīng)徑流；Qg為冰川區(qū)徑流；Q為總徑流。

1.3.3 模型的率定和檢驗(yàn)

選擇2013年和2014年作為率定期，2015年和2016年作為檢驗(yàn)期。采用MOSCEM-UA算法來進(jìn)行參數(shù)優(yōu)化［52］，模型參數(shù)個(gè)數(shù)為17（表1），迭代次數(shù)20 000。

采用三個(gè)目標(biāo)函數(shù)：①Kling-Gupta效率系數(shù)［簡稱KGE，式（15）］檢驗(yàn)徑流深的模擬［53］，KGE越接近1表示二者之間的一致性越高，為符合求最小值的優(yōu)化算法，與1相減作為徑流深模擬的目標(biāo)函數(shù)L1［式（16）］；②計(jì)算δ18O的均方誤差（L2），檢驗(yàn)δ18O的模擬［式（17）］；③各高程帶冰川物質(zhì)平衡的體積偏差效率（L3）［18］，檢驗(yàn)冰川物質(zhì)平衡的模擬［式（18）］。均方誤差和體積偏差效率越接近0表示模擬和實(shí)測之間的偏離程度越小。為了從產(chǎn)生的帕累托前沿識(shí)別出有效參數(shù)集，定義可接受閾值的性能指標(biāo)。L1可接受閾值定義為0.20。L2和L3的可接受閾值分別定義為0.56和0.30。只有低于定義閾值的參數(shù)集才會(huì)被判別為有效參數(shù)集。

式中：γ為模擬和實(shí)際徑流深線性相關(guān)系數(shù)；α為模擬和實(shí)際徑流深標(biāo)準(zhǔn)差比值；β為模擬和實(shí)際徑流深均值比值；Sobs（k）和Ssim（k）分別為流域在第k天實(shí)測和模擬的δ18O；m為數(shù)據(jù)總數(shù)；Mlobs（t）和Mlsim（t）分別為冰川在t年海拔帶l觀測和模擬的冰川物質(zhì)平衡；NB為海拔帶總數(shù)；N為總年數(shù)。

為探討不同率定數(shù)據(jù)對(duì)徑流深、同位素、冰川物質(zhì)平衡等過程模擬檢驗(yàn)的影響，設(shè)計(jì)了4種方案（表2）。所有方案關(guān)于同位素過程的參數(shù)，固定為方案2有效參數(shù)集的一組解。

表2 實(shí)驗(yàn)方案設(shè)計(jì)Table 2 Experimental scheme design

2 結(jié)果與分析

2.1 模型率定檢驗(yàn)及不確定性分析

如表3所示，在率定期，方案1徑流深模擬最好（KGE中位數(shù)為0.893），而對(duì)冰川物質(zhì)平衡和δ18O的模擬表現(xiàn)較差。方案2對(duì)δ18O的模擬最好（均方誤差中位數(shù)為0.545‰），在徑流深和冰川物質(zhì)平衡的模擬中略差于方案1。方案3在冰川物質(zhì)平衡的模擬最好（體積偏差效率中位數(shù)為0.235），徑流深模擬效果下降，δ18O的模擬最差（均方誤差中位數(shù)為0.545‰）。方案4除了在徑流深模擬方面（KGE中位數(shù)為0.821）不如其他方案，冰川物質(zhì)平衡（體積偏差效率中位數(shù)為0.265）和δ18O（均方誤差中位數(shù)為0.555‰）的模擬上均有不錯(cuò)的表現(xiàn)。在檢驗(yàn)期，模擬效果往往隨著率定數(shù)據(jù)的增加而提高。方案1在徑流深和δ18O的模擬上表現(xiàn)較差，但冰川物質(zhì)平衡的模擬效果最好（體積偏差效率中位數(shù)為0.324）。方案2在δ18O和冰川物質(zhì)平衡的模擬中表現(xiàn)較好，而在徑流深的模擬上（KGE中位數(shù)為0.654）最差。方案3徑流深的模擬效果最好（KGE中位數(shù)為0.818），而δ18O和冰川物質(zhì)平衡的模擬均不如其他方案。方案4在各方面的模擬效果均較好（KGE中位數(shù)0.813；均方誤差中位數(shù)為0.608‰；體積偏差效率中位數(shù)為0.369）。對(duì)于4種方案有效參數(shù)集產(chǎn)生的模擬范圍（圖2和圖3），在率定期，方案1在三個(gè)方面產(chǎn)生的不確定性范圍較大，而方案2和方案3在各個(gè)方面模擬的不確定性范圍均減小方案4相比于方案3產(chǎn)生更小的不確定性范圍，說明了多種數(shù)據(jù)集的率定具有減少不確定性的優(yōu)勢，以及額外加入水同位素?cái)?shù)據(jù)對(duì)限制不確定性有重要作用。在檢驗(yàn)期，結(jié)果與率定期一致。

圖2 在率定期4種方案有效參數(shù)集模擬的日徑流深、徑流δ18O、各海拔高度帶冰川年物質(zhì)平衡的不確定性范圍Fig.2 Uncertainty envelopes of daily runoff depth，runoff δ18O，and annual glacier mass balance in different elevation bands，produced by the valid parameter sets of the four schemes in calibration period

圖3 在檢驗(yàn)期4種方案有效參數(shù)集模擬的日徑流深、徑流δ18O、各海拔高度帶冰川年物質(zhì)平衡的不確定性范圍Fig.3 Uncertainty envelopes of daily runoff depth，runoff δ18O，and annual glacier mass balance in different elevation bands，produced by the valid parameter sets of the four schemes in evaluation period

表3 由4種方案有效參數(shù)集生成的徑流深、δ18O和冰川物質(zhì)平衡的評(píng)價(jià)指標(biāo)Table 3 The evaluation index of discharge，δ18O and glacier material balance generated by the valid parameter set of the four schemes

對(duì)比4種方案，發(fā)現(xiàn)將水穩(wěn)定同位素用于模型的率定可以提升模擬效果。方案1和方案2的比較表明，將徑流δ18O用于模型率定，徑流和冰川物質(zhì)平衡的整體模擬性能均有提高，不確定性范圍均有所減小，顯著提高了真實(shí)性和穩(wěn)定性。這表明在沒有地面測量的冰川物質(zhì)平衡數(shù)據(jù)的冰川流域水文模擬中，水同位素?cái)?shù)據(jù)能提高徑流模擬效果和減小不確定性。方案3的模擬結(jié)果也表明冰川物質(zhì)平衡數(shù)據(jù)也有助于提高徑流δ18O的模擬。這進(jìn)一步表明，在消融季，徑流水同位素組成能夠代表冰雪消融過程。方案3與方案4對(duì)比，方案4在多源數(shù)據(jù)率定中額外加入徑流δ18O能進(jìn)一步減少模擬不確定性，徑流、冰川物質(zhì)平衡和徑流δ18O都得到了最好的整體模擬效果。這突出了水穩(wěn)定同位素能提供徑流過程內(nèi)部信息，使模擬更加接近真實(shí)。

綜上所述，同位素信息的加入提供了一個(gè)新的驗(yàn)證指標(biāo)，提升了內(nèi)部變量的模擬效果，有助于得到正確的模擬結(jié)果，能為模擬提供有價(jià)值的信息。這與德國科學(xué)家在中亞吉爾吉斯斯坦得到的結(jié)果一致［18］。

2.2 參數(shù)不確定性與辨識(shí)能力

圖4 比較了4種方案得到的有效參數(shù)集的不確定性范圍。方案1產(chǎn)生了4種方案中參數(shù)最大的不確定性范圍，說明數(shù)據(jù)參數(shù)約束不足。方案2相較于方案1，參數(shù)產(chǎn)生了更小的不確定性范圍（除了Ks、Ca），說明了水同位素?cái)?shù)據(jù)約束模型參數(shù)的能力。方案2相較于方案3，大部分參數(shù)（β、Ce、D、Tt、Fdd、Cg、Frr、Kf，g）產(chǎn)生了更小的不確定性范圍，說明了水同位素?cái)?shù)據(jù)相比冰川物質(zhì)平衡數(shù)據(jù)約束參數(shù)的能力更強(qiáng)。方案4相比方案3，產(chǎn)生了更小的參數(shù)不確定性范圍（除了Wm、Ks、Ca、Cwh、Frr），說明了水同位素和物質(zhì)平衡數(shù)據(jù)結(jié)合能進(jìn)一步約束參數(shù)。

圖4 4種方案產(chǎn)生的有效參數(shù)集的不確定性范圍Fig.4 Uncertainty ranges of the valid parameter sets produced by the four schemes

圖5 展示了水文模塊各參數(shù)對(duì)δ18O模擬效果敏感性測試的結(jié)果。在每個(gè)子圖中，每一個(gè)參數(shù)先驗(yàn)范圍內(nèi)展示對(duì)δ18O的模擬效果，同時(shí)其他參數(shù)固定，橫坐標(biāo)為歸一化的參數(shù)先驗(yàn)范圍。由圖5可知，對(duì)δ18O的模擬顯示出對(duì)Tt有很強(qiáng)的敏感性，其次為Fdd、Cwh和Kf，g。徑流δ18O的模擬對(duì)冰川地區(qū)產(chǎn)流過程的參數(shù)（Tt、Kf，g）和雪產(chǎn)流過程融化（Tt、Fdd和Cwh）的參數(shù)顯示出顯著的敏感性，說明了δ18O變化受融雪和融冰的主導(dǎo)。數(shù)據(jù)表明，δ18O數(shù)據(jù)對(duì)模型參數(shù)的有效范圍具有不同程度的約束能力，對(duì)Tt、Fdd、Cwh和Kf，g等與冰雪消融相關(guān)的參數(shù)識(shí)別能力最強(qiáng)。

圖5 徑流δ18O模擬對(duì)各參數(shù)值的敏感性測試Fig.5 Runoff δ18O simulation sensitivity test to each parameter value

2.3 徑流成分的定量分割

2013—2016年1號(hào)冰川斷面各組分對(duì)徑流的貢獻(xiàn)率如下（圖6，模擬范圍來自方案4）：年均32%～34%來自融雪，48%～51%來自融冰，0%～7%來自地下水，12%～15%來自降雨徑流。季節(jié)尺度下（5月至10月為消融季，11月至次年4月為非消融季），融雪在5月至7月占主導(dǎo)（41%～47%），在8月至10月貢獻(xiàn)為23%～27%，在非消融季占主導(dǎo)（50%～51%）。融冰在5月至7月（38%～41%）占比較大，在8月至10月（57%～63%）占主導(dǎo)。降雨徑流的在消融季貢獻(xiàn)為10%～18%，在非消融季貢獻(xiàn)小于5%。地下水貢獻(xiàn)在消融季小于8%，而在非消融季占主導(dǎo)（49%～50%）。在消融季，融冰是第一大徑流成分，而融雪是第二大徑流成分，降雨徑流次之，最后是地下水。在非消融季，融雪和地下水為主導(dǎo)，融冰和降雨徑流貢獻(xiàn)低。

由圖6也可知，4種不同方案之間不確定性范圍有明顯差異。方案1產(chǎn)生了最大的不確定性區(qū)間，以及更高的融冰貢獻(xiàn)和更小的地下水和融雪貢獻(xiàn)。例如，與方案4相比，在非消融期地下水相差最大可達(dá)45%，融雪消融期相差最大可達(dá)13%。很明顯，低約束模型（方案1）試圖通過更高的融冰貢獻(xiàn)，來彌補(bǔ)其他偏低貢獻(xiàn)。另外三種方案均傾向于產(chǎn)生較低融冰貢獻(xiàn)。方案2傾向于得到更高地下水貢獻(xiàn)比例（1%～15%），方案3傾向于得到更高的融雪貢獻(xiàn)比例（30%～38%）。方案4得到了最低融冰貢獻(xiàn)比例、最高的融雪貢獻(xiàn)比例和最小的不確定范圍。

圖6 4種方案年季節(jié)尺度下徑流組分的貢獻(xiàn)率Fig.6 Contributions of runoff components quantified by four schemes at annual and seasonal scales

對(duì)比4種方案，發(fā)現(xiàn)將穩(wěn)定同位素用于模型的率定會(huì)使徑流組分的妥協(xié)情況減少，得到更為準(zhǔn)確的徑流分割。方案1和方案2之間的比較表明，將水同位素?cái)?shù)據(jù)用于模型校準(zhǔn)，融雪和地下水的貢獻(xiàn)比例均有所提升，降水占比變化不大，融冰比例占比下降約7%，這表明同位素信息可顯著減少冰雪相關(guān)過程不確定性。對(duì)比方案3和方案4，在多源數(shù)據(jù)多目標(biāo)率定的基礎(chǔ)上，額外考慮同位素信息后，對(duì)徑流內(nèi)部妥協(xié)進(jìn)一步限制。融雪和地下水的貢獻(xiàn)比例有所上升（約2%和4%），降雨和融冰的貢獻(xiàn)有所下降。這表明同位素信息也可以增強(qiáng)控制地下水產(chǎn)生相關(guān)的參數(shù)（主要為參數(shù)D）的識(shí)別能力，從而改善對(duì)徑流分割模擬的準(zhǔn)確性。由上述兩組方案對(duì)比也可知，在水文模型校準(zhǔn)過程中加入水同位素?cái)?shù)據(jù)還能進(jìn)一步減小徑流成分貢獻(xiàn)率的不確定性范圍。

綜上所述，缺乏內(nèi)部變量驗(yàn)證的模擬結(jié)果（方案1），雖然率定期得到了最佳的徑流模擬效果，但是在檢驗(yàn)期模擬不理想，這說明徑流各組分間存在較大互相妥協(xié)，造成異參同效，影響了模擬真實(shí)性。本研究中，低約束模型試圖調(diào)整參數(shù)得到更高的融冰貢獻(xiàn)比例，來提升率定期的模擬。同位素信息的加入通過提高對(duì)過程和參數(shù)的辨識(shí)能力，使徑流成分之間的互相妥協(xié)明顯降低。這表明水同位素組成具有提供有關(guān)徑流成分內(nèi)部分配信息的能力。這與Birkel等［54］和Capell等［24］在溫帶濕潤區(qū)得到的結(jié)果一致。

3 討論

寒區(qū)示蹤水文模型FLEXG-iso基于寒區(qū)水文過程模型FLEXG作為框架，以基于不同水源的氫氧穩(wěn)定同位素特征差異為前提，不僅要考慮不同端元混合過程，還要考慮積雪消融升華、土壤蒸發(fā)等同位素分餾效應(yīng)?？紤]了水同位素分餾過程，包括土壤水蒸發(fā)和融雪過程伴隨的分餾［18］。對(duì)混合過程也進(jìn)行了詳細(xì)的考慮。有研究發(fā)現(xiàn)在地下水水庫中設(shè)定更大的混合體積可以獲得更好的同位素模擬效果，這說明參與混合的地下水庫蓄量遠(yuǎn)比單純考慮徑流過程設(shè)置的蓄量要大［24，30］，這與本研究結(jié)論一致。另外，考慮到冰川區(qū)表面存在小型冰面湖［34］，并且冰川區(qū)至觀測站之間有部分非冰川區(qū)域，本研究在冰川區(qū)也添加了額外儲(chǔ)量。在該研究流域中，由于冬季植被覆蓋率低，植被冠層對(duì)雪的攔截很小。因此，沒有考慮升華對(duì)雪同位素的影響。本研究結(jié)果表明，所提出的耦合模型能夠重現(xiàn)該流域徑流中同位素時(shí)間序列。本研究的另一個(gè)出發(fā)點(diǎn)在于探討依靠有限的采樣數(shù)據(jù)是否能得到較好的同位素和水文模擬效果。結(jié)果表明，在冰川流域開展水穩(wěn)定同位素的示蹤輔助水文模擬是可行的，這為后續(xù)在類似流域收集更詳細(xì)觀測數(shù)據(jù)，開展更準(zhǔn)確水文模擬工作提供了基礎(chǔ)。

本研究也存在一定局限性。首先是同位素?cái)?shù)據(jù)的代表性。由于數(shù)據(jù)限制，本研究未考慮同位素組成的高度效應(yīng)，采用站點(diǎn)的測量值代表整個(gè)流域，同位素分布的異質(zhì)性可能會(huì)帶來不確定性。未來研究需加密觀測網(wǎng)絡(luò)，基于實(shí)測得到同位素海拔梯度變化以估算高海拔地區(qū)同位素組成的空間異質(zhì)性，進(jìn)一步改善模擬效果［18，55］。本研究假定降雪和降雨的δ18O相同，并且基于閾值溫度來劃分降雪和降雨，未考慮大氣中的分餾過程和缺乏降雨和降雪的驗(yàn)證，分量取值可能存在偏差。地下水的同位素?cái)?shù)據(jù)由于變化幅度小，采樣的時(shí)間分辨率較粗，可能會(huì)遺漏一些短期的變化波動(dòng)信息。另外，模擬本身存在一定的不確定性。徑流同位素組成的變化取決于多方面因素的影響，例如春季融雪的速率和持續(xù)時(shí)間等［56］，在日尺度模擬中很難充分考慮這些細(xì)節(jié)。完全混合和瞬時(shí)混合的假設(shè)也值得商榷，更復(fù)雜的內(nèi)在機(jī)制仍需進(jìn)一步探討［24，51］。例如，夏季冰雪表面的融水與雪坑中冬季形成的細(xì)粒雪和中粒雪之間進(jìn)行了大量的物質(zhì)交換［57］，瞬時(shí)完全混合還不能充分模擬出這些細(xì)節(jié)。

模型對(duì)冰川物質(zhì)平衡損失的模擬存在低估，可能由于未考慮該區(qū)域的粉塵等導(dǎo)致冰川反射率下降或者冰川覆蓋地區(qū)風(fēng)吹雪造成的積雪重新分配［58-59］。而且，在該模型中未將雪冰川相關(guān)參數(shù)的時(shí)間變化考慮在內(nèi)，這會(huì)導(dǎo)致模擬的偏差［60］。另外，雖然FLEXG-iso在徑流、冰川物質(zhì)平衡和徑流同位素的模型上均得到了較好的結(jié)果，也與Gao等［3］的結(jié)果一致，但融冰貢獻(xiàn)比例仍需進(jìn)一步嚴(yán)格驗(yàn)證。例如，賈玉峰等［14］、Jia等［61］根據(jù)水量平衡模型計(jì)算的26%融冰貢獻(xiàn)比例；Thiel等［11］的研究中2013—2016年冰川融水貢獻(xiàn)率具有較大的年際變化（20%～60%）。耦合遙感等多源數(shù)據(jù)在率定水文模型中的價(jià)值也值得進(jìn)一步探討［62］。

4 結(jié)論

（1）本研究構(gòu)建了水同位素示蹤的冰川流域水文模型FLEXG-iso，并成功應(yīng)用于2013—2016年烏魯木齊河源1號(hào)冰川流域徑流過程模擬，并重現(xiàn)了徑流水穩(wěn)定同位素和冰川物質(zhì)平衡。將水穩(wěn)定同位素?cái)?shù)據(jù)用于模型率定檢驗(yàn)，在率定期和檢驗(yàn)期內(nèi)提升了模擬結(jié)果的可靠性，降低了不確定性。

（2）利用水穩(wěn)定同位素有效約束了模型的參數(shù)范圍，證明利用短期的水同位素?cái)?shù)據(jù)能有效約束控制關(guān)于雪和冰川產(chǎn)流過程的參數(shù)（Tt、Fdd、Cwh和Kf，g），從而提高對(duì)模型參數(shù)的識(shí)別能力，減少因參數(shù)識(shí)別誤差而導(dǎo)致的不確定性。

（3）徑流分割的結(jié)果表明，2013—2016年烏源1號(hào)冰川斷面徑流約32%～34%來自融雪，48%～51%來自融冰，0%～7%來自地下水，12%～15%來自降雨徑流。穩(wěn)定同位素作為輔助信息減少了過多補(bǔ)償融冰貢獻(xiàn)帶來的不確定性。徑流組分之間相互妥協(xié)的情況減少，融冰所占比重相比原模型下降約7%。

通過構(gòu)建穩(wěn)定氫氧同位素示蹤的冰川流域水文模型（FLEXG-iso），為實(shí)驗(yàn)和模型科學(xué)家的溝通建立橋梁，將推動(dòng)寒區(qū)水文學(xué)理論和方法的發(fā)展。在此基礎(chǔ)上也將為未來寒區(qū)水資源預(yù)估、水資源優(yōu)化配置、水災(zāi)害治理，寒區(qū)生態(tài)環(huán)境保護(hù)等提供科學(xué)依據(jù)。