CN113435630B - 一种产流模式自适应的流域水文预报方法及系统 - Google Patents

一种产流模式自适应的流域水文预报方法及系统 Download PDF

Info

Publication number
CN113435630B
CN113435630B CN202110621866.4A CN202110621866A CN113435630B CN 113435630 B CN113435630 B CN 113435630B CN 202110621866 A CN202110621866 A CN 202110621866A CN 113435630 B CN113435630 B CN 113435630B
Authority
CN
China
Prior art keywords
basin
sub
runoff
target
watershed
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202110621866.4A
Other languages
English (en)
Other versions
CN113435630A (zh
Inventor
陈璐
罗东
江平
张寒松
易彬
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Meizhou Hydrological Branch Of Guangdong Provincial Bureau Of Hydrology
Shantou Hydrological Branch Of Guangdong Hydrological Bureau
Huazhong University of Science and Technology
Original Assignee
Meizhou Hydrological Branch Of Guangdong Provincial Bureau Of Hydrology
Shantou Hydrological Branch Of Guangdong Hydrological Bureau
Huazhong University of Science and Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Meizhou Hydrological Branch Of Guangdong Provincial Bureau Of Hydrology, Shantou Hydrological Branch Of Guangdong Hydrological Bureau, Huazhong University of Science and Technology filed Critical Meizhou Hydrological Branch Of Guangdong Provincial Bureau Of Hydrology
Priority to CN202110621866.4A priority Critical patent/CN113435630B/zh
Publication of CN113435630A publication Critical patent/CN113435630A/zh
Application granted granted Critical
Publication of CN113435630B publication Critical patent/CN113435630B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q10/00—Administration; Management
    • G06Q10/04—Forecasting or optimisation specially adapted for administrative or management purposes, e.g. linear programming or "cutting stock problem"
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q10/00—Administration; Management
    • G06Q10/06—Resources, workflows, human or project management; Enterprise or organisation planning; Enterprise or organisation modelling
    • G06Q10/063—Operations research, analysis or management
    • G06Q10/0639—Performance analysis of employees; Performance analysis of enterprise or organisation operations
    • G06Q10/06393—Score-carding, benchmarking or key performance indicator [KPI] analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q50/00—Systems or methods specially adapted for specific business sectors, e.g. utilities or tourism
    • G06Q50/10—Services
    • G06Q50/26—Government or public services
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A10/00—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE at coastal zones; at river basins
    • Y02A10/40—Controlling or monitoring, e.g. of flood or hurricane; Forecasting, e.g. risk assessment or mapping

Abstract

本发明公开了一种产流模式自适应的流域水文预报方法及系统，属于水文学中的水文预报领域。本发明针对流域不同前期土壤含水量特征及降水过程特征，构造了基于蓄满产流、超渗产流和混合产流于一体的流域产流模式判别方法；进一步地，依据无因次单位线原理，根据流域历史洪水过程，推导出了基于小流域地形、地貌特征的无因次综合汇流单位线；并以洪峰目标、峰现时间目标、确定性目标作为评价指标，对各子流域确定的最优产流模式涉及的参数进行优化率定，进而形成完整预报体系。本发明能够针对流域不同土壤含水量特征和降雨特征，确定流域最优下渗模式，并根据流域地形地貌特征，确定流域汇流单位线，具有计算简单、适用范围广、拟合精度高的优点。

Description

一种产流模式自适应的流域水文预报方法及系统
技术领域
本发明属于水文学中的水文预报领域，更具体地，涉及一种产流模式自适应的流域水文预报方法及系统。
背景技术
我国水文预报工作在过去的几十年里得到很大的发展，积累了丰富的经验，大量洪水预报方法已取得显著成效，但由于我国地域辽阔，地貌类型复杂多样，许多中小流域仍存在预报不够精准的情况。究其原因，主要是由于中小流域下垫面条件复杂，降雨变化多端，导致其产流机理复杂多样，传统的单一产流理论根本无法解析其下渗过程，而产流模拟是流域径流模拟的第一步，若产流模拟效果不好，则整个径流过程基本上无法达到满意的精度。
由此可见，现有水文预报技术存在产流方法单一、流域下渗机理无法解析的技术问题，进而，对于中小流域的径流预报，无法给出高精度的预报。
发明内容
针对现有技术的以上缺陷或改进需求，本发明提供了一种产流模式自适应的流域水文预报方法及系统，其目的在于解决水文预报技术存在产流方法单一、流域下渗机理无法解析的技术问题，进而，对于中小流域的径流预报，无法给出高精度的预报。
为实现上述目的，按照本发明的一个方面，提供了一种产流模式自适应的流域水文预报方法，包括如下步骤：
(1)获取目标流域洪水预报相关资料，并将所述目标流域划分为多个子流域；
(2)基于各子流域前期土壤含水量因子和降雨过程因子与相应阈值的大小关系，确定各子流域最优产流模式并拟合相应的产流过程；
(3)将所述产流过程代入对应的汇流单位线，并模拟所述目标流域的初始洪水过程；
(4)以洪峰目标、峰现时间目标、确定性目标作为评价指标，对各子流域确定的最优产流模式涉及的参数进行优化率定；将参数优化后拟合的产流过程代入对应的汇流单位线，得到所述目标流域的最终洪水过程。
进一步地，所述步骤(1)中，
所述目标流域洪水预报相关资料包括：雨量站逐小时降雨数据、蒸发数据、各站点经纬度数据、前期M天日降雨量数据、以及流量站对应时间点的连续逐小时径流数据；
所述将所述目标流域划分为多个子流域包括：利用泰森多边形法将所述目标流域划分为多个子流域。
进一步地，所述步骤(2)具体包括如下子步骤：
(2-1)计算各子流域前期土壤含水量：
Figure BDA0003100148380000021
式中，Wi,j为第i个子流域在第j日对应的场次洪水发生初始阶段的土壤含水量；Pj为j日的日降雨量；Pa,j-1为j-1日前期影响雨量；n为影响本次径流的前期降雨天数；k为土壤消退系数；当Pa,j＞Wm时，Pa,j＝Wm，其中，Wm为流域平均最大蓄水量；
(2-2)计算各子流域x小时和y小时连续最大降雨量：
Figure BDA0003100148380000022
式中，HPi,x为第i个子流域连续x小时最大降雨量；HPi,y为第i个子流域连续y小时最大降雨量；Pi,k为第i个子流域第k个小时的降雨量；n为连续最大降雨起始时刻；
(2-3)预先设置各子流域前期土壤含水量阈值α1、x小时连续最大降雨量阈值α2和y小时连续最大降雨量阈值α3，并确定各子流域最优产流模式：
a.当Wi,j＜α1、HPi,6＜α2且HPi,12＜α3，或者Wi,j＜α1、HPi,6＜α2且HPi,12＞α3，或者Wi,j＞α1、HPi,6＞α2且HPi,12＜α3时，流域发生混合产流，产流公式如下：
ft＝a-bWt-cTt
式中，Tt为洪水t时刻从上涨后开始计数的时段数，初始值默认为1，之后的每一个时段加1，直到上涨中断；a,b,c为常系数，每个目标流域对应不同的参数；
b.当Wi,j＜α1、HPi,6＞α2且HPi,12＜α3时，流域发生超渗产流，产流公式如下：
Figure BDA0003100148380000031
式中，ft为t时刻流域下渗能力；Wt为t时刻流域土壤含水量；C为稳定入渗率，D为土壤吸渗率，均为常参数；
c.当Wi,j＜α1、HPi,6＞α2且HPi,12＞α3时，流域发生超渗产流，产流公式如下：
ft＝m(WM-Wt)n+fc
式中，WM为田间持水量；m、n为常参数，与流域土壤性质有关；fc为流域稳定下渗能力；
d.当Wi,j＞α1、HPi,6＜α2且HPi,12＜α3，或者Wi,j＞α1、HPi,6＜α2且HPi,12＞α3时，流域发生蓄满产流，产流公式如下：
Figure BDA0003100148380000041
式中，R为径流量；PE为除去蒸发的降雨量；W0为初始土壤含水量，WM为田间持水量；W'mm为最大土壤含水量；A，B为蓄水容量曲线系数；
e.当Wi,j＞α1、HPi,6＞α2且HPi,12＞α3，流域发生超渗产流，产流公式如下：
Figure BDA0003100148380000042
式中，KF为渗透系数。
进一步地，所述步骤(3)中，各子流域的汇流单位线由以下方式确定：
(3-1)根据流域地形地貌特征，确定各子流域无因次综合单位线：
Figure BDA0003100148380000043
式中，Li为第i个子流域的河长，Ji为第i个子流域的坡降，θi为第i个子流域特征参数值，mi为第i个子流域滞时参数；通过θi与mi，查表确定各子流域对应无因次综合单位线ui～xi；
(3-2)通过统计流域历史洪水过程，确定各子流域起涨历时：
Figure BDA0003100148380000044
式中，vi为第i个子流域滞时；
Figure BDA0003100148380000046
为第i个子流域起涨历时；Ki为第i个子流域单位线参数；Δt为流域单位线计算时段；
(3-3)建立无因次综合单位线与时段单位线转化关系，确定各子流域的汇流单位线：
Figure BDA0003100148380000045
式中，Fi为第i个子流域面积；qi～ti为第i个子流域的汇流单位线。
进一步地，所述步骤(3)具体包括：
将所述产流过程代入对应的汇流单位线，得到各子流域的径流过程；将所述各子流域的径流过程对应的流量过程进行错位累加，得到所述目标流域的初始洪水过程。
进一步地，所述步骤(4)中，
所述洪峰目标Obj1表示为：
Figure BDA0003100148380000051
所述峰现时间目标Obj2表示为：
Figure BDA0003100148380000052
所述确定性目标Obj3表示为：
Figure BDA0003100148380000053
式中，Qobs,i为流量实测值；Qsim,i为流量预测值；Tobs,i为峰现实测值；Tsim,,i为峰现预测值；
Figure BDA0003100148380000054
为实测流量均值；Q'obs,i为实测场次洪水洪峰值；Q'sim,i为场次洪水预测值；N为场次洪水数。
进一步地，所述步骤(4)中，对各子流域确定的最优产流模式涉及的参数进行优化率定包括：
以洪峰目标、峰现时间目标、确定性目标作为评价指标，采用SCEUA优化算法对各子流域确定的最优产流模式涉及的参数进行优化率定。
本发明另一方面提供了一种产流模式自适应的流域水文预报系统，包括：
获取与划分模块，用于获取目标流域洪水预报相关资料，并将所述目标流域划分为多个子流域；
产流模式判断模块，用于基于各子流域前期土壤含水量因子和降雨过程因子与相应阈值的大小关系，确定各子流域最优产流模式并拟合相应的产流过程；
初始洪水过程模拟模块，用于将所述产流过程代入对应的汇流单位线，并模拟所述目标流域的初始洪水过程；
最终洪水过程模拟模块，用于以洪峰目标、峰现时间目标、确定性目标作为评价指标，对各子流域确定的最优产流模式涉及的参数进行优化率定；将参数优化后拟合的产流过程代入对应的汇流单位线，得到所述目标流域的最终洪水过程。
总体而言，通过本发明所构思的以上技术方案与现有技术相比，能够取得下列有益效果：
(1)本发明根据流域不同前期土壤含水量特征及降雨过程特征，构造了基于蓄满产流、超渗产流和混合产流于一体的流域产流模式判别方法，对中小流域产流过程进行解析。从而，提高了中小流域径流预报精度。
(2)本发明根据各子流域地形、地貌特征，确定了流域无因次单位线方程，并通过统计流域历史洪水过程，确定各子流域起涨历时，进而建立无因次单位线与时段单位线转化关系，确定各子流域汇流单位线。对中小流域汇流过程进行解析，提高了中小流域径流预报精度。
附图说明
图1是本发明实施例提供的一种产流模式自适应的流域水文预报方法的流程图；
图2是本发明实施例研究区域的站点分布及其河流分布图；
图3是本发明实施案例所选洪水过程及其降雨过程图；
图4(a)和图4(b)是各下渗模型计算的流域净雨过程图；
图5是汇流单位线过程图；
图6是各下渗模型模拟洪水与实测洪水对比图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。
如图1所示，为本发明实施例提供的一种产流模式自适应的流域水文预报方法的流程图，包括如下步骤：
(1)获取目标流域洪水预报相关资料，并将所述目标流域划分为多个子流域；
具体的，目标流域洪水预报相关资料包括：雨量站逐小时降雨数据、蒸发数据、各站点经纬度数据、前期M天日降雨量数据、以及流量站对应时间点的连续逐小时径流数据；
所述将所述目标流域划分为多个子流域包括：利用泰森多边形法将所述目标流域划分为多个子流域。
(2)基于各子流域前期土壤含水量因子和降雨过程因子与相应阈值的大小关系，确定各子流域最优产流模式并拟合相应的产流过程；
具体的，所述步骤(2)具体包括如下子步骤：
(2-1)计算各子流域前期土壤含水量：
Figure BDA0003100148380000071
式中，Wi,j为第i个子流域在第j日对应的场次洪水发生初始阶段的土壤含水量；Pj为j日的日降雨量；Pa,j-1为j-1日前期影响雨量；n为影响本次径流的前期降雨天数；k为土壤消退系数；当Pa,j＞Wm时，Pa,j＝Wm，其中，Wm为流域平均最大蓄水量；
(2-2)计算各子流域x小时和y小时连续最大降雨量：
Figure BDA0003100148380000081
式中，HPi,x为第i个子流域连续x小时最大降雨量；HPi,y为第i个子流域连续y小时最大降雨量；Pi,k为第i个子流域第k个小时的降雨量；n为连续最大降雨起始时刻；
(2-3)根据流域历史洪水过程以及下垫面特征，预先设置各子流域前期土壤含水量阈值α1、x小时连续最大降雨量阈值α2和y小时连续最大降雨量阈值α3，并确定各子流域最优产流模式：
a.当Wi,j＜α1、HPi,6＜α2且HPi,12＜α3，或者Wi,j＜α1、HPi,6＜α2且HPi,12＞α3，或者Wi,j＞α1、HPi,6＞α2且HPi,12＜α3时，流域发生混合产流，产流公式如下：
ft＝a-bWt-cTt
式中，Tt为洪水t时刻从上涨后开始计数的时段数，初始值默认为1，之后的每一个时段加1，直到上涨中断；a,b,c为常系数，每个目标流域对应不同的参数；
b.当Wi,j＜α1、HPi,6＞α2且HPi,12＜α3时，流域发生超渗产流，产流公式如下：
Figure BDA0003100148380000082
式中，ft为t时刻流域下渗能力；Wt为t时刻流域土壤含水量；C为稳定入渗率，D为土壤吸渗率，均为常参数；
c.当Wi,j＜α1、HPi,6＞α2且HPi,12＞α3时，流域发生超渗产流，产流公式如下：
ft＝m(WM-Wt)n+fc
式中，WM为田间持水量；m、n为常参数，与流域土壤性质有关；fc为流域稳定下渗能力；
d.当Wi,j＞α1、HPi,6＜α2且HPi,12＜α3，或者Wi,j＞α1、HPi,6＜α2且HPi,12＞α3时，流域发生蓄满产流，产流公式如下：
Figure BDA0003100148380000091
式中，R为径流量；PE为除去蒸发的降雨量；W0为初始土壤含水量，WM为田间持水量；Wm'm为最大土壤含水量；A，B为蓄水容量曲线系数；
e.当Wi,j＞α1、HPi,6＞α2且HPi,12＞α3，流域发生超渗产流，产流公式如下：
Figure BDA0003100148380000092
式中，KF为渗透系数。
(3)将所述产流过程代入对应的汇流单位线，并模拟所述目标流域的初始洪水过程；
具体的，各子流域的汇流单位线由以下方式确定：
(3-1)根据流域地形地貌特征，确定各子流域无因次综合单位线：
Figure BDA0003100148380000093
式中，Li为第i个子流域的河长，Ji为第i个子流域的坡降，θi为第i个子流域特征参数值，mi为第i个子流域滞时参数；通过θi与mi，查表确定各子流域对应无因次综合单位线ui～xi；
(3-2)通过统计流域历史洪水过程，确定各子流域起涨历时：
Figure BDA0003100148380000094
式中，vi为第i个子流域滞时；
Figure BDA0003100148380000105
为第i个子流域起涨历时；Ki为第i个子流域单位线参数；Δt为流域单位线计算时段；
(3-3)建立无因次综合单位线与时段单位线转化关系，确定各子流域的汇流单位线：
Figure BDA0003100148380000101
式中，Fi为第i个子流域面积；qi～ti为第i个子流域的汇流单位线。
进一步地，在各子流域汇流单位线确定之后，将所述产流过程代入对应的汇流单位线，得到各子流域的径流过程；将所述各子流域的径流过程对应的流量过程进行错位累加，得到所述目标流域的初始洪水过程。
(4)以洪峰目标、峰现时间目标、确定性目标作为评价指标，对各子流域确定的最优产流模式涉及的参数进行优化率定；将参数优化后拟合的产流过程代入对应的汇流单位线，得到所述目标流域的最终洪水过程。
具体的，将步骤(1)中各子流域的降雨、蒸发、径流数据输入步骤(2)和(3)中，并分别以洪峰目标、峰现时间目标、确定性目标作为评价指标，采用SCEUA算法对各产流模型参数进行优化率定，其中目标函数公式如下：
洪峰目标Obj1表示为：
Figure BDA0003100148380000102
峰现时间目标Obj2表示为：
Figure BDA0003100148380000103
确定性目标Obj3表示为：
Figure BDA0003100148380000104
式中，Qobs,i为流量实测值；Qsim,i为流量预测值；Tobs,i为峰现实测值；Tsim,,i为峰现预测值；
Figure BDA0003100148380000111
为实测流量均值；Q'obs,i为实测场次洪水洪峰值；Q'sim,i为场次洪水预测值；N为场次洪水数。
对洪峰目标、峰现时间目标、确定性目标函数进行权重分配，得到综合目标：
Obj＝mObj1+nObj2+pObj3
其中，m、n、p为权重系数，且m+n+p＝1；优选的，Obj＝0.4Obj1+0.4Obj2+0.2Obj3。
具体SCEUA算法计算步骤如下：
(4-1)初始化，假定待优化参数个数为n，种群数位s,分区数位p；
(4-2)在可行空间内，随机产生一组样本：X1,X2,…,Xs，并计算对应的目标函数值Fi(i＝1,2,…,s)；
(4-3)将样本(X1,X2,…,Xs)按相对应的目标函数值以升序排列，并将其存入数组D(Xi,Fi)(i＝1,2,…,s)；
(4-4)将D按顺序依次放入p个分区中，第k个分区中存入样本点为：Ak:{Xp(i-1)+k}(i＝1,2,…,m)，其中
Figure BDA0003100148380000112
(4-5)对每个分区分别运用CCE算法独立进行演算；
(4-6)将进化演算后的样本再依次以目标函数值升序排列存入数组D，并以(4-4)中的方法将样本点再依次放入p个分区；
(4-7)检查样本点是否满足收敛条件，如果满足则计算结束，否则继续重复(4-4)～(4-7)。
另一方面，本发明还提供了一种产流模式自适应的流域水文预报系统，包括：
获取与划分模块，用于获取目标流域洪水预报相关资料，并将所述目标流域划分为多个子流域；
产流模式判断模块，用于基于各子流域前期土壤含水量因子和降雨过程因子与相应阈值的大小关系，确定各子流域最优产流模式并拟合相应的产流过程；
初始洪水过程模拟模块，用于将所述产流过程代入对应的汇流单位线，并模拟所述目标流域的初始洪水过程；
最终洪水过程模拟模块，用于以洪峰目标、峰现时间目标、确定性目标作为评价指标，对各子流域确定的最优产流模式涉及的参数进行优化率定；将参数优化后拟合的产流过程代入对应的汇流单位线，得到所述目标流域的最终洪水过程。
上述产流模式自适应的流域水文预报系统中各个模块的划分仅用于举例说明，在其他实施例中，可将产流模式自适应的流域水文预报系统按照需要划分为不同的模块，以完成上述产流模式自适应的流域水文预报系统的全部或部分功能。
为了更清晰的展现本发明的目的和结构及技术方案，以下以梅州流域结合附图对本发明进一步详细说明，具体实施步骤包括：
(1)收集梅州流域内尖山站点以上各水文站点的洪水过程摘录资料，利用泰森多边形法将目标流域划分为多个子流域；
具体地，梅州流域尖山站点以上区域面积1578km2，共有21个雨量站点和一个流量站点，研究区域的站点分布及其河流分布如图2所示，选取2013年8月17日～8月19日的实测洪水为例进行实施例研究，场次降雨及洪水过程如图3所示。
(2)根据梅州流域历史洪水过程及下垫面条件，给定流域前期土壤含水量阈值及连续6h、12h最大降水量阈值分别为45、40、70；计算流域前期土壤含水量和连续6h、12h最大降水量，确定各子流域最优产流模式并拟合相应的产流过程；
具体地，流域前期土壤含水量和连续6h、12h最大降水量见下表1；各产流模型计算的流域净雨过程线如图4(a)、图4(b)所示。
(3)基于子流域的地形地貌，采用无因次单位线计算各流域的汇流单位线；将所述产流过程代入对应的汇流单位线，并模拟所述目标流域的初始洪水过程；
具体地，流域汇流单位线如图5所示，各模型模拟的洪水过程如图6所示。
(4)分别以洪峰相对误差、峰现时间误差和确定性系数为目标函数，采用SCEUA优化算法对各子流域确定的最优产流模式涉及的参数进行优化率定；将参数优化后拟合的产流过程代入对应的汇流单位线，得到所述目标流域的最终洪水过程；
具体地，各模型最优参数见下表2，所计算的洪峰相对误差、峰现时间误差和确定性系数见下表3。
表1流域产流模式影响因子
Figure BDA0003100148380000131
表2各下渗模型的率定最优参数
Figure BDA0003100148380000132
Figure BDA0003100148380000141
表3各下渗模型的目标函数计算值
Figure BDA0003100148380000142
根据表1数据可得，按照本方法所提产流判别式，本次洪水符合上述情形g，应该选用混合产流模式。根据图6及表2的结果显示，针对洪峰相对误差指标，垂向-混合产流模型效果最好，其洪峰相对误差为1.93％，为所有模型中误差最小，其它产流模型如：菲利浦、霍尔坦、格林-安普特以及蓄满产流模型洪峰相对误差分别为：-4.82％、18.23％、3.27％以及-17.63％；针对峰现时间误差指标，垂向-混合与格林-安普特产流模型效果最好，其峰现时间误差皆为0h，其它产流模型如：菲利浦、霍尔坦以及蓄满产流模型峰现时间误差分别为：-1h、-1h以及-5h，都比实际洪峰出现要偏早；针对确定性系数指标，垂向-混合产流模型效果最好，其确定性系数为0.974，为所有模型中最大的。其它产流模型如：菲利浦、霍尔坦、格林-安普特以及蓄满产流模型的确定性系数分别为：0.813、0.678、0.912以及0.624。综上，采用垂向-混合产流模型模拟出的洪水过程效果最好，各项指标误差最小。
综上所述，本发明公开了一种产流模式自适应的流域水文预报方法及系统，针对流域不同前期土壤含水量特征及降水过程特征，构造了基于蓄满产流、超渗产流和混合产流于一体的流域产流模式判别方法；进一步的，依据无因次单位线原理，推导出了基于小流域地形、地貌特征的无因次综合汇流单位线；并以洪峰目标、峰现时间目标、确定性目标作为评价指标，对各子流域确定的最优产流模式涉及的参数进行优化率定，进而形成完整预报体系。本发明能够针对流域不同土壤含水量特征和降雨特征，确定流域最优下渗模式，并根据流域地形地貌特征，确定流域汇流单位线，具有计算简单、适用范围广、拟合精度高的优点。
本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (7)

1.一种产流模式自适应的流域水文预报方法，其特征在于，包括如下步骤：
(1)获取目标流域洪水预报相关资料，并将所述目标流域划分为多个子流域；
(2)基于各子流域前期土壤含水量因子和降雨过程因子与相应阈值的大小关系，确定各子流域最优产流模式并拟合相应的产流过程；
所述步骤(2)具体包括如下子步骤：
(2-1)计算各子流域前期土壤含水量：
Figure FDA0003657802200000011
式中，Wi,j为第i个子流域在第j日对应的场次洪水发生初始阶段的土壤含水量；Pj为j日的日降雨量；Pa,j-1为j-1日前期影响雨量；h为影响本次径流的前期降雨天数；g为土壤消退系数；当Pa,j＞Wm时，Pa,j＝Wm，其中，Wm为流域平均最大蓄水量；
(2-2)计算各子流域x小时和y小时连续最大降雨量：
Figure FDA0003657802200000012
式中，HPi,x为第i个子流域连续x小时最大降雨量；HPi,y为第i个子流域连续y小时最大降雨量；Pi,k为第i个子流域第k个小时的降雨量；n为连续最大降雨起始时刻；
(2-3)预先设置各子流域前期土壤含水量阈值α1、x小时连续最大降雨量阈值α2和y小时连续最大降雨量阈值α3，并确定各子流域最优产流模式：
a.当Wi,j＜α1、HPi,6＜α2且HPi,12＜α3，或者Wi,j＜α1、HPi,6＜α2且HPi,12＞α3，或者Wi,j＞α1、HPi,6＞α2且HPi,12＜α3时，流域发生混合产流，产流公式如下：
ft＝a-bWt-cTt
式中，Tt为洪水t时刻从上涨后开始计数的时段数，初始值默认为1，之后的每一个时段加1，直到上涨中断；a,b,c为常系数，每个目标流域对应不同的参数；
b.当Wi,j＜α1、HPi,6＞α2且HPi,12＜α3时，流域发生超渗产流，产流公式如下：
Figure FDA0003657802200000021
式中，ft为t时刻流域下渗能力；Wt为t时刻流域土壤含水量；C为稳定入渗率，D为土壤吸渗率，均为常参数；
c.当Wi,j＜α1、HPi,6＞α2且HPi,12＞α3时，流域发生超渗产流，产流公式如下：
ft＝m(WM-Wt)r+fc
式中，WM为田间持水量；m、r为常参数，与流域土壤性质有关；fc为流域稳定下渗能力；
d.当Wi,j＞α1、HPi,6＜α2且HPi,12＜α3，或者Wi,j＞α1、HPi,6＜α2且HPi,12＞α3时，流域发生蓄满产流，产流公式如下：
Figure FDA0003657802200000022
式中，R为径流量；PE为除去蒸发的降雨量；W0为初始土壤含水量，WM为田间持水量；Wm'm为最大土壤含水量；A，B为蓄水容量曲线系数；
e.当Wi,j＞α1、HPi,6＞α2且HPi,12＞α3，流域发生超渗产流，产流公式如下：
Figure FDA0003657802200000031
式中，KF为渗透系数；
(3)将所述产流过程代入对应的汇流单位线，并模拟所述目标流域的初始洪水过程；
(4)以洪峰目标、峰现时间目标、确定性目标作为评价指标，对各子流域确定的最优产流模式涉及的参数进行优化率定；将参数优化后拟合的产流过程代入对应的汇流单位线，得到所述目标流域的最终洪水过程。
2.如权利要求1所述的方法，其特征在于，所述步骤(1)中，
所述目标流域洪水预报相关资料包括：雨量站逐小时降雨数据、蒸发数据、各站点经纬度数据、前期M天日降雨量数据、以及流量站对应时间点的连续逐小时径流数据；
所述将所述目标流域划分为多个子流域包括：利用泰森多边形法将所述目标流域划分为多个子流域。
3.如权利要求1所述的方法，其特征在于，所述步骤(3)中，各子流域的汇流单位线由以下方式确定：
(3-1)根据流域地形地貌特征，确定各子流域无因次综合单位线：
Figure FDA0003657802200000032
式中，Li为第i个子流域的河长，Ji为第i个子流域的坡降，θi为第i个子流域特征参数值，mi为第i个子流域滞时参数；通过θi与mi，查表确定各子流域对应无因次综合单位线ui～xi；
(3-2)通过统计流域历史洪水过程，确定各子流域起涨历时：
Figure FDA0003657802200000041
式中，vi为第i个子流域滞时；
Figure FDA0003657802200000047
为第i个子流域起涨历时；Ki为第i个子流域单位线参数；Δt为流域单位线计算时段；
(3-3)建立无因次综合单位线与时段单位线转化关系，确定各子流域的汇流单位线：
Figure FDA0003657802200000042
式中，Fi为第i个子流域面积；qi～ti为第i个子流域的汇流单位线。
4.如权利要求1或3所述的方法，其特征在于，所述步骤(3)具体包括：
将所述产流过程代入对应的汇流单位线，得到各子流域的径流过程；将所述各子流域的径流过程对应的流量过程进行错位累加，得到所述目标流域的初始洪水过程。
5.如权利要求1所述的方法，其特征在于，所述步骤(4)中，
所述洪峰目标Obj1表示为：
Figure FDA0003657802200000043
所述峰现时间目标Obj2表示为：
Figure FDA0003657802200000044
所述确定性目标Obj3表示为：
Figure FDA0003657802200000045
式中，Qobs,i为流量实测值；Qsim,i为流量预测值；Tobs,i为峰现实测值；Tsim,,i为峰现预测值；
Figure FDA0003657802200000046
为实测流量均值；Q'obs,i为实测场次洪水洪峰值；Q'sim,i为场次洪水预测值；N为场次洪水数。
6.如权利要求1或5所述的方法，其特征在于，所述步骤(4)中，对各子流域确定的最优产流模式涉及的参数进行优化率定包括：
以洪峰目标、峰现时间目标、确定性目标作为评价指标，采用SCEUA优化算法对各子流域确定的最优产流模式涉及的参数进行优化率定。
7.一种产流模式自适应的流域水文预报系统，其特征在于，包括：
获取与划分模块，用于获取目标流域洪水预报相关资料，并将所述目标流域划分为多个子流域；
产流模式判断模块，用于基于各子流域前期土壤含水量因子和降雨过程因子与相应阈值的大小关系，确定各子流域最优产流模式并拟合相应的产流过程；
所述产流模式判断模块具体用于执行如下子步骤：
(2-1)计算各子流域前期土壤含水量：
Figure FDA0003657802200000051
式中，Wi,j为第i个子流域在第j日对应的场次洪水发生初始阶段的土壤含水量；Pj为j日的日降雨量；Pa,j-1为j-1日前期影响雨量；h为影响本次径流的前期降雨天数；g为土壤消退系数；当Pa,j＞Wm时，Pa,j＝Wm，其中，Wm为流域平均最大蓄水量；
(2-2)计算各子流域x小时和y小时连续最大降雨量：
Figure FDA0003657802200000052
式中，HPi,x为第i个子流域连续x小时最大降雨量；HPi,y为第i个子流域连续y小时最大降雨量；Pi,k为第i个子流域第k个小时的降雨量；n为连续最大降雨起始时刻；
(2-3)预先设置各子流域前期土壤含水量阈值α1、x小时连续最大降雨量阈值α2和y小时连续最大降雨量阈值α3，并确定各子流域最优产流模式：
a.当Wi,j＜α1、HPi,6＜α2且HPi,12＜α3，或者Wi,j＜α1、HPi,6＜α2且HPi,12＞α3，或者Wi,j＞α1、HPi,6＞α2且HPi,12＜α3时，流域发生混合产流，产流公式如下：
ft＝a-bWt-cTt
式中，Tt为洪水t时刻从上涨后开始计数的时段数，初始值默认为1，之后的每一个时段加1，直到上涨中断；a,b,c为常系数，每个目标流域对应不同的参数；
b.当Wi,j＜α1、HPi,6＞α2且HPi,12＜α3时，流域发生超渗产流，产流公式如下：
Figure FDA0003657802200000061
式中，ft为t时刻流域下渗能力；Wt为t时刻流域土壤含水量；C为稳定入渗率，D为土壤吸渗率，均为常参数；
c.当Wi,j＜α1、HPi,6＞α2且HPi,12＞α3时，流域发生超渗产流，产流公式如下：
ft＝m(WM-Wt)r+fc
式中，WM为田间持水量；m、r为常参数，与流域土壤性质有关；fc为流域稳定下渗能力；
d.当Wi,j＞α1、HPi,6＜α2且HPi,12＜α3，或者Wi,j＞α1、HPi,6＜α2且HPi,12＞α3时，流域发生蓄满产流，产流公式如下：
Figure FDA0003657802200000062
式中，R为径流量；PE为除去蒸发的降雨量；W0为初始土壤含水量，WM为田间持水量；W'mm为最大土壤含水量；A，B为蓄水容量曲线系数；
e.当Wi,j＞α1、HPi,6＞α2且HPi,12＞α3，流域发生超渗产流，产流公式如下：
Figure FDA0003657802200000071
式中，KF为渗透系数；
初始洪水过程模拟模块，用于将所述产流过程代入对应的汇流单位线，并模拟所述目标流域的初始洪水过程；
最终洪水过程模拟模块，用于以洪峰目标、峰现时间目标、确定性目标作为评价指标，对各子流域确定的最优产流模式涉及的参数进行优化率定；将参数优化后拟合的产流过程代入对应的汇流单位线，得到所述目标流域的最终洪水过程。
CN202110621866.4A 2021-06-04 2021-06-04 一种产流模式自适应的流域水文预报方法及系统 Active CN113435630B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110621866.4A CN113435630B (zh) 2021-06-04 2021-06-04 一种产流模式自适应的流域水文预报方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110621866.4A CN113435630B (zh) 2021-06-04 2021-06-04 一种产流模式自适应的流域水文预报方法及系统

Publications (2)

Publication Number Publication Date
CN113435630A CN113435630A (zh) 2021-09-24
CN113435630B true CN113435630B (zh) 2022-07-15

Family

ID=77803725

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110621866.4A Active CN113435630B (zh) 2021-06-04 2021-06-04 一种产流模式自适应的流域水文预报方法及系统

Country Status (1)

Country Link
CN (1) CN113435630B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114943051B (zh) * 2022-05-16 2023-05-23 广东省水文局梅州水文分局 一种基于综合降雨因子的洪峰预报方法

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106529176B (zh) * 2016-11-11 2018-04-06 中国水利水电科学研究院 一种双核双驱洪水预报方法
CN107729695A (zh) * 2017-11-21 2018-02-23 中国水利水电科学研究院 一种用于小流域次洪模拟的水文模型率定方法
CN110928965B (zh) * 2018-12-13 2022-11-04 河海大学 一种基于流域精细分类的多模型灵活架构的模拟方法
US11428078B2 (en) * 2019-07-11 2022-08-30 Halliburton Energy Services, Inc. Systems and methods for forecasting well productivity
CN110781259B (zh) * 2019-09-18 2021-02-09 河海大学 一种基于地貌单位线的水文模型的构建方法
CN111861014B (zh) * 2020-07-23 2021-03-16 广东省水文局梅州水文分局 一种分布式单位线洪水预报方法

Also Published As

Publication number Publication date
CN113435630A (zh) 2021-09-24

Similar Documents

Publication Publication Date Title
CN111651885A (zh) 一种智慧型海绵城市洪涝预报方法
CN111665575B (zh) 一种基于统计动力的中长期降雨分级耦合预报方法及系统
CN110232471B (zh) 一种降水传感网节点布局优化方法及装置
CN108021773B (zh) 基于dss数据库的分布式水文模型多场次洪水参数率定方法
CN112686426B (zh) 一种基于水电站流域关键点的来水量预警方法及系统
CN114429053B (zh) 流域尺度wefe系统适配性模拟优化方法
CN114741987B (zh) 考虑洪水预报模型绝对误差拟合残差分布的洪水概率预报模型
CN110134907B (zh) 一种降雨缺失数据填补方法、系统及电子设备
CN111737853A (zh) 一种基于swmm模型的低影响开发多目标区间优化配置方法
CN115238947A (zh) 气候变化下旱涝急转事件的社会经济暴露度预估方法
CN113807545A (zh) 一种基于深度学习和物理模型的河湖生态流量预报预警方法
CN110930282B (zh) 一种基于机器学习的局地降雨雨型分析方法
CN113435630B (zh) 一种产流模式自适应的流域水文预报方法及系统
CN115759445A (zh) 一种基于机器学习和云模型的分类洪水随机预报方法
CN111611541A (zh) 基于Copula函数的无资料地区降水数据推求方法及系统
CN110580377A (zh) 一种基于响应曲面法的植草沟设计参数优化方法
CN114969655A (zh) 流域输沙量的模拟估算方法
CN116305933B (zh) 一种基于dem数据的坡面产汇流简易计算方法及装置
CN112116176B (zh) 基于气象因素的地区绝缘子群污闪风险评估方法和系统
CN115965243A (zh) 一种用于土地利用变化下洪涝风险动态评价的模型链
CN114491978B (zh) 基于时变参数水文不确定性处理器的日模型实时预报方法
CN113094876B (zh) 未来情景下基于arima模型的水库防洪风险计算方法和系统
CN114611290B (zh) 一种基于量变参数水文不确定性处理器的场次洪水水文模型实时预报方法
CN117350180B (zh) 一种复杂工况下洪水快速预测模型集合构建方法
CN114528708B (zh) 一种流域侵蚀产沙指数模拟方法和系统

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant