分布式地貌瞬时单位线在闽江流域的应用探究

董丰成 石 朋1， 陆美霞 鲍庆煜 赵兰兰 陈颖冰 瞿思敏 崔彦萍

(1.河海大学 水文水资源与水利工程科学国家重点实验室， 南京 210098；2.河海大学 水文水资源学院，南京 210098；3.江苏省水文水资源勘测局南京分局， 南京 210008；4.水利部信息中心， 北京 100053；5.江苏省水文水资源勘测局， 南京 210029)

在流域汇流计算中，传统单位线方法已愈发不能满足社会生产作业的需求.传统单位线需通过丰富的水文实测资料分析得到.然而在许多流域，由于水文监测站建设时间短、布设密度不足等原因导致实测降水径流资料缺失，这无疑限制了传统单位线方法的实现和应用.当用系统分析观点解释流域汇流时，单位线就是流域的水文响应[1].一直以来，水文学家试图将流域的水文响应同地形地貌特征建立联系，最早在此方向取得突破性研究成果的是Rodriguez-Iturbe[2]、Valdes[3]等人关于地貌瞬时单位线(GIUH)理论的建立和应用.鉴于不依赖实测径流资料、计算便捷等优点，越来越多的水文工作者开展对GIUH的研究[1]，并将其应用到无资料地区的汇流计算中[4-5].

流域下垫面地形地貌特征复杂多变，流域的汇流特征在空间上存在显著差异[6-7].为此，郭良[8]、Hao[9]等按水系的自然分水岭划分单元，推求出单元上的地貌瞬时单位线以区分流域空间上的汇流特性；Jin[10]基于河流等级划分单元，考虑了汇流性质在空间上的变异性，通过三参数伽马分布地貌瞬时单位线构建了流域的综合地貌瞬时单位线.这些学者在划分单元时均使用了数字高程模型(DEM)生成流域自然水系的结果.值得考虑的是，在利用ArcGIS提取整个流域水系时，所设定的面积阈值对最终结果影响很大.阈值越小，生成的水系越复杂，这也势必会对河流分级结果造成影响[11-12]，其中1级水系的真实性就有待商榷.因此，本文在现有划分单元的基础上，重新制定单元划分方案，进而构建流域分布式地貌瞬时单位线，并将其与蓄满产流模型进行耦合，在闽江建阳流域进行应用和检验.

1 基本方法介绍

1.1 单元流域划分

考虑到径流形成过程的完整性以及面积阈值对生成水系结果的影响，基于流域数字高程模型(DEM)，利用ArcGIS平台自动提取流域水系，并按照汇入最高级河道的河网集水区为基准划分出单元流域.

1.2 地貌瞬时单位线

分布式地貌瞬时单位线的构建采用伽马分布的地貌瞬时单位线[13]，其具体表达形式如下：

式中，v为平均流速(m/s)；Γ(x)为x的伽马函数；RB为霍顿分叉比；RL为霍顿河长比；RA为霍顿面积比；LΩ为最高级河道长度(m)；t为计算时段长.

利用单元流域的DEM数据提取计算地貌瞬时单位线时所需的地形地貌参数(RA、RB、RL和LΩ，一般地，RA取值范围在3.0～6.0，RB取值范围在2.5～5.0，RL取值范围在1.5～3.5)，单元流域的平均流速利用文献[14]中介绍的Eagleson-Bras公式估算得到.

1.3 产流计算

利用河海大学赵人俊教授提出的蓄满产流模型[15]计算单元流域的产流量并划分出不同的径流成分(地表径流、壤中流和地下径流).在蓄满产流模型中，假设在流域上任意一点，土壤含水量在未达到田间持水量前，降水量均用于补给土壤缺水量，而在蓄满后的降水将全部产生径流.其利用了指数分布形式(式4)表达了土壤缺水量在空间分布的异质性：

式中，α为产流面积比；b为反映流域包气带蓄水容量分布特征的参数；W'm为流域各点包气带蓄水量值；Wmm为其中的最大值(mm).

1.4 汇流计算

为了区分不同径流成分的汇流性质，将地表径流直接汇入河网，而壤中流和地下径流经过线性水库(KKSS，KKG)的调蓄后再汇入河网，然后利用流域的地貌瞬时单位线模拟水体从进入河网到单元流域出口的汇集过程.

1.5 河道洪水演算

河道洪水演算利用马斯京根法(式5～7)，将各单元的演算结果进行叠加，即得到整个流域的出流过程.

式中，C0、C1、C2均为流量演算系数；I1、I2分别为计算时段始末的河段入流量(m3/s)；Q1、Q2分别为计算时段始末的河段出流量(m3/s)；k为蓄量常数；x为流量比重因数；Δt为计算时段长.

2 应用实例

2.1 流域概况

