基于改进新安江模型的岩溶地区径流过程模拟——以庙沟岩溶流域为例

李玉坤，燕子琪，王纪元，陈 丽，周 宏

(中国地质大学(武汉)地质调查研究院，湖北 武汉 430074)

流域水文模型是在流域产汇流结构概化的基础上，以一个实体或数学模型来模拟流域降雨-径流的形成过程，是一种常用的径流过程模拟方法。流域水文模型通常分为集总式水文模型如新安江模型、SAC模型和水箱模型，分布式水文模型如SWAT模型、SHE模型，半分布式水文模型如TOPMODEL模型等。

岩溶地区含水介质内部结构复杂多样，实测资料获取难度较大，使得流域水文模型在该地区的应用成为难点。分布式水文模型基于水动力学机制来描述和模拟流域内的水文过程，但对资料要求较高，因而在缺乏资料地区的应用受到限制；黑箱模型又称经验模型，应用门槛低、适用性强，但模型不能反映含水系统内部的水文过程，也无法用于评价和管理人类活动对泉流量的影响[1]。新安江模型是一种概念模型，由于采用张力水和自由水蓄水容量曲线来描述张力水和自由水蓄水容量分布的空间异质性，因此已成为一种具有分布意义的概念性流域水文模型。近年来，随着GIS和RS技术以及计算机的发展，新安江模型参数的精度和计算分区的精度已得到显著提升，如李致家等[2]将新安江模型应用于栅格之上，与高精度雷达降雨资料耦合进行了洪水预报；井立阳等[3]根据地形因子分布曲线推求了自由水蓄水容量;Yuan等[4]应用双源蒸散发模型改进了新安江模型结构，模拟了在植被参与的水循环过程中的径流过程；Zhao等[5]将新安江模型应用于包括土壤雨水侵蚀、土地资源利用等土壤资源评价中。目前新安江模型已由功能单一的洪水预报模型发展为适用于多领域、多功能的流域水文模型。

在岩溶地区，新安江模型的应用也得到了发展。如庄一鸰等[6]在总结南方岩溶地区水文特性的基础上，对新安江模型的稳定下渗、水源划分和汇流模块进行了改进，并取得了较好的模拟效果；郝庆庆等[7]在新安江模型原有自由水库的基础上，再增加一个线性水库用来模拟深层地下水的出流过程，使模拟精度得到了显著提高。但目前新安江模型在岩溶地区的应用研究仍较为有限。为此，本文选取庙沟岩溶流域为典型岩溶流域，以降雨入渗为切入点，对新安江模型进行了改进，以期提高新安江模型对岩溶地区径流过程模拟的精度。

1 研究区概况

庙沟岩溶流域位于湖北省宜昌市兴山县黄粮镇高阳镇境内,为高岚河流域的子流域，整体上处于中国大陆二级阶地向三级阶地过渡地带。该流域面积约为25 km2，平均海拔为1 100 m，地貌类型为溶蚀侵蚀中山沟谷地貌。该流域西侧为斜坡，岩性以志留系粉砂岩为主；东侧为岩溶洼地地貌，岩性以寒武-奥陶系碳酸盐岩为主，区内落水洞、岩溶漏斗等岩溶形态十分发育,见图1。

图1 庙沟岩溶流域水系概况图Fig.1 Drainage map of Miaogou karst basin

2 研究方法

2.1 新安江模型及预模拟

新安江模型以蓄满产流为理论基础，该模型结构为分散性结构。如图2所示，降雨信号输入模型后转化为张力水，此时一部分水量耗于蒸散发，采用三层蒸发模型进行蒸散发计算，另一部分水量转化为自由水，并将水源划分为地表径流、壤中流和地下径流三部分；最后采用马斯京根法进行河道汇流计算。

图2 新安江模型结构简图Fig.2 Structure diagram of Xin’anjiang model

为了验证原始新安江模型对岩溶流域径流过程的模拟效果，采用原始新安江模型对庙沟岩溶流域进行了一次预模拟，选取的时间序列为2014年3月2日到2014年9月2日，期间共有4次小洪峰，1次大洪峰。原始新安江模型预模拟采用遗传算法对模型参数进行率定，预模拟参数率定值见表1，预模拟的效果图见图3。

