产流误差的动态系统响应曲线修正方法应用

阙家骏 包为民，2

（1.河海大学 水文水资源学院，南京 210098；2.河海大学 水文水资源与水利工程科学国家重点实验室，南京 210098）

洪水预报是防洪非工程措施中的一项重要内容，准确及时的洪水预报对防洪减灾、水库调度兴利、工农业安全生产等意义重大［1］.在实际流域中，由于存在许多误差因素，如设备故障，导致资料缺测或不合理的观测数据；水利工程、农田蓄放水误差；流域水文规律的变化；水文规律简化误差，即模型结构误差，故需利用实时系统能获得的观测信息和一切能利用的其他信息对预报误差进行实时修正，以弥补流域水文模型的不足［2］.洪水预报实时修正方法有很多种［3－6］，如AR模型、动态系统响应曲线、卡尔曼滤波、递推最小二乘估计等.AR模型结构简单、计算简便，但前提是流量误差的前后时段具有较强相关性，否则修正效果不甚理想.卡尔曼滤波较为复杂，可容纳的信息量大，但实时洪水预报常难以提供修正所需的足够信息，这就大大限制了方法的使用，导致其修正效果与AR模型效果类似.本文根据东圳流域的实际情况建立概念性水文模型，分别采用DSRC和二阶AR模型对预报结果进行实时修正，并将修正结果与预报结果进行比较分析.

1 流域概况

东圳水库位于木兰溪支流延寿溪中游，是一座以灌溉为主，结合防洪、发电、供水、养殖等综合利用的多年调节的大型水库，流域面积321km2，流域内植被较好，平均年降水量1 740mm，多年平均入库水量3.02亿m3，水库供应莆田市的生活用水和工农业用水，特别是湄洲湾北岸工业发展用水，是莆田市的生命之源［7］.

流域属亚热带海洋性季风气候，10月～次年3月为旱季，4～9月为雨季，雨量充沛，累积雨量约占全年雨量的80%.其中4月下旬～6月下旬为梅雨季，常出现持续连绵的阴雨，洪水特点是历时长、洪量大、但洪峰一般较小，累积总雨量约占全年降雨量的1／3，称前汛；7～9月为台风季节，台风夹带降雨一般历时短，雨量集中、强度大、洪水特点是洪峰较大，称后汛.洪水陡涨陡落，汇流较快，洪水过程一般为3d左右.

2 东圳流域概念性水文模型

根据东圳流域的地理位置、降水、气候、下垫面条件和洪水过程线分析，属于典型的湿润地区，拟选用三水源新安江模型.模型分为蒸散发、产流、分水源和汇流计算4个层次结构，其中蒸散发计算为3层蒸散发模式；产流计算为蓄满产流模式；分水源计算为自由水蓄水库结构；坡面汇流计算为单位线或线性水库；河道汇流计算为马斯京根法.模型将流域分成若干个单元流域，对划分好的每块单元流域分别进行蒸散发、产流、水源划分和汇流计算，得到单元流域出口的流量过程；对单元流域出口的流量过程进行出口以下的河道汇流计算，得到该单元流域在全流域出口的流量过程；将每块单元流域在全流域出口的流量过程线性叠加，即为全流域出口总的流量过程，模型结构如图1所示［2］.图中下标i为单元站点；P为降水量；IMP为各单元流域不透水面积比；Ew为蒸散发能力；EU、EL、ED分别为上层、下层、深层蒸散发量；E为实际蒸散发量；WU、WL、WD分别为上层、下层、深层张力水蓄量；R为蓄满产流计算的总径流量；RS、RI、RG分别为地表径流、壤中流、地下径流；I和Q分别为马斯京根河道演算的入流和出流.

图1 三水源新安江模型结构图

3 实时修正方法

3.1 产流误差的动态系统响应曲线修正方法

产流误差的动态系统响应曲线修正方法是将计算值与实测值之间的误差（流量误差）归咎于产流误差，将图1中蓄满产流之后的部分概化为一个系统，假定流域出口断面的流量过程变化是由于产流量R的变化所引起的系统响应，即