研究流域选取了位于闽江上游的建阳流域(如图1所示)，该流域毗邻中国的东南沿海，属于典型的亚热带季风气候，温和潮湿，属于高雨区，多年平均降水量在1 800～2 200 mm之间，雨季多发生在每年的4～9月.流域面积4 848 km2，海拔在45～2 154 m范围之间，是典型的山丘区地形.河床普遍呈蜿蜒曲折的峡谷相间形态，滩多流急，多大石顽礁，水位暴涨暴落，具有流程短等特点.流域内有岭阳、黄坑等17个雨量站，建阳水文站为流域控制站，具有1988～1999年间丰富的日径流资料及洪水相关资料.

图1 研究区域位置及站网示意图

2.2 计算结果及分析

利用建阳流域的30 m×30 m的SRTM DEM数据，以50 km2的面积阈值生成自然水系，划分出14个单元流域(如图2所示)，单元面积在60～673 km2之间.利用划分出各单元的DEM，通过River Tools工具提取各单元的地形地貌参数(见表1)并计算出地貌瞬时单位线.

图2 流域分区示意图

表1 单元流域部分地貌参数及平均流速数值

如图3所示，单元流域的地貌瞬时单位线间存在明显的差异，这也反映了流域汇流特性在空间上的异质性.其中②单元流域形状及水系分布为扇形结构，同样情况下形成的洪峰流量相对较大，并且洪水历时相对较短；④和⑦单元流域是受流域上、中游主河道控制的单元流域，因坡降较大的原因也易形成短历时洪水.

对建阳水文站1988～1999年的25场洪水采用上述的方法进行模拟计算，其模拟特征值统计见表2.

图3 单元流域的地貌瞬时单位线

表2 洪水模拟结果统计

洪峰和径流深精度的评定采用百分比误差形式，绝对值越小表示模拟精度越高.根据计算结果，径流深相对误差在-29.4%至15%之间，误差绝对值小于20%的有23场，均值为8.7%.模拟的洪峰流量与实测值比较接近，相对误差较小，其值分布在-4.7%至10.3%范围内，平均值为3.5%，并且这25场洪水的计算误差均在许可误差的范围内.从图4上和表2中可以看出，相对于径流深的模拟效果，洪峰流量的模拟效果表现更优，洪峰流量的计算值普遍比实测值大，而径流深的模拟值相对于实测值偏小，说明在应用分布式地貌瞬时单位线模拟洪水时，计算的洪水过程涨、退形势相对迅急.洪峰滞时代表了模拟洪峰相对于实测值的出现时间，“+”表示模拟滞后于实测值，而“-”表示提前.对于大部分场次洪水实测与计算洪峰的出现时间相差不大，其中的20场次相差在2 h以内，计算洪峰平均滞后实测约为1.1 h.确定性系数表示模拟和实测洪水过程的拟合优度，数值越接近1表明模拟的结果越接近实际洪水过程.根据计算结果，确定性系数最大值可以达到0.962，最小值为0.537，其平均值为0.844；超过0.85共有14场，超过0.8有20场，占总数的80%.

图4 洪水模拟结果统计箱线图

19880228和19880620场次洪水的径流深相对误差超过了许可误差范围.对比图4中的洪水过程，19880228场次洪水的主峰模拟结果较好；在主峰形成前出现了次峰，但其计算与实测的峰现时差较大，可能由于同时期降水主要分布在中上游地区，而靠近流域出口的计算单元内的降水量很少或者没有.19880620场次在主峰消退后又出现了次峰，但计算与实测的洪峰量级差别较大；观察次峰附近的洪水过程可知，其实测过程涨洪迅猛，很可能是水利工程泄蓄作业造成的影响.19920704场次洪水滞后时间较长，但确定性系数仍高于平均值，为0.945；由图4显示，此场次洪水过程中出现了两个量级差别不大的洪峰；在第二次涨洪阶段的前期降水比较丰富，但峰现前的洪水过程增长趋势突然减缓，这可能是当季的农业或其他需求引水导致的结果.

图5 部分场次洪水过程模拟对比(Q c和Q o分别表示计算和实测流量)

3 结 论

本文基于现有的单元流域划分方法，以汇入最高级河道的河网集水区为基准重新制定划分单元方案；提取了各单元内的地形地貌参数，从而构建分布式地貌瞬时单位线.将其与蓄满产流模型进行耦合，在闽江建阳流域进行应用检验.由各单元流域内地貌瞬时单位线的趋势表明，单元间的汇流特性差异显著.根据洪水模拟计算的结果，模型精度较高，确定性系数平均值达到0.844，并且洪峰模拟的表现最好.据该流域应用效果显示，此方法能够较好地为无资料或资料匮乏地区的洪水预报提供借鉴.

地貌瞬时单位线发展到今天，其理论正确性已得到了证明，但其推广应用仍存在一些难点.譬如：地形资料的精度对于提取到的地形地貌参数具有重要影响，进而影响最终得到的GIUH[12]，在今后的研究中，需要进一步加强对这个问题的探究.