表1 原始新安江模型预模拟参数率定值[8]Table 1 Pre-simulation parameters of the original Xin’anjiang model[8]

图3 原始新安江模型预模拟的效果图Fig.3 Pre-simulation results of the original Xin’anjiang model

从原始新安江模型预模拟的结果可以发现以下现象：

(1) 根据贾云飞等[9]对金沙江流域张力水蓄水容量的研究，一般流域平均张力水蓄水容量WM取值范围在90～150 mm之间，而预模拟的WM值为31.88 mm，与经验值相比偏小。

(2) 模拟流量对降雨的响应过于敏感，多大于实际流量，即使无效降雨，模拟值也会有所响应。

(3) 模拟洪峰流量较为滞后，起峰过程长，退水过程短，峰形较宽缓，拖尾现象不明显。

结合庙沟岩溶流域的水文地质条件，对上述现象分析如下：

非岩溶地区降雨降至地面之后，经过缓慢入渗至地下水，这一过程在原始新安江模型中体现为三层蒸散发模型和表层、下层、深层3个串联张力水水库的计算，因此该模型对降雨信号的过滤能力强，较好地模拟了非岩溶地区含水介质对降雨的调蓄功能。而在岩溶地区，降雨经过地表产流通过岩溶管道直接汇入地下，相对于非岩溶区，岩溶区降雨在地表和浅层土壤滞留的时间更短，汇入地下转化为地下水的速度更快，效率更高，所以岩溶地区地下水对降雨的响应更为迅速，退水也更为迅速。因此，将新安江模型直接应用于岩溶地区径流过程的模拟会造成以下两个结果：①加长了洪峰的响应时间，导致模拟洪峰流量较为滞后；②减小了流域平均张力水蓄水容量WM值和自由水蓄水容量SM值，减小了产流阈值，缩短了模拟值的响应时间。这将导致模拟结果对降雨的响应过于敏感、拖尾现象不明显。

2.2 新安江模型的改进

Hartmann等[10]认为，流域的岩溶发育程度越高，渗透性能越好。而岩溶地区降雨主要通过洼地、落水洞等岩溶地貌直接灌入并补给地下水，因此等效于降雨中的一部分受包气带调蓄，参与蒸发过程，另一部分直接补给地下水(见图4)。故本文根据新安江模型预模拟的结论，将直接补给地下水的一部分降雨量定义为P1，总降雨量为P2，则参与蒸发过程的降雨量为P2-P1，并定义了降雨分配系数φ：

(1)

郭绪磊[11]、罗明明[12]的研究表明，南方岩溶地区地下水降雨入渗系数多在0.07～0.7之间，故本文中φ的取值范围为0～0.7。

图4 降雨-自由水转换结构示意图Fig.4 Schematic diagram of rainfall-free water transition

另外，在地下径流(RG)出流之前，增加一个线性水库用以模拟深层溶洞水(快速径流)与裂隙水(慢速径流)的水量交换过程(见图5)，进入快速流水库和慢速流水库的水量由岩溶管道和裂隙的发育情况决定[7]。

图5 溶洞快速和裂隙慢速径流结构示意图Fig.5 Schematic diagram of fast runoff in karst caves and slow runoff in fracture

溶洞快速径流流量和裂隙慢速径流流量的计算公式如下:

(2)

式中：Q出，1、Q出，2分别代表快速径流和慢速径流的出库流量(m3/s)；Q入，1、Q入，2分别代表快速径流和慢速径流的入库流量(m3/s)；α1、α2分别为溶洞水和裂隙水的地下水衰退系数。

2.3 模型参数的敏感性分析

原始新安江模型改进后共涉及17个参数，其中φ、α1和α2为新增参数。本文首先采用One-factor-At-a-Time(OAT)法对模型参数进行局部敏感性分析。参数的敏感系数Si表达式为

(3)

式中:x为参数变化量;y为因变量模拟流量；Si为参数的敏感系数即敏感度，其大小代表敏感性的强弱。

