1. 引言
降水是水循环中不可或缺的一部分,是径流形成的物质基础。近几年来,受全球气候变化和副热带高压的共同影响下,我国的降水逐渐呈现出年内时空分配不均、汛期强度大的特点,形成湍急的地表径流 [1] ,由此带来的洪涝灾害也逐渐呈现出普遍且复杂的特征。相对于江南、华南等地区,我国西北地区河流海拔较高,属于典型的高寒山区,常年积雪覆盖,因此其地表径流大部分来自于冰川融雪。近年来,在全球气候变暖的背景下,我国西北地区冰川融化及融雪速度正处于快速上升趋势,冰川融雪正加速流失,给洪水的泛滥造成潜在的威胁。此外,西北地区的地形往往山高坡陡,一旦发生较大洪水,则严重影响着当地居民的正常生活,并且防洪形势往往比江南、华南等地区更为严峻,特别是在中小河流发生的洪水灾害,往往具有发生频率高、形成过程复杂的特征,造成极大的破坏。
研究流域水文状态的变化规律并摸清洪水的形成机理是降低洪涝灾害损失的必要条件之一,而水文模型是研究流域径流特性在不同时空尺度变化规律的重要技术手段。近几十年来,国内外很多学者相继开发了很多水文模型,而分布式模型由于具有能够考虑地形和下垫面分布不均等因素的特征,成为国内外水文研究者的热点课题 [2] [3] 。我国最常用的模型则是由河海大学赵人俊教授研制的新安江模型。它虽然属于集总式模型,但其具有功能齐全并能够客观地反映流域水文过程的每一个环节的特征,在我国湿润及半湿润地区得到了广泛的应用 [4] [5] 。而新安江模型在处理冰川及融雪径流方面相对较弱,而对于有积雪覆盖的流域进行水文预报,积雪消融非常重要,因此在类似我国西北地区具有高纬度、高寒气候等特征的地区并不适用。新安江模型在运用于该地区的水文预报之前,需要对模型进行一定的改进,并考虑融雪径流因素,才能够拓宽其适用性。本文以黑河上游为例,依托流域内实测的降雨径流数据以及气温数据,在传统的新安江模型的基础上,加入融雪径流模块,将流域分为低山丘陵带和高山积雪带,并对高山积雪带的产汇流机制进行了改进,同时对所涉及到的参数进行率定验证,评估其在黑河流域洪水模拟的精度,为西北地区以融雪洪水为主的洪水预报业务提供技术支撑,具有很好的实践意义。
2. 研究区概况及基础数据
2.1. 研究区概况
黑河为我国第二大内陆河,其源头发源于祁连山北麓,莺落峡水文站为黑河流域上游控制站,其控制面积约10,000 km2,位于37˚24'N~39˚06'N,98˚33'~101˚09'E,见图1。研究区总体呈现西高东低的态势,以高山积雪带和低山丘陵带为主,其中高山积雪带约占整个研究区的22%。研究区海拔在1675~5073 m之间,具有明显的大陆性气候特征。研究区降水量相对充足且蒸发少,随高程增加而增加,属于高寒半干旱气候带,平均气温低于4℃,年日照时长达2600 h。研究区水系属于典型融雪补给河流,积雪持续时间最长可达200d左右,在春夏季融雪期,融雪径流对河流的补给量可达到约70%以上,降雨径流量仅占全年总量的19% [6] [7] 。研究区域基本没有人类活动影响,因此其水文数据不做还原即可直接使用,具有天然性和一致性的特征。据莺落峡水文站统计数据显示,黑河上游多年平均水资源量约为16.41亿m3,径流量年内分布较为不均匀,6~8月径流量约占全年的53%。由于研究区山高坡陡,河流汇水速度快导致洪水比较集中,因此融雪洪水灾害频发。
![](//html.hanspub.org/file/11-2411266x8_hanspub.png?20230927101034083)
Figure 1. Water system map of the upper reaches of the Heihe River basin and an overview of the underlying surface
图1. 黑河流域上游流域水系图和下垫面概况
2.2. 基础数据
基础数据主要包括DEM数据、气象数据、水文数据以及下垫面数据。其中,DEM数据来源于地理空间数据云,分辨率为90 m × 90 m;气象数据主要包括流域内气象站点的雨量数据和日平均气温数据,来源于中国气象数据共享网(http://data.cma.cn),时间序列为2010~2019年;水文数据主要包括水文站的实测数据,主要包括流域内札马什克、祁连、莺落峡三个水文站2010~2019年的逐日流量系列,可从水文年鉴中获得;下垫面数据来源于新一代的EOS/MODIS卫星数据,考虑到云量较多,采用MODIS的八日合成积雪覆盖产品MOD10A1,并选用春季融雪季节(3~6月)所涉及到的影响作为遥感数据元。
3. 模型原理
传统的新安江模型分为四个模块,即蒸散发模块、产流模块、分水源模块、汇流模块。而在耦合融雪径流模块中是将传统新安江模型的部分模块进行了改进,也就是将融雪径流模块耦合至新安江模型的部分模块中,并将下垫面分为低山丘陵带和高山积雪带两种类型分别计算,详细介绍如下。
3.1. 蒸散发模块
低山丘陵带的蒸散发仍沿用传统新安江模型采用三层蒸散发模式;高山积雪带可以认为常年温度较低,不考虑蒸散发因素。
3.2. 产流模块
3.2.1. 降水输入处理
雨雪划分是改进新安江模型中非常重要的一部分,主要是通过比较当日的日均气温与融雪临界温度来划分。当日均温度大于临界温度时,该日的所有降水被处理为降雨;当日均温度小于临界温度时,该日的所有降水被处理为降雪。用下式表示:
(1)
式中:Pt为第t日降水,Prain,t为第t日以降雨形式的降水,Psnow,t为第t日以降雪形式的降水,Tt为第t天的日均温度、Tc为雨雪划分的临界温度,一般为0℃。
除此之外,还要判断是降在高山积雪带或者是低山丘陵带,所谓低山丘陵带即土壤植被覆盖的区域,它们的产流机制各自不同,详见3.2.2和3.2.3。
3.2.2. 融雪径流模块
对于融雪径流的计算,度日因子法是应用比较广泛的方法,主要是基于冰川或积雪消融与气温之间的经验关系,因此只需要气温值,则可以构建度日模型,按照如下公式计算 [8] 。其中,SWEM (Snow Water Equivalent Maximum)为最大融雪水当量,即积雪完全融化所形成液态水的总量;a为融雪度日因子(mm/(˚C∙d)),即温度每上升1℃所融化的积雪深度;Tmelt为积雪融化的阈值温度。积雪持水能力决定了融雪产流过程,将在下节详述。
(2)
3.2.3. 耦合融雪径流模块的产流量计算
黑河上游由于处在高海拔地段,其下垫面不仅包含传统意义上的土壤及植被(即低山丘陵带),还包括积雪区(即高山积雪带),不同区域的产流机制也不尽相同。并且,不同的降水形式(降雨或降雪)对不同的下垫面的产流机制也不相同,因此需要进一步区分其产流机制。
1) 低山丘陵带产流
根据研究区域DEM和土地利用分布图显示,下垫面类型若位于海拔相对较低的地区,地表温度则相对高山积雪带较高,当出现降雪则积雪融化速度相对较快。根据研究区多年气象资料显示,在低山丘陵带内,降水主要以降雨为主,降雪事件并不多,且单次降雪事件的降雪量一般较少,因此其产流机制可以假定与蓄满产流的原理一致。
2) 高山积雪带产流 [9]
该下垫面类型主要为海拔相对较高且积雪常年不化的积雪区,其产流机制涉及到积雪最大融雪水当量、实际融雪水当量、积雪存储量、雪被持水能力等因素,并且降雨和降雪在高山积雪带的产流形式是不同的。需要分为降水是以降雨或者降雪的形式分别讨论。
第一种情况:降水的形式为降雨。此时产流的来源一部分来自降雨填满积雪孔隙之后剩余的水,另一部分来自积雪自身的融化。也就是说,降雨先补充积雪孔隙,当蓄水量达到积雪持水能力,剩余的水将在重力的作用下排出。积雪持水能力和土壤持水能力一样,也具有分布的不均匀性,我们用M表示积雪持水能力的不均匀系数。积雪自身的融化也就是实际融雪水当量,通过与前一天的积雪实际融雪当量进行比较得出。这样具体产流过程可用下式表示:
(3)
其中 [8] :
(4)
(5)
式中:Rsnow,t为高山积雪带总产流量;SWE为实际融雪水当量,SN为积雪覆盖量,PE为有效降雨量;SW为积雪含水量;SC为积雪平均持水容量,也称积雪平均持水能力,SCM为其最大值;A为积雪持水容量曲线对应的纵坐标值。这里的融雪产流量源于降水(PE)与积雪消融(SWE)两个部分的贡献,SWE是度日因子和气温的函数,与降水无关。
积雪覆盖量SN的变化可认为与温度的变化有关,具体计算规则如下:
a) 若相邻两天温度均大于融雪临界温度Tmelt,则
(6)
b) 若相邻两天温度均小于融雪临界温度Tmelt,则
(7)
c) 若相邻两天温度的范围包含临界温度Tmelt,则
(8)
第二种情况:降水的形式为降雪。此时,降雪作为积雪存储量储存在雪盖上,若有产流,则全部来源于积雪消融,用公式表示则为:
(9)
第三种情况:无降水。此时产流量就仅来源于融雪水当量,是温度的和前一天融雪水当量的函数,按照式(4)进行融雪径流计算。
3.3. 分水源模块
分水源模块仍采用传统新安江三水源模型的思路,分为地表径流、壤中流和地下径流。而在耦合融雪径流模块的新安江模型中,低山丘陵带和高山积雪带产流机制具有差异性,需要分别计算。
对于低山丘陵带产流,仅考虑降雨径流,产流面积比仅考虑降雨产生的径流量,按照下式计算:
(10)
对于高山积雪带产流,其主要来源是地表积雪融化和以降雪形式的降水,按下式计算。
(11)
其中,PEsnow为降雪形式的降水。相应的,低山丘陵带的地表径流RS则按照下式计算:
(12)
壤中流RI和地下径流RG仍利用出流系数进行计算:
(13)
(14)
类比低山丘陵带的计算方法,高山积雪带的RS、RI、RG分别用下式计算:
(15)
(16)
(17)
3.4. 汇流模块
汇流模块与传统新安江模型原理相同,分为坡地汇流、河网汇流和河道汇流三个阶段。只不过在考虑融雪径流的新安江模型中,是将低山丘陵带和高山积雪带的RS、RI、RG合并起来之后进行计算,即:
(18)
(19)
(20)
在坡地汇流阶段,地表径流QS、壤中流QI和地下径流QG分别按照下式计算:
(21)
(22)
(23)
在河网汇流阶段,利用滞后演算法计算:
(24)
(25)
河道汇流应用马斯京根法进行计算:
(26)
式中:U为单位换算系数,CI、CG、CS分别为壤中流消退系数、地下径流消退系数和地表径流消退系数,L为河网汇流滞时。在马斯京根法中,I为上游流量,Q为下游流量,C0、C1、C2为马斯京根计算系数。
4. 模型的率定与检验
4.1. 需要率定的参数清单
在传统的新安江模型中,需要率定的参数包括流域蓄水容量WM、蒸散发折算系数K、流域蓄水容量分布曲线指数B、深层蒸散发系数C、不透水面积比例IM、自由水容量SM、自由水容量分布曲线指数EX、地下水出流系数KG、壤中流出流系数KI、地下水消退系数CG、壤中流消退系数CG、河网水流消退系数CS、河网水流滞时L、马斯京根法参数K和X。而在考虑融雪径流的新安江模型中,还需要把积水持水能力SC、积雪持水能力的不均匀系数M、雨雪划分临界温度Tc、积雪融化阈值温度Tmelt和度日因子a作为需要率定的参数。
同时,传统新安江模型的初始状态变量中,产流和分水源模块主要包括上层、下层、深层张力水含量WU、WL、WP以及自由水水深S、产流面积比FR;汇流模块主要包括马法演算初始流量QX。而考虑融雪模块的新安江模型中,还需要将积雪储量SN、积雪含水量SW作为模型的初始状态,其余仍与传统新安江模型一致。
4.2. 模型参数范围及率定检验方法
4.2.1. 参数范围
传统新安江模型的参数范围在此处不再赘述,这里指介绍与融雪有关并较为敏感的参数和状态变量。
1) 积雪融化阈值温度Tmelt
可以参考相关经验,在融雪季节开始时取临界气温为4.5℃,在结束时0.75℃,使用时根据流域特点做相应调整,本文对以往学者研究采用的经验值调整,将该参数设置在0.75℃~4.5℃之间。
2) 度日因子a
根据前人的研究及经验 [10] [11] ,度日因子的计算方法有两种,一种是通过野外观测得到,另一种是通过经验公式,具体计算方法为:
(27)
其中,
和
分别指的是雪和水的密度,而水的密度是定值,因此度日因子的值直接取决于雪的密度。据相关学者 [12] [13] 从新疆天山东侧、巩乃斯河谷等地区的观测数据显示,积雪密度基本稳定在0.06~0.24 g/cm3之间,若海拔更高一些,则积雪更稳定,密度相应也就更高,将有0.3~0.35 g/cm3左右。本文参考已有的学术成果,将度日因子取值范围设定在0.6~3.5之间。
3) 积雪储量SN
积雪储量是状态变量,冬春季节温度较低,积雪覆盖一般较厚,取值更高,一般能达到30~50 mm,夏季温度较高,但由于海拔也较高,高山积雪带的积雪一般不会完全融化,一般也会有10~20 mm左右。在模拟的时候,视模拟具体时间设置SN的取值在10~50之间。
4.2.2. 参数率定及检验方法
参数率定采用SCE-UA算法,该方法能有效解决非线性约束的最优化问题,可以快速有效的获取到水文模型参数的全局最优解,在水文模型参数优选中得到了广泛的应用,具体原理可参见文献 [14] 。根据水文预报评估标准,本文采用确定性系数(R2)和纳什效率系数(NSE)对水文模型的模拟效果进行评价。计算公式分别为:
(28)
(29)
5. 结果与讨论
5.1. 参数率定成果
选取2010~2016年3~8月的资料作为模型的率定期,2017~2019年3~8月的资料作为检验期,利用上述方法对流域内札马什克和莺落峡两个水文站及其以上区间的参数进行率定,计算时间步长Δt = 6 h,参数的初始范围和率定结果如表1和表2所示。
![](Images/Table_Tmp.jpg)
Table 1. Calibration results of production and confluence parameters for SMS_3 and LAG_3 of hydrological stations in the basin
表1. 流域内水文站三水源产汇流参数率定成果
![](Images/Table_Tmp.jpg)
Table 2. Parameter calibration results of Muskingum Method at Yingluoxia Station
表2. 莺落峡站马斯京根法参数率定成果
从表1和表2结果可以看出,不管是在传统新安江模型模块,还是融雪径流模型模块,参数率定的结果均在模型的初始范围内。
5.2. 模拟结果及分析
根据所率定的参数结果,用2017~2019年的资料进行模拟验证,同时为了进一步说明考虑融雪模块的新安江模型的模拟效果,本文利用传统的新安江模型,利用相同的资料进行率定验证,以进一步验证本文所用方法的有效性,评价因子成果如表3所示。
![](Images/Table_Tmp.jpg)
Table 3. Calibration and inspection evaluation factor results of hydrological stations
表3. 流域内水文站率定和检验评价因子成果
从两个站点的整体模拟结果来看,耦合融雪的新安江模型要比传统新安江模型的模拟精度要高,表现为同等条件下的R2和NSE值均比传统模型要高。具体来看,札马什克站用传统新安江模型模拟的R2值在0.4左右,NSE值在0.2左右,模拟效果不太理想;莺落峡站的R2虽然达到了0.7左右,但从NSE的值来看,还没有超过0.6,模拟效果基本满意。两个站用耦合融雪的新安江模型来模拟,札马什克的R2接近0.8,NSE也接近了0.6,
![](//html.hanspub.org/file/11-2411266x40_hanspub.png?20230927101034083)
Figure 2. Comparison of flood simulation results using different methods at Zamashk Station in 2019
图2. 札马什克站2019年用不同方法进行洪水模拟结果对比
![](//html.hanspub.org/file/11-2411266x41_hanspub.png?20230927101034083)
Figure 3. Comparison of flood simulation results using different methods at Yingluoxia Station in 2019
图3. 莺落峡站2019年用不同方法进行洪水模拟结果对比
能基本还原天然的径流过程;莺落峡站的R2和NSE均到了0.8左右,模拟效果较好。而不管用什么模型,莺落峡的模拟效果均比札马什克要好,究其原因,可能是由于区间模拟效果本身要劣于河道模拟效果,因为区间模拟所涉及到的参数要多一些,并且融雪径流都来自于区间。莺落峡的流量除了区间之外,还包括上游札马什克和祁连的河道来水,而札马什克本身的流量仅来源于区间,从而导致札马什克的模拟效果要相对差一些。
为了说明洪水过程的模拟效果,以2019年为例,利用传统新安江模型的洪水过程模拟结果与本方法计算出来的结果进行对比,得到札马什克和莺落峡水文站的模拟结果分别如图2和图3所示,同时洪峰和洪量的对比如表4所示。从两个站点的洪水过程来看,利用耦合融雪径流的新安江模型能够较好地模拟出6月14日、6月26日和7月5日的洪水过程。以莺落峡站为例,6月14日的这场洪水的洪峰流量大约在450 m3/s左右,用传统新安江模型模拟出来的洪峰大概在300 m3/s左右,而用耦合融雪的新安江模型模拟的结果能够很好的模拟出洪峰的值(约488 m3/s),说明用传统模型模拟出来的洪峰与实测洪峰的差值基本来自于融雪径流。从洪量结果的模拟分析来看,两个模型的模拟结果均存在一定的误差,但利用改进的新安江模型模拟的洪量结果与实测结果相比,误差明显要小一些;而误差产生的原因可能是由于西部地区普遍属于资料稀缺地区,从而在洪量的计算上可能存在一定的误差,此外该模型还属于改进试验,模型构建方面仍然基于一些假设,因此还不能100%还原洪
![](Images/Table_Tmp.jpg)
Table 4. Comparison of flood peak and flood volume simulation using different methods at Zamashke and Yingluoxia Stations
表4. 札马什克和莺落峡站用不同方法进行洪峰和洪量模拟的对比
水过程的总洪量。从气温数据分析来看,6月9日至6月14日有一轮气温上升的过程,上升幅度约为5℃,从而进一步确定流量的差值是来源于融雪径流。但需要指出的是,模型状态变量SN的设定是很关键的,因为在西北地区的融雪径流的来源大部分都来自于SN。仍以莺落峡站为例,7月5日的洪水没有前两次大,通过实测资料分析,其原因在于对应的降雨量相对较少,但气温相对前两天较高,若SN设置的过大,则7月5日的模拟径流也将会很大,因为模型最初积雪储量较大,到了7月5日这天可能也有一定的积雪储量,温度升高雪会融化,产生融雪径流。而实际通过降雨和水文资料分析,这一场洪水主要来自区间的降雨以及上游来水,并且到了7月温度会升高,下垫面的积雪储量一般不会很多。而当SN设置过小时,则不能模拟出融雪洪水过程,因为此时下垫面已经没有雪被,径流仍以区间降雨为主。
6. 结语
本文在传统新安江模型的基础上,加入了融雪径流模块并将研究区分为低山丘陵带和高山积雪带,分别对其产流机制进行了研究,从而考虑了融雪对径流的影响,在我国西北地区黑河流域上游进行了洪水过程模拟,并利用相关指标因子进行了评价,主要得到以下结论:
1) 传统的新安江模型在札马什克和莺落峡水文站的纳什效率系数仅为0.2和0.5左右,而考虑融雪径流模块的新安江模型进行模拟时,纳什效率系数就达到了0.6和0.8左右,并且在洪水过程的模拟中,能够较好的再现场次洪水过程,并且模拟出来的洪峰流量与实测洪峰相差不大。这表明该模型能够在黑河流域推广应用,并对西北地区以融雪径流为主的洪水预报提供了技术支撑,具有很好的实践意义。
2) 模型在模拟的过程中,由于融雪模块主要体现在区间产流上,而决定区间产流量的因素较多,且融雪径流可能还受到雪水下渗、蒸散发等因素,因此当模型输入仅有区间输入时,一定程度上限制了模型的预报精度,表现在札马什克的模拟效果要稍微次于莺落峡的模拟效果。在今后的研究中,还要从雪被的水量平衡以及产汇流特征等方面做进一步的研究,以进一步提高洪水预报精度。