冰川变化的定量评估

冰川变化的定量评估

冰川是自然界的重要淡水资源。世界上3.4%以上的淡水资源存储在冷水（主要是北极冰盖和格陵冰盖）中。中国是世界上中、低纬度冰川最为发育的国家。据第一次冰川编目数据,中国有冰川46377条,总面积59406km2,冰川储量5590km3,相当于长江入海口年径流量的5倍,仅次于美国、俄罗斯、加拿大。随着全球气候变暖,世界各地冰川不断萎缩,冰川变化对干旱区水资源变化、海平面升降以及全球水文循环等的影响越来越明显。冰川水文主要研究冰川与气候相互作用、热量和物质收支、消融与径流过程、冰川对河流的补给、冰川洪水、冰川泥石流、冰面湖以及冰川对气候变化的响应规律。冰川是天然的固体水库,中国每年平均冰川融水量约为620亿m3,与黄河多年平均入海径流量相当。在西北内陆干旱区,冰川融水的重要性尤其突出,塔里木河各源流区冰川融水补给比例多在30%～80%。实际上正是由于冰川和积雪的存在,才使得中国深居内陆腹地的干旱区形成了许多人类赖以生存的绿洲,也使得中国干旱区有别于世界上其它地带性干旱区。这种冰川积雪-绿洲景观及其相关的水文和生态系统稳定和持续存在的核心是冰雪,没有冰雪就没有绿洲。冰川融水同时也是西欧和北美一些国家水利发电的重要资源;此外,冰川又是自然灾害的发源地,如突发性冰川洪水、冰川泥石流及雪崩等自然灾害都会对当地和中下游的人民生命和财产造成严重威胁。通过定量模拟冰川区产流(消融)和汇流过程,可为冰川水资源变化和管理、冰川洪水预防和冰川水能利用等提供理论依据,此外,对冰川内水流运动的理解,对冰川内部压力的分布、冰川滑动和跃动及冰湖溃决等冰川动力学研究至关重要。本文从冰川水文的主要过程——产流(消融)和汇流两方面对冰川水文模型发展进行汇总和分析,以期为冰冻圈动态过程及其对气候、生态、水资源和环境影响的综合研究,特别是冰川水文过程研究,提供水文模型方面的基础信息。1-文化条件:传统冰面结论,包括5.冰川区由于气温较低,当降水发生时,冰川表面并不一定会直接产流,而当气温较高时,即使没有降水发生,冰川表面也会出现产流现象。即冰川消融与降水关系不密切,而主要取决于冰川表面的气象条件如气温、风速、相对湿度、太阳辐射等和冰川本身物理特征如粒雪粒径、粒雪的新旧程度、表碛特征以及坡度、坡向等。北半球冰川消融时间一般为6～8月,而南半球冰川则在12～2月。传统的冰川消融模型就是构建冰面消融与气象因子(如气温)之间的相关关系,或从冰面消融的物理机制角度出发,研究冰面消融和辐射平衡之间的关系。随着3S技术发展和计算机技术的提高,基于能量平衡基础的分布式水文物理模型不断得到改进,成为冰川消融模型发展的必然趋势。因此,冰川消融模型大体上可以分为两类:即基于气象因子的统计模型和基于物理机制的能量平衡模型。1.1冰土温度对冰酒冰流速的影响影响冰面消融的气象因素主要有气温、风速、相对湿度、降水以及地形等,且不同性质和规模的冰川由于其下垫面和气候特征条件的不同,其热力学特性也是不同的。由于冰川区一般缺少或没有常规气象特别是辐射平衡资料,而气温是反映辐射平衡、湍流热交换等状况的综合指标,且资料相对容易获得,所以基于气温的冰川消融模型应用最广,且结构简单、易于推广。相关研究表明,冰川或积雪的消融与气温之间存在着密切的关系。基于气温的统计模型主要有冰川平衡线法和度日因子法。(1)平衡线上的年监狱警察网冰川平衡线又称零平衡线,在此高度上冰川的年累积量等于年消融量。1924年Ahlmann首先提出平衡线高度处的消融量与气温之间存在着密切的关系,Khodakov和Krenke则首先用经验公式计算了平衡线高度处的年消融量。Kotlyakov也把一个类似的经验公式应用到阿尔卑斯山。模型的基本形式为A=a(Τ+b)c(1)A=a(T+b)c(1)式中A为平衡线上消融或积累水当量,mmw.e/a;T为消融期平衡线上的平均气温,℃;a、b、c为经验参数。为公式单位表示一致性,本文将消融、积累、降雪量和雪深的单位统一折算为消融单位,即mmw.e.或mw.e.刘潮海等根据上述原理,提出了中国天山平衡线上的冰川年消融量经验公式:A=0.78(Τ+9.0)3.09(2)A=0.78(T+9.0)3.09(2)(2)日气形势模型度日因子模型也是基于冰川消融量与气温关系的统计模型,它与冰川平衡线法的主要区别在于该模型是基于正积温的,而冰川平衡线法则是基于月平均气温的。模型基本形式如下:Μ=DDF×ΡDD(3)M=DDF×PDD(3)式中M为某一时段冰川或雪的消融水当量,mmw.e.;DDF为冰川或雪的度日因子,mm/(℃·d);PDD为某一时段的正积温,℃。度日因子是Finsterwalder和Schunk在1887年对阿尔卑斯山冰川的研究中提出的,在随后的一个世纪内,度日因子模型被广泛应用到冰岛、挪威等地的冰川研究中。最近,Hock为了提高模型的时空精度,进一步将模型改进为Μ=(DDF+a?ΙGsG0)Τ(4)M=(DDF+a?IGsG0)T(4)式中a为辐射系数;I为晴天太阳直接辐射,W/m2;Gs为太阳总辐射,W/m2;G0为晴天太阳总辐射,W/m2;T为气温,℃。最近,张勇等分析了天山南坡科其卡尔巴契冰川度日因子的空间变化规律。基于气象因子的统计,模型由于结构简单,输入参数较少且易于获取,模拟结果较为理想等特点在中国目前被广泛应用到冰川径流估算、冰川物质平衡对气候变化规律等研究中。但由于统计模型只是在一定数据区域内对水文物理意义的统计分析,无法精确表征冰川消融的实际物理过程,由此造成模型不易空间推广,且其外延预报结果缺乏可信度,特别是在全球变暖的背景下。1.2冰的能量平衡模型随着典型冰川区各种自动水文、气象仪器的布设、消融和辐射平衡观测数据的积累以及冰川消融理论的发展,具有物理机制的能量平衡模型成为研究冰川消融的另一个选择。能量平衡模型形式众多,从其发展历程来看,主要分为单点、基于海拔高度梯度以及分布式能量平衡模型三类。(1)能量平衡计算方法QΜc=Q*+QΗ+QE+QG+QΡ(5)QMc=Q∗+QH+QE+QG+QP(5)式中QMc为融化热,W/m2;Q*为净辐射,W/m2;QH为感热,W/m2;QE为潜热,W/m2;QG为地热通量,W/m2;QP为高温降水产生的热通量,W/m2。Munro在其模型中还提出了不同类型冰川能量通量的简化计算方法,在实际计算中忽略了地热和降水热量,简化了感热和潜热通量的估算方法,这些方法均不同于非冰川区的常规算法,如增加了大气稳定度校正系数等,同时将模拟结果和电子自动消融仪观测值进行比较,发现两者的结果是很接近的。此外,还指出由于受地形、风等的作用,积雪分布不均匀,从而导致粒雪密度不均,使电子自动消融仪观测值并不十分可靠,而用能量平衡计算的冰雪消融量是比较可信的。在国内,康尔泗构建了一个基于常规气象观测要素的参数化能量平衡模型,对中国乌鲁木齐河源1号冰川的消融过程进行了模拟,计算和讨论了消融期冰川表面辐射平衡和能量平衡,为研究高山冰川和周边大气之间能水交换过程,冰川对气候变化的响应以及冰川融水径流的形成提供了一个便于实际应用并具有明确物理意义的计算方法。(2)平衡方程的结构张寅生等将一个基于平衡线高程对气候响应的物理气候学模型应用到我国大陆型山地冰川,也指出能量平衡方程是随平衡线高度变化而变化的。其模型结构如下:QR(Ζ)+QS(Ζ)+QL(Ζ)+Qm(Ζ)=0(6)QR(Z)+QS(Z)+QL(Z)+Qm(Z)=0(6)式中QR为净辐射,W/m2;QS为感热,W/m2;QL为潜热,W/m2;Qm为消融热,W/m2;Z为高程,m。(3)最大波长辐射通量计算其实,分布式模型的理论雏形在早期的冰川产流(消融)模型中就有所体现。1982年Munro就分析了地形对短波辐射的影响,在能量平衡方程中引入地形参数:cosi=cosΖsinh+sinΖcoshcos(A-A1)(7)cosi=cosZsinh+sinZcoshcos(A−A1)(7)式中i为太阳射线和坡面法线之间的夹角,rad;Z为坡度,rad;h为太阳高度角,rad;A为太阳方位角,rad;A1为坡向,rad。白重瑗也从坡度、坡向和周围山体等地形因子对辐射量的影响出发,对我国天山乌鲁木齐河源冰川与空冰斗辐射平衡进行过比较计算。Arnold将一个分布式模型应用到瑞士的一个山谷冰川流域中,系统地介绍了地形对净短波辐射和反射率的影响。模型输入参数有流域DEM、太阳高度和方位角、积雪深度初始分布面积图以及常规的气象观测数据等。模型把单元格分为受太阳直射和被地形阴影遮盖两类,分别计算单元格净短波辐射通量。计算方法如下:(1)太阳直射单元格净短波辐射通量Q*=(1-α)L/ρw×Q1sinh×cosi(8)Q∗=(1−α)L/ρw×Q1sinh×cosi(8)式中α为反射率;L为融化热,J/kg;ρw为水密度,kg/m3;Q1为实测短波辐射值,W/m2。(2)遮蔽单元格净短波辐射通量:Q*=(1-α)L/ρw×[0.2Q1cos2(Ζ/2)+αmQ1sin2(Ζ/2)](9)Q∗=(1−α)L/ρw×[0.2Q1cos2(Z/2)+αmQ1sin2(Z/2)](9)式中αm为流域平均反射率。Arnold同时还根据Oerlemans的思路对反射率的时空变化进行了模拟,认为反射率是海拔、雪深、表面融化、冻结作用和新雪等综合作用的结果,并分别针对上述变量提出了相应修正公式:反射率随海拔高度的变化:αb=a1arctan[(h-E+300)/200]+a2(10)αb=a1arctan[(h−E+300)/200]+a2(10)式中αb为随海拔高度变化的反射率;h为高程,m;E为平衡线高度,m;a1、a2为经验参数。雪深、表面融化及再冻结作用对反射率的作用:αs=αos-(αos-αb)ea3d-a4Μ(11)αs=αos−(αos−αb)ea3d−a4M(11)式中αs为受雪深、表面融化及再冻结作用影响的反射率;αos为旧雪的标准反射率;M为累积融水当量,mw.e.;d为雪深,mw.e.;a3、a4为经验参数。新雪对反射率的作用:αs(n)=αns-(αns-αs)e-5?000p(12)αs(n)=αns−(αns−αs)e−5?000p(12)式中αs(n)为受新雪作用影响的反射率;αns为新雪的标准反射率;p为新雪深度,mw.e.。分布式能量平衡模型虽然从物理机制上详细地揭示冰川的消融过程,但由于输入参数较多,理论和结构比较复杂,加上冰川监测存在一定的困难,因此,分布式能量平衡模型在中国实际应用中还存在一定困难,目前我国还没有比较成熟的分布式能量平衡模型。最近,Konya把一个基于经验关系的分布式辐射和气温指数模型(RTModel)应用到俄罗斯的Koryto冰川,模型结构简单,所需参数少,仅需气温和总辐射作为模型输入参数。其结构如下:Μ=αR+βΤ+γ(13)M=αR+βT+γ(13)式中R为总辐射,W/m2;T为气温,℃;α、β、γ为回归系数。R和T由整个冰川的DEM插值所得,α、β、γ由观测点的M、R和T进行多元回归确定。2冰内基因的汇流形成计算在冰川消融期,融水和液态降水在冰面的水道网中流动(该过程类似于普通的坡面汇流过程),对中大型冰川尤其是海洋性冰川来说,其中大部分的融水沿冰川表面的裂隙和垂直通道(冰内竖井)进入冰川内部,少部分沿冰面两侧和末端直接流出冰川体,而进入冰川内的水又汇流成一条或几条巨大水流从冰川末端下面流出来。该过程是一个水、热以及冰川运动相互作用的过程,冰内水道系统在三者的作用下,其管道形态、大小和长度都是动态变化的。从20世纪60、70年代人们就开始用河流流量记录、示踪物和冰内凿洞等方法对冰川内水道系统进行研究。Fountain在前人的研究基础上,系统总结了冰川的汇流方式,从粒雪的水文特性、冰川内部和下部的排泄方式、冰川水储量、冰湖溃决及冰川动力学等方面对冰川汇流进行了详细描述,同时,他认为冰内管道系统对时间和输水量变化的响应研究是很重要的。但由于目前对冰川内部实际排泄方式,水热交换及冰川运动对冰川内部管道系统的作用,管道系统对外界环境变化的响应以及对不同规模、类型和地域冰川的汇流方式等缺乏详细的认识,基于物理过程来量化冰川汇流还存在很大的困难。本文从冰川表面汇流、内部汇流和冰川下部汇流3个方面,对现有冰川水文汇流的研究成果进行概述,旨在为今后冰川融水汇流模型构建提供一些借鉴。2.1裂隙或垂直通道冰川表面径流的汇流方向是由冰川表面的地形坡度决定的,按最陡坡度原则即D8(Deterministiceight-neighbors)原则,确定各单元的水流方向和路径长度,由此获得单元格到裂隙或垂直通道的时间和距离。但由于裂隙或垂直通道位置的测量是一件十分危险和耗时的事,Flowers等在对加拿大Trapridge冰川的研究中发现,裂隙或垂直通道的位置可以通过分析冰川表面曲率并利用卫星遥感影像进行校正获得,即冰川凸起部位一般可以作为冰川表面张力的指示,而在张力、太阳辐射和水流等共同作用下,这些部位可能沿粒雪纹理演变成裂隙并最终可能发展为垂直通道。积雪覆盖与否对冰川表面汇流也有很大影响,由于冰的相对隔水性,当表面有积雪时,水在重力的作用下先从表面沿晶体间缝隙向积雪底部运动,在冰面上形成一个饱和水层,然后随地形坡度向下流动,所以积雪对冰川表面汇流有滞后效应。(1)饱和水流运动速率水从积雪表面到达底部的时间可用Colbeck提出的公式计算。D=ked(3ρwgμ-1)1/3k1/3q2/3(14)D=ked(3ρwgμ−1)1/3k1/3q2/3(14)式中D为水流到达积雪底部的时间,s;ke为雪层有效孔隙度;d为雪深,m;ρw为水的密度,kg/m3;g为重力加速度,m/s2;μ为水的粘滞系数,Pa·s;k为雪的渗透比(量纲一);q为水通量,m3/s。积雪层下部饱和水流运动速率用如下公式表示:CS=ρwgμ-1k?θksnow(15)CS=ρwgμ−1k?θksnow(15)式中Cs为饱和积雪层下部的水流速率,m/s;θ为坡度,rad;ksnow为雪的孔隙度。(2)水流速度、水力半径及曼宁模拟Arnold采用曼宁公式计算每个单元格冰川融水的流速:Ci=R2/3θ1/2/n(16)Ci=R2/3θ1/2/n(16)式中Ci为水流速度,m/s;R为冰川表面水道的水力半径,m;n为曼宁糙率系数,m-1/3·s。由上所述,根据每个单元格水流方向可以得到该单元格到下一个单元格的距离,除以式(15)(或式(16))得到的水流速度,即可以得到该单元格到下一个单元格的汇流时间,最后根据每个单元格到裂隙或垂直通道的汇流路径,从而计算每个单元格到达裂隙或垂直通道的总汇流时间。2.2冰内管道方程冰川融水进入冰川内部的主要通道是裂隙或垂直通道,但冰川内部的管道一般很细小,且在任何特定时间内,均有两种相反的作用决定管道的大小和长度,即Shreve提出的管道理论:①水在管道中流动会融化管壁上的冰而使管道扩大;②若上覆冰的压力超过水压,则冰流进入管道并减少管道的直径。该理论对冰床处的冰内管道是同样适用的。Shreve认为水在总势能(重力势和冰的压力势)的作用下会沿最陡的方向形成一个向上的树枝状的排泄结构,除坡度大于表面坡度一个量级的管道内的情况外,水流方向基本由表面坡度控制。Nye推导了冰川管道(假设管道中充满了水)中非稳定状态水流的微分方程。Flowers等为了简化冰川内部的汇流过程,把冰川内部管道系统理想化为3种结构:垂直通道(即冰内竖井,假设可到达冰川的底部)、表面裂隙通道和底部裂隙通道。冰川内部的水流连续方程可以用下面方程描述:∂he∂t+∂Qej∂xj=ϕr‚e-ϕe‚s(17)∂he∂t+∂Qej∂xj=ϕr‚e−ϕe‚s(17)式中he为冰川内水的平均体积,m3;Qejej为单宽流量,m3/s;ϕr,e为冰川内部和冰川表面水流交换量,m3/s;ϕe,s为冰川内部和冰川下部水流交换量,m3/s;t为时间,s;x为方向差;j为水流方向。2.3冰的研究过程冰川下部一般存在管道式(channelized,快速)和分布式(distributed,慢速)两种排泄系统。管道式排泄系统一般由R管道组成,而分布式排泄系统则由冰下空穴、水膜、N水道、冰下沉积物等组成,且对于不同的冰川其分布式排泄系统的组分也是不同的。在消融期初期,分布式排泄系统随融水输入量的增加可转变为管道式排泄系统,而在消融期末期,随着冰川运动及来水量的减少,管道式排泄系统则会转变成分布式排泄系统。现代冰川下部汇流过程研究可以追溯到20世纪60年代,Mathews对加拿大SouthLeduc冰川的研究结果。SouthLeduc冰川下部水压观测结果表明,夏季水压要明显低于冬季水压,而当大量融水或大降水过程发生时,水压会急剧增加,这反映了冰川表面和冰川下部存在着明显的水力联系。随后Hodge第一次用钻孔技术对冰川下部汇流进行了系统的研究。目前,冰川下部汇流模型主要有黑箱模型、概念性集总模型和分布式模型3种。(1)钻孔间的水压响应关系Murray和Clarke把一个黑箱模型应用于加拿大Trapridge冰川,用来研究与冰川排泄系统相联系的钻孔和无联系的钻孔之间的水压响应关系。研究发现钻孔中水压存在着3种明显不同的水力联系:显著性相关、显著性差异和介于二者之间的关系。Murray和Clarke指出这三种关系可用一个低阶的非线性差分方程来很好的描述。随后Corn等用计算智能方法(遗传算法、神经网络、模糊数学、自组织映射等)改进了Murray和Clarke所提出的模型算法。(2)冰的排泄系统的结构分析方法Clarke认为冰川下部的排泄系统可以用一种类似电路分析的方法进行概念性集总模拟,模型对冰川下部排泄系统R管道水力特性、排泄系统形态的转变、冰川储水量对汇流的非线性作用以及管道系统的季节性变化都能做出很好的解释。Fountain等认为这一方法在研究冰川下部排泄系统的结构特征的方法中有很大的发展潜力。模型把冰川内部水压差和流量输入输出比作电路中的电压降、管道系统当作电阻、冰川内部的储水结构看作电容、管道形态变化认为是电源开关,各种排泄方式可由电路的串连和并联组合所得。(3)模型形式的建立Arnold把美国环保局研制的雨洪管理模型SWMM(stormwatermanagementmodel)中的管网汇流子系统引入分布式冰川模型中,裂隙或垂直通道(冰内竖井)代表SWMM中的进水或排水口,分布式的和渠道式的两种排泄系统的组分代表SWMM中的管道,连接进水口和排泄口。雨洪管理模型SWMM中的管网汇流子系统水流运动过程采用圣维南方程组求解。方程形式为∂A∂t+∂Q∂x=q(18)vg×∂v∂x+1g×∂v∂t+∂h∂x+Sf-S0=0(19)∂A∂t+∂Q∂x=q(18)vg×∂v∂x+1g×∂v∂t+∂h∂x+Sf−S0=0(19)式中A为管道面积,m2;Q为流量,m3/s;q为管道侧向单宽流量,m3/s;t为时间,s;x为长度,m;h为管内水深,m;g为重力加速度,m/s2;v为流速,m/s;Sf为摩阻坡度;S0为底坡;Sf、S0均为量纲一。上式中A是一个随融水输入量Q和冰川流变而变化的量,Arnold根据管道内水流速度和冰川下部管道中的水压Pw及上覆冰的压力Pi以小时为步长对模型中管道面积及其长度进行调节。管道融化速率如下:W′=[(πA)1/2ρw(fr?v3/4)]×1L(20)式中W′为管道的融化速率,kg/(m·s);A为管道面积,m2;ρw为水密度,kg/m3;fr为摩擦系数;v为流速,m/s;L为水的融化潜热,J/kg。管道长度的变化速率如下:Μ′=-(pi-pw)|pi-pw|m-1×2(1mB)mA(21)式中M′为由于冰川蠕变造成的管道面积的变化速率,m2/s;pi为上覆冰的压力,Pa;pw为管道中的水压,Pa;m为冰的格伦定律中的指数;B为冰的格伦定律中的阿列纽斯参数,N·m-2·s1/m。如前所述,冰川的下部排泄系统可由不同的元素组成,从而模型的形式也会随之改变。在加拿大Trapridge冰川下面存在一个面积很广,厚度很薄的透水层,Flowers等在前人研究的基础上,把一个分布式的与冰川表面、内部和地下水耦合模型应用到Trapridge冰川。模型假设这个透水层是饱和的,且孔隙率随水量的变化而变化,水流运动方式用达西定律描述为∂hs∂t+∂Qsj∂xj=bs-ϕr‚s+ϕe‚s-ϕs‚a(22)式中hs为透水层的平均水体积,m3;Qsj为流量,m3/s;ϕr,s为冰川下部与冰川边缘的交换量,m3/s;ϕe,s为冰川下部与冰川内部的交换量,m3/s;ϕs,a为冰川下部与地下水的交换量,m3/s;bs为地热和冰川滑动所造成的融水量,m3/s;t为时间,s;x为方向差;j为水流方向。3第三,文化研究和模型研究有待进一步加强随着全球变暖,冰川不断萎缩,冰川水文物理模型已经成为国内外冰川学和全球变化研究领域的热点之一。近百年的发展过程中,冰川水文模型得到了较大的发展,特别是基于物理机制的分布式冰川水文模型的应用,为真实揭示冰川径流的形成,冰川物质波动对气候的响应,冰川消融加剧和萎缩,突发性冰川洪水等机理的研究提供了有力的工具。但由于冰川产汇流过程的复杂性以及流域尺度空间变异性,使得在描述和模拟冰川径流形成方面与实际物理过程相差较远。其中冰川产流(消融)模型发展已相对较为成熟,基于气象因子的统计消融模型日趋成熟,但其无法应对全球不断变化的气候;分布式消融物理模型能够高精度模拟和预测全球变化