改进的新安江模型各参数取值同表1，模型参数的敏感度Si计算结果及其排序见表2[13]，模型参数敏感性分类见表3。

表2 改进的新安江模型参数敏感度及其排序Table 2 Sensitivity of parameters of the improved Xin’anjiang model and the ranking

注：敏感度较低的模型参数未列出。

表3 改进的新安江模型参数敏感性分类Table 3 Classification of the sensitivity of parameters of the improved Xin’anjiang model

结合表2和表3可知，改进的新安江模型参数SM、WM、KG和φ为极敏感参数，KC、KI和CI为敏感参数。其中，张力水蓄水容量WM和自由水蓄水容量SM均属于模型中线型水库的出流阈值，WM控制产流，SM控制自由水水量，影响水源划分，两者直接影响着模型中水量的转换，参数的敏感性较高；地下水出流系数KG控制着自由水S中地下水的出流比例，对于Nash-Sutcliffe效率系数DC有着十分显著的影响[14]，因而参数的敏感性较高；降雨分配系数φ描述了通过岩溶管道直接灌入地下转化为自由水的部分降雨，实质上该参数减少了水分蒸发量，并将减少的水分直接转化为自由水，减小了模型对降雨的调蓄作用，使得模拟径流对降雨的响应更加迅速，而这一点恰恰与岩溶地区的产流特点相吻合。

2.4 φ值对岩流流域径流过程模拟的影响分析

为了进一步探究参数φ值对岩溶流域径流过程模拟效果的影响，本文针对庙沟岩溶流域140831号洪峰流量，利用改进的新安江模型对其进行了模拟，分析不同φ值对岩溶流域径流过程模拟效果的影响，并与参数SM值和WM值进行比较，其模拟结果见图6。其中图6(a)、(b)、(c)分别为只改变φ、SM、WM参数时，庙沟岩溶流域洪峰流量模拟曲线与实际流量曲线的对比图。

由图6(b)、(c)可见，随着参数SM、WM的增大，模拟洪峰流量整体上逐渐减小，强降雨和弱降雨的模拟流量的“平均流量水平”逐渐降低。

表4为不同敏感参数取值时模拟得到的庙沟岩溶流域洪水信息统计表。

由表4可知，当WM值逐渐增大时，平均流量、洪峰流量和确定性系数DC值逐渐减小，流量误差逐渐增大；当SM逐渐增大时，平均流量、洪峰流量和DC值逐渐增大，流量误差逐渐减小，但当SM取值大于80 mm时，DC值和流量误差分别保持在96.4%和19.78%；在参数φ值逐渐增大的过程中，各项指标均符合先增大后减小的规律，并且当参数φ取值为0.5时，岩溶流域径流过程的模拟效果达到最优。

图6 不同敏感参数取值时庙沟岩溶流域洪峰流量模拟曲线与实际流量曲线的对比图Fig.6 Comparison between simulated and actual curves of runoff process in Miaogou karst basin with different sensitive parameter values

表4 不同敏感参数取值时模拟得到的庙沟岩溶流域洪水信息统计表Table 4 Statistical table of simulated flood information of Miaogou karst basin with different sensitivity parameter values

由此可见，WM值控制模拟流量的“平均流量水平”，影响模拟流量的纵向“高度”，决定着模型模拟效果的整体水平；SM为自由水蓄水容量，向下由地下水日出流系数KG控制地下水径流，壤中流出流系数RI控制壤中流，当自由水量大于SM值时，出流成为地表径流RS，而岩溶地区地表径流成分较小，因此当自由水量大于SM值时，SM值的大小不会改变模型的模拟效果，DC值和流量误差保持稳定。

结合图6(a)可以发现，当参数φ值较小时，洪水峰值更多落在图中Ⅰ区，较小洪峰流量的模拟效果较好，较大洪峰流量的模拟效果较差，改进新安江模型对强降雨的响应效率较低，模拟效果更接近于原始新安江模型；随着参数φ值的增大，弱降雨模拟洪峰流量逐渐减小，强降雨模拟洪峰流量逐渐增大，洪水峰值逐渐由Ⅰ区向Ⅱ区过渡；当参数φ值较大时，洪水峰值主要落在Ⅱ区，较小洪峰流量的模拟效果较差，较大洪峰流量的模拟效果较好，说明改进新安江模型对强降雨的响应效率更高。整体上模拟值的“平均流量水平”变化幅度相对较小。表4中数据也表明，参数φ值在整个取值域中流量误差的变化范围在3%以内，而SM值和WM值在整个取值范围内误差的最大变化范围分别为7.5%和48.5%。