式中，Q（t）为流域出口的总径流过程；R＝［r1，r2，…，rm］T为产流量系列；θ为模型参数；t为时间.式（1）的微分表达式为

式中，Rc＝［rc1，rc2，…，rcm］T为计算得到的产流量初始系列值.

假设样本系列长度为m，则实测流量为Q（t）＝［Q1，Q2，…，Qm］T；计算流量为Q（Rc，θ，t）＝［Q（Rc，θ，1），Q（Rc，θ，2），…，Q（Rc，θ，m）］T，由式（2）可得：

式中，ΔR＝［Δr1，Δr2，…，Δrn］T为产流量误差；E＝［ξ1，ξ2，…，ξm］T为经实时修正后的预报系统残差，一般为白噪声向量；U为动态系统响应曲线矩阵，其表达式为

式（4）U中的每一项都可以用下式进行差分求解：

对式（4）进行最小二乘估计可得：

预报结果的修正式为

式中，为修正后产流量.

3.2 AR模型

AR模型假定各时段的流量误差存在相关性，可用P阶的线性关系来进行模拟.对于误差系列｛ε1，ε2，…，εt，…，εt＋L，…｝，可建立式（8）表达式.

预报结果的修正式为

假设样本系列长度为m，实测流量为Q（t）＝［Q1，Q2，…，Qm］T；计算流量为Qc（t）＝［Qc1，Qc2，…，Qcm］T，可得模型误差系列｛ε1，ε2，…，εm｝，由式（8）可得矩阵形式为

式中：

对式（10）进行最小二乘估计可得：

4 应用检验

在将实时修正方法应用于实际流域后，修正效果用如下3项指标来比较衡量.

1）洪峰相对误差

2）次洪径流深相对误差：

3）确定性系数：

本文选用东圳流域1965～1999年共16场历史洪水的降雨、蒸发和流量资料，利用三水源新安江模型进行计算，模型参数见表1.所有参数在计算和修正过程中保持不变.

表1 东圳流域使用的三水源新安江模型参数

其中KC为蒸散发折算系数；WM为平均张力水容量；UM、LM分别为上层、下层张力水容量；B为张力水蓄水容量曲线方次；C为深层蒸散发折算系数；SM为表层自由水蓄水容量；EX为表层自由水蓄水容量曲线方次；KI、KG分别为表层自由水蓄水库对壤中流、地下水的日出流系数；CS、CI、CG分别为地表水、壤中流、地下水消退系数；KE、XE为马斯京根法演算参数［2］.

表2 新安江模型计算结果、AR模型和DSRC实时修正结果

表3 新安江模型计算结果、AR模型和DSRC实时修正结果比较

5 结 论

准确实时的洪水预报是东圳流域防洪安全的重要保证，实时修正模型能够有效提高概念性水文模型的预报精度.本文分别应用DSRC和AR模型对东圳流域的预报结果进行实时修正，分析结果表明，相比于AR模型，DSRC对于洪峰的修正效果更佳，对径流深的修正效果也较好，这对于水库流域防洪抗灾的意义重大，值得进一步地研究和推广应用.

［1］ 刘国纬.论防洪减灾非工程措施的定义与分类［J］.水科学进展，2003，14（1）：98－103.

［2］ 包为民.水文预报［M］.北京：中国水利水电出版社，2006.

［3］ 包为民，王 浩，赵 超，等.AR模型参数的抗差估计研究［J］.河海大学学报：自然科学版，2006，34（3）：258－261.

［4］ 司 伟，包为民，瞿思敏.洪水预报产流误差的动态系统响应曲线修正方法［J］.水科学进展，2013，24（4）：497－503.

［5］ 郭 磊，赵英林.基于误差自回归的洪水实时预报校正算法的研究［J］.水电能源科学，2002，20（3）：25－27.

［6］ 瞿思敏，包为民.实时洪水预报综合修正方法初探［J］.水科学进展，2003，14（2）：167－171.

［7］ 刘金英.东圳水库生态环境与水资源保护的探讨［J］.水利科技，2003，2（2）：14－15.