3 模型应用与分析

3.1 模型参数率定

南方地区蒸散发折算系数KC一般取0.6，上层初始含水量WU、下层初始含水量WL和深层初始含水量WD是比较重要的指标，其值的大小可以看作是整个流域内蓄水量变化的起点，直接影响后期模型参数率定的准确度以及模型模拟的效果[15]。

流域平均张力水蓄水容量WM可用实测资料分析，即选择前期特别干旱的一次全流域产流过程进行分析，根据庙沟岩溶流域实际降雨资料，确定WM取值为50 mm。利用2014年3月2日——2014年9月15日降雨径流数据作为改进的新安江模型参数率定期，采用遗传算法对模型参数进行率定，最后确定了Ⅰ区和Ⅱ区改进新安江模型的参数，见表5。

表5 改进的新安江模型参数率定值Table 5 Parameter values of the improved Xin’anjiang model

3.2 模拟结果与分析

根据参数率定表确定的改进新安江模型参数，对Ⅰ区和Ⅱ区改进的新安江模型分别赋值并计算2014年3月1日—2014年12月1日的降雨径流，其模拟结果见图7,模拟结果的精度见表6。

图7 模拟径流量与实际径流量的关系图Fig.7 Diagram of the relationship between simulated flow and actual flow

表6 新安江模型模拟精度评定表Table 6 Simulation accuracy evaluation table for Xin’anjiang model

由图7可见，改进的新安江模型的模拟效果较好，所模拟的洪峰起峰时间及峰值均与实际洪峰过程较为吻合。由表6可知，改进后新安江模型模拟结果的确定性系数DC由0.592增长至0.901，平均绝对误差和相对误差均有较大幅度的减小。总体上看，改进的新安江模型模拟值偏小，从而导致水量平衡偏小。

为了探究改进的新安江模型对庙沟岩溶流域洪峰过程的模拟效果，利用改进的新安江模型针对庙沟岩溶流域2014年4月10日、2014年7月14日、2014年8月13日以及2014年8月31日4次小洪峰过程进行模拟，并分析模型参数对改进新安江模型模拟效果的影响，其模拟结果见图8和表7。

图8 庙沟岩溶流域4次洪峰径流过程模拟结果与实际流量的对比图Fig.8 Comparison between the simulated and actual runoff during the runoff process of four flood peak periods in Miaogou karst basin

表7 庙沟岩溶流域4次洪峰径流过程模拟值与实际值信息表Table 7 Information table of simulated runoff and actual runoff in the runoff process of four flood peak periods

注：a为原始新安江模型；b为改进的新安江模型。

由图8和表7可见，改进新安江模型的确定性系数明显优于原始新安江模型，140410号洪峰可能是由于降雨-流量的相关关系较差，故而其模拟的效果较差。

4 结 论

(1) 原始新安江模型直接应用于岩溶地区径流过程的模拟效果较差，本文根据岩溶地区降水灌入并补给地下水的特点，对原始新安江模型入渗模块结构进行了调整，增加了降雨-自由水直接转换环节，改进后新安江模型模拟结果的确定性系数由0.592增为0.901，单次径流过程模拟精度最高可达0.993。

(2) 降雨分配系数φ实质上将部分降雨直接转化为自由水，减少了包气带对降雨的调蓄时间，缩短了模拟径流对降雨的响应时间，提高了新安江模型降雨-径流的转化效率。

(3) 降雨分配系数φ的敏感度为1.991，在新安江模型众多参数中排在第4位，仅次于SM、WM和KG参数，且该参数更侧重于改变确定性系数DC值，对模拟流量的影响较小。

综上所述，改进的新安江模型对庙沟岩溶流域径流过程的模拟效果有显著提升。