CN113190958B - 一种全分布式考虑流域内水库影响的汇流模拟方法 - Google Patents

一种全分布式考虑流域内水库影响的汇流模拟方法 Download PDF

Info

Publication number
CN113190958B
CN113190958B CN202110322305.4A CN202110322305A CN113190958B CN 113190958 B CN113190958 B CN 113190958B CN 202110322305 A CN202110322305 A CN 202110322305A CN 113190958 B CN113190958 B CN 113190958B
Authority
CN
China
Prior art keywords
grid
reservoir
processed
target
unit
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
CN202110322305.4A
Other languages
English (en)
Other versions
CN113190958A (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.)
Hohai University HHU
Original Assignee
Hohai University HHU
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 Hohai University HHU filed Critical Hohai University HHU
Priority to CN202110322305.4A priority Critical patent/CN113190958B/zh
Publication of CN113190958A publication Critical patent/CN113190958A/zh
Application granted granted Critical
Publication of CN113190958B publication Critical patent/CN113190958B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • 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
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00—Details relating to the application field
    • G06F2113/08—Fluids
    • 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

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开了一种全分布式考虑流域内水库影响的汇流模拟方法，方法包括获取目标流域中各栅格单元的坡度、流向和汇流累计值，基于水库位置，将目标流域中的栅格单元划分为水库栅格和非水库栅格，构建水库栅格的蓄泄模型，基于该模型获取各水库栅格的泄流量，基于MCT方法构建非水库栅格单元的汇流模型，基于该模型获取各非水库栅格的泄流量。本发明提供的方法应用目标流域的数字高程数据，数据来源稳定可靠，方法考虑了流域内水库对汇流计算的影响，计算结果精确合理，同时解决了有水库的流域的汇流计算问题，有利于数字水文学及水文模型研究的深入发展。

Description

一种全分布式考虑流域内水库影响的汇流模拟方法
技术领域
本发明涉及水文技术领域，具体涉及一种全分布式考虑流域内水库影响的汇流模拟方法。
背景技术
洪水灾害是最常见的自然灾害，给国民经济和人民生命财产带来巨大损失。近年来，受气候变化影响，由局地强降水造成的中小河流突发性洪水频繁发生，已成为造成人员伤亡的主要灾种，对中小流域突发性洪水进行快速预报成为亟待解决的重要问题。减轻洪水灾害的方法包括堤防建设等工程措施及洪水预警预报等非工程措施。
随着遥感、地理信息以及数字流域等技术的发展，基于栅格数字高程模型(DEM，Digital Elevation Model)的分布式水文模型以其充分考虑降雨和下垫面条件空间变化的特点，现已有成为非工程措施中洪水预报的主要工具的趋势，发展分布式水文模型，对应对洪水灾害有重要意义。在构建分布式的水文模型时，通常将流域划分为若干个正交的网格，以每一个网格作为计算单元。在每一个计算单元中通过产流模型计算产生的径流，然后按照一定的汇流顺序逐单元进行汇流计算。
水库作为防洪减灾工程措施的重要组成部分，可以拦蓄洪水，削减进入下游河道的洪峰流量，达到减免洪水灾害的目的。但是，大规模的水库建设，显著改变了流域尺度的水文循环和流域的自然汇流过程，给流域洪水预报带来了新的挑战。目前洪水预报中考虑库塘的影响，主要是对预报洪量进行修正或者是在集总式模型进行改进，在分布式水文模型的汇流过程中考虑水库的影响，是洪水预报必然的发展趋势。一般的水库均有设计时的特征库容资料，包括死库容、兴利库容、防洪库容等等，但在众多中小流域中，存在大量的无资料小水库和塘坝，这些小库塘对洪水过程影响显著，却通常仅有总库容资料，缺乏特征库容资料。如何在汇流中考虑这一部分库塘的影响，是分布式水文模型建模的重点和难点之一。
为了解决以上难题，进一步促进分布式水文模型中汇流计算的发展，需要更深入理解水库对于汇流计算的影响，研究基于水库真实位置，将水库蓄泄过程加入基于正交网格的汇流演算中，并且同时考虑了有资料的水库和无资料的水库，从而构建一种全分布式考虑流域内水库影响的汇流模拟方法，正是发明人需要解决的问题。
发明内容
本发明的目的：提供一种基于水库真实位置，计及水库蓄泄洪过程的汇流模拟方法。
技术方案：本发明提供的一种全分布式考虑流域内水库影响的汇流模拟方法，用于获取目标流域内各栅格单元在目标时刻t时的泄流量，方法包括如下步骤：
步骤1：基于目标流域的数字高程数据，获取目标流域中各栅格单元的坡度、流向、汇流累计值；基于各栅格单元的流向和汇流累计值，计算各栅格单元的演算次序；然后进入步骤2；
步骤2：根据目标流域内各水库的位置，将各栅格单元区分为水库栅格和非水库栅格；然后进入步骤3；
步骤3：分别针对各水库栅格i，将其作为第一待处理栅格，依次执行步骤3.1.1至步骤3.1.2：
步骤3.1.1：判断第一待处理栅格i在目标时刻t时的库容ResStroagei,t是否大于或等于其总库容TotalStoragei，是则根据如下公式获取第一待处理栅格i在目标时刻t时的泄流量ResDrainagei,t；否则进入步骤3.1.2；
ResDrainagei,t＝(TotalStoragei-ComputeStoragei)*RDC
+(ResStroagei,t-TotalStoragei)
其中：ComputeStoragei为第一待处理栅格i的运行库容；RDC为目标流域的水库泄流系数；
步骤3.1.2：判断第一待处理栅格i在目标时刻t时的库容ResStroagei,t是否大于或等于其运行库容ComputeStoragei，是则跟据如下公式获取第一待处理栅格i在目标时刻t时的泄流量ResDrainagei,t，否则取第一待处理栅格i在目标时刻t时的泄流量ResDrainagei,t＝0；
ResDrainagei,t＝(ResStroagei,t-ComputeStoragei)*RDC
通过执行步骤3.1.1至3.1.2，获取第一待处理栅格对应目标时刻t的泄流量，进而获取各水库栅格分别对应目标时刻t的泄流量；
分别针对演算次序等于1的各非水库栅格，执行如下操作：根据产流模型，获取各非水库栅格在目标时刻t时的产流量，将该产流量作为对应非水库栅格的泄流量；
针对演算次序不等于1的各非水库栅格，按照演算次序由小到大的顺序依次获取各非水库栅格的泄流量，获取对应非水库栅格的泄流量的方法包括如下步骤：
分别针对演算次序不等于1的各非水库栅格，将其作为第二待处理栅格i'，根据如下公式：
Oi′,t＝C1Ii′,t+C2Ii′,t-1+C3Oi′,t-1
获取第二待处理栅格目标时刻t的泄流量Oi′,t，进而获取各非水库栅格分别对应目标时刻t的泄流量；
其中，C1、C2、C3为和第二待处理栅格i'的参考流量以及流速相关的系数，Ii',t为第二待处理栅格i'在目标时刻t时的入流量，Ii',t-1为第二待处理栅格i'在时刻t-1时的入流量,Oi′,t-1为第二待处理栅格i'时刻t-1时的泄流量。
作为本发明的一种优选方案，在步骤1中，获取演算次序的方法包括如下步骤：
步骤1.1：分别针对汇流累计值等于1的各栅格单元，将其演算次序赋值为1；分别针对汇流累计值不等于1的各栅格单元，初始化其演算次序n＝0；然后进入步骤1.2；
步骤1.2：分别针对各演算次序n＝0的栅格单元，将其作为待处理单元；同时分别针对各待处理单元，执行如下操作，然后进入步骤1.3；
根据与待处理单元相邻各栅格单元的流向，获取水流流向待处理单元的各相邻栅格单元：
若水流流向待处理单元的各栅格单元的演算次序均不等于零，则对水流流向待处理单元的各相邻栅格单元中的演算次序最大值加1获取一累计值，将该累计值作为待处理单元的演算次序；
若至少一个水流流向待处理单元的相邻栅格单元的演算次序等于零，则不改变待处理单元的演算次序；
步骤1.3：判断是否存在演算次序n＝0的栅格单元，是则返回步骤1.2；否则停止循环。
作为本发明的一种优选方案，在步骤1中，根据D8流向法获取各栅格单元的流向。
作为本发明的一种优选方案，在步骤3.1.1中：在初始时刻时，第一待处理栅格i的库容ResStroagei,0＝computeStoragei；
在初始时刻之后的目标时刻t时，根据如下公式：
ResStroagei,t＝ResStroage′i,t-1+Ii,t
获取第一待处理栅格i在目标时刻t时的库容ResStroagei,t；
其中，Ii,t为第一待处理栅格i在目标时刻t时的入流量，ResStroage′i,t-1第一待处理栅格i在目标时刻t-1时的更新库容。
作为本发明的一种优选方案，根据如下公式：
ResStroage′i,t-1＝ResStroagei,t-1-ResDrainagei,t-1
获取第一待处理栅格i在目标时刻t-1时的更新库容ResStroage′i,t-1；
其中，ResDrainagei,t-1为第一待处理栅格i在目标时刻t-1时的泄流量，ResStroagei,t-1为第一待处理栅格i在目标时刻t-1时的库容。
作为本发明的一种优选方案，根据如下公式：
Il,t-1＝Ol′,t-1+Rl,t-1
Il,t＝Ol′,t+Rl,t
获取目标流域中各栅格单元l在目标时刻t时的入流量ll,t，以及各栅格单元l在时刻t-1时的入流量ll,t-1；
其中，Ol′,t-1和Ol′,t分别为与各栅格单元l相邻的上游栅格单元l′在t-1时刻和t时刻的泄流量，如l′为水库栅格，即为ResDrainagel′,t-1和ResDrainagel′,t；Rl,t-1和Rl,t分别为根据产流模型计算得到的各栅格单元l在t-1时刻和t时刻的产流量。
作为本发明的一种优选方案，在步骤3.1.1中，根据如下方法获取第一待处理栅格i的运行库容ComputeStoragei：
若和第一待处理栅格i所对应水库的相关的资料中记载了其兴利库容ActiveStoragei和死库容DeadStoragei，根据如下公式：
ComputeStoragei＝ActiveStoragei+DeadStoragei
获取第一待处理栅格i的运行库容ComputeStoragei；
若和第一待处理栅格i所对应水库相关的资料中未记载其兴利库容ActiveStoragei和死库容DeadStoragei，根据如下公式：
Figure BDA0002993358280000051
获取第一待处理栅格i的运行库容ComputeStoragei；
其中，根据如下公式：
Figure BDA0002993358280000052
获取目标流域内的平均库容系数
Figure BDA0002993358280000053
其中，N为目标流域内兴利库容和死库容均已知的水库的总数量，CSCn为目标流域内兴利库容和死库容均已知的第n个水库的库容系数，ActiveStoragen为目标流域内兴利库容和死库容均已知的第n个水库的兴利库容；DeadStoragen为与目标流域内兴利库容和死库容均已知的第n个水库的死库容；TotalStoragen为与目标流域内兴利库容和死库容均已知的第n个水库的总库容。
作为本发明的一种优选方案，在步骤3.1.2中，根据如下公式：
Figure BDA0002993358280000054
Figure BDA0002993358280000055
Figure BDA0002993358280000056
获取和第二待处理栅格i'的参考流量以及流速相关系数C1、C2、C3；
其中，CCt-1、CCt、CDt、CDt-1为中间参数，根据如下公式
Figure BDA0002993358280000057
Figure BDA0002993358280000058
Figure BDA0002993358280000059
Figure BDA00029933582800000510
获取中间系数CCt-1、CCt、CDt、CDt-1；
其中，Δt为预设的计算时段长；Δx为预设的计算空间步长，cvt-1和cvt分别为第二待处理栅格i'在t-1时刻和t时刻的流速，βt-1和βt分别为第二待处理栅格i'在t-1时刻和t时刻的预设校正因子；Bt-1和Bt为与第二待处理栅格i'在t-1时刻和t时刻的所对应的河面宽度值，S0为与第二待处理栅格i'所对应的河底坡度；Qt-1和Qt分别为第二待处理栅格i'在t-1时刻和t时刻的参考流量。
有益效果：相对于现有技术，本发明提供的方法，研究基于水库真实位置，将水库蓄泄过程加入基于正交网格的汇流演算中，从而构建一种全分布式考虑流域内水库影响的汇流模拟方法，既保证了计算结果的精度与可靠性，同时解决了无资料水库较多流域的汇流计算问题。且本方法主要应用流域数字高程模型，数据来源稳定可靠，方法中变量之间的函数关系明确，有利于在计算机上编程实现，保证了结果的客观合理性，可以进一步促进数字水文学以及分布式模型的深入发展。
附图说明
图1是根据本发明实施例提供的汇流模拟方法流程框图；
图2是根据本发明实施例提供的流向栅格示意图；
图3是根据本发明实施例提供的汇流累计栅格示意图；
图4是根据本发明实施例提供的栅格演算次序计算方法示意图；
图5是根据本发明实施例提供的演算次序栅格示意图；
图6是根据本发明实施例提供的屯溪流域水库栅格示意图；
图7是根据本发明实施例提供的屯溪流域水库蓄泄量示意图；
图8是根据本发明实施例提供的屯溪流域汇流模拟结果示意图。
具体实施方式
下面结合附图对本发明作进一步描述。以下实施例仅用于更加清楚地说明本发明的技术方案，而不能以此来限制本发明的保护范围。
本发明提供的一种全分布式考虑流域内水库影响的汇流模拟方法，用于获取目标流域内各栅格单元在目标时刻t时的泄流量。
方法包括如下步骤：
步骤1：基于目标流域的数字高程数据，计算目标流域中各栅格单元的流向、坡度、汇流累计值；并基于各栅格单元的流向和汇流累计值，计算各栅格单元的演算次序；根据计算出的各参数，构建对应的栅格文件，然后进入步骤2。
获取各栅格的汇流累计值的方法为：
步骤a、初始化汇流累计栅格，将流域内所有栅格单元的汇流累积值设置为1。
步骤b、根据各栅格单元的高程值，以栅格单元C为中心，获取与栅格单元C相邻的各栅格单元中比栅格单元C的低的各栅格单元，提取其中高程最低的栅格单元CD，并计算C和CD之间的高程差ΔE和距离Distance，基于ΔE和Distance计算栅格单元C的坡度Slope；
Slope＝ΔE/Distance
步骤c、将栅格单元CD作为入流栅格，将栅格单元C作为出流栅格，如果入流栅格CD的汇流累计值为1，则栅格CD的汇流累积值为栅格C的汇流累计值加1；如果CD的汇流累积值不为1，则栅格CD的汇流累积值更新为CD原先的汇流累积值加上栅格C的汇流累计值；
步骤d、基于D8流向法获取栅格单元C中水流的流向：通过栅格单元CD和栅格单元C的空间位置关系，确定栅格单元C中水流的流向，将水流分为八个方向：若CD在C的右方，则C中的流向设置为1；若CD在C的右下方，则C中的流向设置为2；若CD在C的下方，则C中的流向设置为3；若CD在C的左下方，则C中的流向设置为4；若CD在C的左方，则C中的流向设置为5；若CD在C的左上方，则C中的流向设置为6；若CD在C的上方，则C中的流向设置为7；若CD在C的右上方，则C中的流向设置为0；
步骤e、逐栅格循环重复步骤b至步骤d，计算得到每一个栅格单元的坡度、汇流累计值，以及各栅格的流向，进而获取如图3所示的汇流累计栅格和如图2所示的流向栅格。
获取各栅格的演算次序的方法包括如下步骤：
步骤1.1：分别针对汇流累计值等于1的各栅格单元，将其演算次序赋值为1；分别针对汇流累计值不等于1的各栅格单元，初始化其演算次序n＝0；然后进入步骤1.2；
步骤1.2：分别针对各演算次序n＝0的栅格单元，将其作为待处理单元；同时分别针对各待处理单元，执行如下操作，然后进入步骤1.3；
根据与待处理单元相邻各栅格单元的流向，获取水流流向待处理单元的各相邻栅格单元：
若水流流向待处理单元的各栅格单元的演算次序均不等于零，则对水流流向待处理单元的各相邻栅格单元中的演算次序最大值加1获取一累计值，将该累计值作为待处理单元的演算次序；
若至少一个水流流向待处理单元的相邻栅格单元的演算次序等于零，则不改变待处理单元的演算次序；
步骤1.3：判断是否存在演算次序n＝0的栅格单元，是则返回步骤1.2；否则停止循环。
以如图4所示的局部的数个栅格单元为例说明演算次序的计算方法：
图中A、B、C、D、E分别表示栅格单元，其中栅格单元A、D的汇流累计值等于1，则将A、D的演算次序赋值为1，对于其他汇流累计值不等于1的栅格单元，将其演算次序的初始值赋值为0，执行如下操作：
第一轮循环：
对于B：水流流向B的相邻栅格为A和E，E的演算次序等于0，所以B的演算次序不变，仍为0；
对于C：水流流向C的相邻栅格分别为B，B的演算次序是0，所以C的演算次序不变，仍为0；
对于E：水流流向E的相邻栅格只有D，D的演算次序是1，对D的演算次序加1，将加1后的结果赋值给E的演算次序，E的演算次序变为2；
第一轮循环结束，其中存在演算次序等于0的栅格，继续第二轮循环：
第二轮循环：
对于B：水流流向B的相邻栅格为A和E，A、E的演算次序均不等于0，所以B的演算次序为A、E中演算次序最大的值加1，将加1后的结果赋值给B的演算次序，B的演算次序变为3；
对于C：水流流向C的相邻栅格为B，B的演算次序是0，所以C的演算次序不变，仍为0；
第二轮循环结束，其中存在演算次序等于0的栅格，继续第三轮循环：
第三轮循环：
对于C：水流流向C的相邻栅格只有B，B的演算次序不等于0，所以C的演算次序为B的值加1，将加1后的结果赋值给C的演算次序，C的演算次序变为4；
第三轮循环结束，A、B、C、D、E的演算次序分别为1、3、4、1、2，不存在演算次序等于0的栅格单元，循环结束。
获取各栅格单元的演算次序后，进而获取如图5所示的演算次序栅格。
步骤2：根据目标流域内各水库的位置，将各栅格单元区分为水库栅格和非水库栅格；然后进入步骤3。
具体的，通过如下方法获取水库栅格，将目标流域中除水库栅格外的其他栅格作为非水库栅格：
步骤f：针对各栅格单元，将各栅格单元均赋值为0；
步骤g：根据目标流域中各水库的经纬度，逐水库计算得到在水库栅格矩阵Rreservoir_Raster中各水库所在栅格单元的行数ReservoirRow和列数ReservoirCol：
ReservoirColi＝Int(((ReservoirLonm-XllCorner)/DEMPrecision))+1
ReservoirRowi＝TotalRow-Int(((ReservoirLatI-YllCorner)/DEMPrecision))
式中，下标m表示第m个水库；Int(…)为取整计算；ReservoirLonm为第m个水库的经度坐标；ReservoirLatm为第m个水库的维度坐标；XllCorner为栅格矩阵左下角栅格单元经度；YllCorner为栅格矩阵左下角栅格单元纬度；DEMPrecision为栅格单元的边长。
步骤h：逐栅格单元循环，将行列数与步骤g中计算结果相匹配的栅格单元赋值为1，否则仍然为0，赋值为1的栅格单元为水库栅格，其他的栅格单元为非水库栅格。
步骤3：基于蓄泄洪模型获取各水库栅格的泄洪量，具体方法如下。
分别针对各水库栅格i，将其作为第一待处理栅格，依次执行步骤3.1.1至步骤3.1.2：
步骤3.1.1：判断第一待处理栅格i在目标时刻t时的库容ResStroagei,t是否大于或等于其总库容TotalStoragei，是则根据如下公式获取第一待处理栅格i在目标时刻t时的泄流量ResDrainagei,t；否则进入步骤3.1.2；
ResDrainagei,t＝(TotalStoragei-ComputeStoragei)*RDC+(ResStroagei,t-TotalStoragei)
其中：ComputeStoragei为第一待处理栅格i的运行库容；RDC为目标流域的水库泄流系数。其中总库容是水库建成时便确定的库容，相当于该水库兴建的总规模，每个不同水库的总库容和水库的结构相关，且不会随着时间的变化而变化；水库在某时刻的库容指的是该水库在该时刻的蓄水量，是和时间相关的量。
根据如下方法获取第一待处理栅格i的运行库容ComputeStoragei：
完整的特征库容资料是水库建成时就构建好的，包括水库的总库容、兴利库容、死库容等参数资料。当水库有完整的特征库容资料时，利用水库的死库容和兴利库容计算蓄泄洪模型中的运行库容，即：若和第一待处理栅格i所对应水库的相关的资料中记载了其兴利库容ActiveStoragei和死库容DeadStoragei，根据如下公式：
ComputeStoragei＝ActiveStoragei+DeadStoragei
获取第一待处理栅格i的运行库容ComputeStoragei；
当水库仅有总库容资料，缺乏其他特征库容资料时，利用总库容和有详细库容资料的水库的平均库容系数
Figure BDA0002993358280000101
计算蓄泄洪模型的运行库容，即：若和第一待处理栅格i所对应水库相关的资料中未记载其兴利库容ActiveStoragei和死库容DeadStoragei，根据如下公式：
Figure BDA0002993358280000102
获取第一待处理栅格i的运行库容ComputeStoragei；
其中，根据如下公式：
Figure BDA0002993358280000103
Figure BDA0002993358280000104
获取目标流域内的平均库容系数
Figure BDA0002993358280000105
其中，N为目标流域内兴利库容和死库容均已知的水库的总数量，CSCn为目标流域内兴利库容和死库容均已知的第n个水库的库容系数，ActiveStoragen为目标流域内兴利库容和死库容均已知的第n个水库的兴利库容；DeadStoragen为与目标流域内兴利库容和死库容均已知的第n个水库的死库容；TotalStoragen为与目标流域内兴利库容和死库容均已知的第n个水库的总库容。以屯溪流域为目标流域，计算得到的屯溪流域的平均库容系数
Figure BDA0002993358280000106
等于0.79。
根据如下方法获取第一待处理栅格i的t时刻的库容ResStroagei,t：
在初始时刻时，第一待处理栅格i的库容ResStroagei,0＝computeStoragei；初始时刻指的是汇流演算开始的第一个时刻。
在初始时刻之后的目标时刻t时，根据如下公式：
ResStroagei,t＝ResStroage′i,t-1+Ii,t
获取第一待处理栅格i在目标时刻t时的库容ResStroagei,t；
其中，Ii,t为第一待处理栅格i在目标时刻t时的入流量，ResStroage′i,t-1第一待处理栅格i在目标时刻t-1时的更新库容。
作为本发明的一种优选方案，根据如下公式：
ResStroage′i,t-1＝ResStroagei,t-1-ResDrainagei,t-1
获取第一待处理栅格i在目标时刻t-1时的更新库容ResStroage′i,t-1；
其中，ResDrainagei,t-1为第一待处理栅格i在目标时刻t-1时的泄流量，ResStroagei,t-1为第一待处理栅格i在目标时刻t-1时的库容。
步骤3.1.2：判断第一待处理栅格i在目标时刻t时的库容ResStroagei,t是否大于或等于其运行库容ComputeStoragei，是则跟据如下公式获取第一待处理栅格i在目标时刻t时的泄流量ResDrainagei,t，否则说明第一待处理栅格i在目标时刻t时，其上游来水全部拦蓄，则取第一待处理栅格i在目标时刻t时的泄流量ResDrainagei,t＝0；
ResDrainagei,t＝(ResStroagei,t-ComputeStoragei)*RDC
通过执行步骤3.1.1至3.1.2，获取第一待处理栅格对应目标时刻t的泄流量，进而获取各水库栅格分别对应目标时刻t的泄流量。
以屯溪流域为目标流域，根据专家经验和试验获取屯溪流域的RDC等于0.75。以2012年4月23日8时开始的一场洪水为例，获取的洪水起止时刻水库蓄量的空间变化量如图7所示。
在计算完第一待处理栅格i在目标时刻t时的泄流量后，根据如下公式：
ResStroage′i,t＝ResStroagei,t-ResDrainagei,t
获取第一待处理栅格i在目标时刻t时的更新库容ResStroage′i,t；
其中，ResDrainagei,t为第一待处理栅格i在目标时刻t时的泄流量，ResStroagei,t为第一待处理栅格i在目标时刻t时的库容。
基于MCT汇流原理，构建水库栅格及非水库栅格的汇流模型。
对于各非水库栅格，根据如下方法获取其泄流量：
分别针对演算次序等于1的各非水库栅格，执行如下操作：根据蓄满产流模型，获取各非水库栅格在目标时刻t时的产流量，将该产流量作为对应非水库栅格的泄流量；
针对演算次序不等于1的各非水库栅格，按照演算次序由小到大的顺序依次获取各非水库栅格的泄流量，即：先对演算次序等于2的各个栅格单元进行演算，再对演算次序等于3的各个栅格单元进行演算，以此类推直至对最大演算次序的栅格单元进行演算完毕，进而各演算次序不等于1的非水库栅格的泄流量。计算各演算次序不等于1的非水库栅格的泄流量的方法包括如下步骤：
分别针对演算次序不等于1的各非水库栅格，将其作为第二待处理栅格i'，根据如下公式：
Oi′,t＝C1Ii′,t+C2Ii′,t-1+C3Oi′,t-1
获取第二待处理栅格目标时刻t的泄流量Oi′,t，进而获取各非水库栅格分别对应目标时刻t的泄流量；
其中，C1、C2、C3为和第二待处理栅格i'的参考流量以及流速相关的系数，Ii',t为第二待处理栅格i'在目标时刻t时的入流量，Ii',t-1为第二待处理栅格i'在时刻t-1时的入流量,Oi′,t-1为第二待处理栅格i'时刻t-1时的泄流量。
根据如下公式：
Figure BDA0002993358280000121
Figure BDA0002993358280000122
Figure BDA0002993358280000123
获取和第二待处理栅格i'的参考流量以及流速相关系数C1、C2、C3；
其中，CCt-1、CCt、CDt、CDt-1为中间参数，根据如下公式：
Figure BDA0002993358280000124
Figure BDA0002993358280000125
Figure BDA0002993358280000126
Figure BDA0002993358280000127
获取中间系数CCt-1、CCt、CDt、CDt-1；
其中，Δt为预设的计算时段长；Δx为预设的计算空间步长，cvt-1和cvt分别为第二待处理栅格i'在t-1时刻和t时刻的流速，βt-1和βt分别为第二待处理栅格i'在t-1时刻和t时刻的预设校正因子；Bt-1和Bt为第二待处理栅格i'在t-1时刻和t时刻的所对应的河面宽度值，为方便举例，假设河道为矩形，河面宽度不随时间而变化，则Bt-1＝Bt，S0为与第二待处理栅格i'所对应的河底坡度；Qt-1和Qt分别为第二待处理栅格i'在t-1时刻和t时刻的参考流量。在分析中使用了虚拟河道的概念，河底坡度即为坡面栅格的坡度，对于非河道的坡面栅格，其河底坡度等于该栅格单元的坡度；坡面栅格的虚拟河道的河道宽度值是通过经验公式获取的，等于该栅格的汇流累计值的β次方与设定系数
Figure BDA0002993358280000134
的乘积，其中，β也是设定系数。
Figure BDA0002993358280000135
和β通过将河道栅格已知的河道宽度和汇流累计值建立回归关系得到。
根据如下公式：
Il,t-1＝Ol′,t-1+Rl,t-1
Il,t＝Ol′,t+Rl,t
获取目标流域中各栅格单元l在目标时刻t时的入流量ll,t，以及各栅格单元l在时刻t-1时的入流量ll,t-1；
其中，Ol′,t-1和Ol′,t分别为与各栅格单元l相邻的上游栅格单元l′在t-1时刻和t时刻的泄流量，如l′为水库栅格，即为ResDrainagel′,t-1和ResDrainagel′,t；Rl,t-1和Rl,t分别为根据蓄满产流模型计算得到的各栅格单元l在t-1时刻和t时刻的产流量。
根据如下公式：
O′t＝Ot-1+(It-It-1)
Qt-1＝(Ot-1+It-1)/2；Qt＝(O′t+It)/2
Figure BDA0002993358280000131
Figure BDA0002993358280000132
Figure BDA0002993358280000133
其中，n为曼宁糙率系数，通过试验得到；α为边坡夹角，假设屯溪流域为矩形河道，边坡夹角为90度；At-1和At为第二待处理栅格i'在t-1时刻和t时刻的过水断面面积，矩形河道时等于该时刻的河道水深乘以河面宽度，可通过牛顿迭代法求解；Pt-1和Pt为第二待处理栅格i'在t-1时刻和t时刻的所对应的河道湿周，约等于该时刻的河面宽度在矩形河道时，Pt-1＝Pt≈Bt-1＝Bt，其他符号意义同前所述。
获取第二待处理栅格i'在t-1时刻和t时刻的参考流量Qt-1和Qt，以及第二待处理栅格i'在t-1时刻和t时刻的流速cvt-1和cvt，以及第二待处理栅格i'在t-1时刻和t时刻的预设校正因子，βt-1和βt。
以屯溪流域为目标流域，以2012年4月23日8时开始的一场洪水为例，洪水模拟结果如图8所示。
本发明基于流域DEM数据计算出流域中所有栅格单元的坡度、流向和汇流累计量；基于水库经纬度计算得到水库栅格和非水库栅格；建立水库栅格的蓄泄模型；基于MCT方法，构建非水库栅格单元的汇流模型。本发明考虑了流域内水库对汇流计算的影响，数据来源稳定可靠，计算结果精确合理，同时解决了有水库的流域的汇流计算问题。且本发明基于水库的真实位置，在汇流演算中分布式的考虑了流域内有资料水库(有各特征库容数据)和无资料水库(仅有总库容数据)的作用，这样既保证了计算结果的精度与可靠性，同时解决了无资料水库较多流域的汇流计算问题。且本方法主要应用流域数字高程模型，数据来源稳定可靠，方法中变量之间的函数关系明确，有利于在计算机上编程实现，保证了结果的客观合理性，可以进一步促进数字水文学以及分布式模型的深入发展。
以上所述仅是本发明的优选实施方式，对于本技术领域的普通技术人员来说，在不脱离本发明技术原理的前提下，还可以做出若干改进和变形，这些改进和变形也应视为本发明的保护范围。

Claims (7)

1.一种全分布式考虑流域内水库影响的汇流模拟方法，其特征在于，用于获取目标流域内各栅格单元在目标时刻t时的泄流量，方法包括如下步骤：
步骤1：基于目标流域的数字高程数据，获取目标流域中各栅格单元的坡度、流向、汇流累计值；基于各栅格单元的流向和汇流累计值，计算各栅格单元的演算次序；然后进入步骤2；
计算演算次序的方法包括如下步骤：
步骤1.1：分别针对汇流累计值等于1的各栅格单元，将其演算次序赋值为1；分别针对汇流累计值不等于1的各栅格单元，初始化其演算次序n＝0；然后进入步骤1.2；
步骤1.2：分别针对各演算次序n＝0的栅格单元，将其作为待处理单元；同时分别针对各待处理单元，执行如下操作，然后进入步骤1.3；
根据与待处理单元相邻各栅格单元的流向，获取水流流向待处理单元的各相邻栅格单元：
若水流流向待处理单元的各栅格单元的演算次序均不等于零，则对水流流向待处理单元的各相邻栅格单元中的演算次序最大值加1获取一累计值，将该累计值作为待处理单元的演算次序；
若至少一个水流流向待处理单元的相邻栅格单元的演算次序等于零，则不改变待处理单元的演算次序；
步骤1.3：判断是否存在演算次序n＝0的栅格单元，是则返回步骤1.2；否则停止循环；
步骤2：根据目标流域内各水库的位置，将各栅格单元区分为水库栅格和非水库栅格；然后进入步骤3；
步骤3：分别针对各水库栅格i，将其作为第一待处理栅格，依次执行步骤3.1.1至步骤3.1.2：
步骤3.1.1：判断第一待处理栅格i在目标时刻t时的库容ResStroagei，t是否大于或等于其总库容
Figure FDA0003726149500000011
是则根据如下公式获取第一待处理栅格i在目标时刻t时的泄流量ResDrainagei，t；否则进入步骤3.1.2；
ResDrainagei，t＝(TotalStoragei-ComputeStoragei)*RDC+(ResStroagei，t-TotalStoragei)
其中：ComputeStoragei为第一待处理栅格i的运行库容；RDC为目标流域的水库泄流系数；
步骤3.1.2：判断第一待处理栅格i在目标时刻t时的库容ResStroagei，t是否大于或等于其运行库容ComputeStoragei，是则跟据如下公式获取第一待处理栅格i在目标时刻t时的泄流量ResDrainagei，t，否则取第一待处理栅格i在目标时刻t时的泄流量ResDrainagei，t＝0；
ResDrainagei，t＝(ResStroagei，t-ComputeStoragei)*RDC
通过执行步骤3.1.1至3.1.2，获取第一待处理栅格对应目标时刻t的泄流量，进而获取各水库栅格分别对应目标时刻t的泄流量；
分别针对演算次序等于1的各非水库栅格，执行如下操作：根据产流模型，获取各非水库栅格在目标时刻t时的产流量，将该产流量作为对应非水库栅格的泄流量；
针对演算次序不等于1的各非水库栅格，按照演算次序由小到大的顺序依次获取各非水库栅格的泄流量，获取对应非水库栅格的泄流量的方法包括如下步骤：
分别针对演算次序不等于1的各非水库栅格，将其作为第二待处理栅格i'，根据如下公式：
Oi′，t＝C1Ii′，t+C2Ii′，t-1+C3Oi′，t-1
获取第二待处理栅格目标时刻t的泄流量Oi′，t，进而获取各非水库栅格分别对应目标时刻t的泄流量；
其中，C1、C2、C3为和第二待处理栅格i'的参考流量以及流速相关的系数，Ii'，t为第二待处理栅格i'目标时刻t时的入流量，Ii'，t-1为第二待处理栅格i'时刻t-1时的入流量，Oi′，t-1为第二待处理栅格i'时刻t-1时的泄流量。
2.根据权利要求1所述的全分布式考虑流域内水库影响的汇流模拟方法，其特征在于，在步骤1中，根据D8流向法获取各栅格单元的流向。
3.根据权利要求1所述的全分布式考虑流域内水库影响的汇流模拟方法，其特征在于，在步骤3.1.1中：在初始时刻时，第一待处理栅格i的库容ResStroagei，0＝computeStoragei；
在初始时刻之后的目标时刻t时，根据如下公式：
ResStroagei，t＝ResStroagei′，t-1+Ii，t
获取第一待处理栅格i在目标时刻t时的库容ResStroagei，t；
其中，Ii，t为第一待处理栅格i在目标时刻t时的入流量，ResStroagei′，t-1第一待处理栅格i在目标时刻t-1时的更新库容。
4.根据权利要求3所述的全分布式考虑流域内水库影响的汇流模拟方法，其特征在于，根据如下公式：
ResStroage′i，t-1＝ResStroagei，t-1-ResDrainagei，t-1
获取第一待处理栅格i在目标时刻t-1时的更新库容ResStroage′i，t-1；
其中，ResDrainagei，t-1为第一待处理栅格i在目标时刻t-1时的泄流量，ResStroagei，t-1为第一待处理栅格i在目标时刻t-1时的库容。
5.根据权利要求3所述的全分布式考虑流域内水库影响的汇流模拟方法，其特征在于，根据如下公式：
Il，t-1＝Ol′，t-1+Rl，t-1
Il，t＝O′l，t+Rl，t
获取目标流域中各栅格单元l在目标时刻t时的入流量ll，t，以及各栅格单元l在时刻t-1时的入流量ll，t-1；
其中，Ol′，t-1和Ol′，t分别为与各栅格单元l相邻的上游栅格单元l′在t-1时刻和t时刻的泄流量，如l′为水库栅格，即为ResDrainagei′，t-1和ResDrainagei′，t；Rl，t-1和Rl，t分别为根据产流模型计算得到的各栅格单元l在t-1时刻和t时刻的产流量。
6.根据权利要求1所述的全分布式考虑流域内水库影响的汇流模拟方法，其特征在于，在步骤3.1.1中，根据如下方法获取第一待处理栅格i的运行库容ComputeStoragei：
若和第一待处理栅格i所对应水库的相关的资料中记载了其兴利库容ActiveStoragei和死库容DeadStoragei，根据如下公式：
ComputeStoragei＝ActiveStoragei+DeadStoragei
获取第一待处理栅格i的运行库容ComputeStoragei；
若和第一待处理栅格i所对应水库相关的资料中未记载其兴利库容ActiveStoragei和死库容DeadStoragei，根据如下公式：
Figure FDA0003726149500000041
获取第一待处理栅格i的运行库容ComputeStoragei；
其中，根据如下公式：
Figure FDA0003726149500000042
Figure FDA0003726149500000043
获取目标流域内的平均库容系数
Figure FDA0003726149500000044
其中，N为目标流域内兴利库容和死库容均已知的水库的总数量，CSCn为目标流域内兴利库容和死库容均已知的第n个水库的库容系数，ActiveStoragen为目标流域内兴利库容和死库容均已知的第n个水库的兴利库容；DeadStoragen为与目标流域内兴利库容和死库容均已知的第n个水库的死库容；TotalStoragen为与目标流域内兴利库容和死库容均已知的第n个水库的总库容。
7.根据权利要求1所述的全分布式考虑流域内水库影响的汇流模拟方法，其特征在于，在步骤3.1.2中，根据如下公式：
Figure FDA0003726149500000045
Figure FDA0003726149500000046
Figure FDA0003726149500000047
获取和第二待处理栅格i′的参考流量以及流速相关系数C1、C2、C3；
其中，CCt-1、CCt、CDt、CDt-1为中间参数，根据如下公式
Figure FDA0003726149500000051
Figure FDA0003726149500000052
Figure FDA0003726149500000053
Figure FDA0003726149500000054
获取中间系数CCt-1、CCt、CDt、CDt-1；
其中，Δt为预设的计算时段长；Δx为预设的计算空间步长，cvt-1和cvt分别为第二待处理栅格i′在t-1时刻和t时刻的流速，βt-1和βt分别为第二待处理栅格i′在t-1时刻和t时刻的预设校正因子；Bt-1和Bt为与第二待处理栅格i′在t-1时刻和t时刻的所对应的河面宽度值，S0为与第二待处理栅格i′所对应的河底坡度；Qt-1和Qt分别为第二待处理栅格i′在t-1时刻和t时刻的参考流量。
CN202110322305.4A 2021-03-25 2021-03-25 一种全分布式考虑流域内水库影响的汇流模拟方法 Active CN113190958B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110322305.4A CN113190958B (zh) 2021-03-25 2021-03-25 一种全分布式考虑流域内水库影响的汇流模拟方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110322305.4A CN113190958B (zh) 2021-03-25 2021-03-25 一种全分布式考虑流域内水库影响的汇流模拟方法

Publications (2)

Publication Number Publication Date
CN113190958A CN113190958A (zh) 2021-07-30
CN113190958B true CN113190958B (zh) 2022-08-12

Family

ID=76973909

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110322305.4A Active CN113190958B (zh) 2021-03-25 2021-03-25 一种全分布式考虑流域内水库影响的汇流模拟方法

Country Status (1)

Country Link
CN (1) CN113190958B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114840989B (zh) * 2022-04-22 2022-11-11 河海大学 一种栅格尺度考虑水利工程拦蓄的河道汇流演算方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108643116A (zh) * 2018-05-08 2018-10-12 河海大学 一种山区性水库库区河道水面宽度的估算方法
CN110096751A (zh) * 2019-04-03 2019-08-06 河海大学 一种估算无资料地区中小型水库蓄量的方法
CN110543692A (zh) * 2019-08-07 2019-12-06 河海大学 一种基于下垫面特征的重配置汇流模拟方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108643116A (zh) * 2018-05-08 2018-10-12 河海大学 一种山区性水库库区河道水面宽度的估算方法
CN110096751A (zh) * 2019-04-03 2019-08-06 河海大学 一种估算无资料地区中小型水库蓄量的方法
CN110543692A (zh) * 2019-08-07 2019-12-06 河海大学 一种基于下垫面特征的重配置汇流模拟方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
Flash Flood Simulation for Ungauged Catchments Based on the Distributed Hydrological Model;Pengfei Jia等;《Water》;20190104;全文 *
Multiple hydrological models comparison and an improved Bayesian model averaging approach for ensemble prediction over semi-humid regions;Wenbo Huo等;《Stochastic Environmental Research and Risk Assessment》;20180829;全文 *
基于分布式水文模型的小流域山洪预报方法与应用;包红军等;《暴雨灾害》;20170430;第36卷(第2期);全文 *

Also Published As

Publication number Publication date
CN113190958A (zh) 2021-07-30

Similar Documents

Publication Publication Date Title
CN106884405B (zh) 一种无资料地区溃堤型山洪灾害分析评价方法
Tayfur Modern optimization methods in water resources planning, engineering and management
CN111046563B (zh) 一种梯级水库群连续溃决洪水模拟方法
CN113204927B (zh) 一种针对城市降雨洪涝过程的不同尺度分级嵌套模拟方法
CN113610264A (zh) 精细化电网台风洪涝灾害预测模型
CN113190958B (zh) 一种全分布式考虑流域内水库影响的汇流模拟方法
CN112052545B (zh) 一种基于元胞自动机的城市地表径流与管网汇流耦合方法
US20220186451A1 (en) Remote sensing-based dynamic estimation method for outflow process of ungauged free overflow reservoir
Xu et al. Four-decades of bed elevation changes in the heavily regulated upper Atchafalaya River, Louisiana, USA
Meema et al. Real-time optimization of a large-scale reservoir operation in Thailand using adaptive inflow prediction with medium-range ensemble precipitation forecasts
Chaudhary et al. Integrated 1D and 2D numerical model simulations for flushing of sediment from reservoirs
CN110377868B (zh) 一种基于实时雨情的动态水系提取方法
Wu et al. Simulation of lake-groundwater interaction under unsteady-state flow using the sloping lakebed method
Lee Implicit discontinuous Galerkin scheme for discontinuous bathymetry in shallow water equations
Idrees et al. Complementary modeling approach for estimating sedimentation and hydraulic flushing parameters using artificial neural networks and RESCON2 model
CN112926786B (zh) 一种基于关联规则模型和数值模拟的浅水湖泊目标水位逆向预测方法及系统
Vedula et al. Real-time multipurpose reservoir operation: a case study
Chen et al. The development of a real-time flooding operation model in the Tseng-Wen Reservoir
Javed et al. Sediment flushing strategy for reservoir of proposed Bhasha Dam, Pakistan
Huang et al. Optimal design of interception for flood control: an integrated simulation approach
Konan Urban Flood modelling and Floodplain Mapping using ArcGIS, HEC-HMS and HEC-RAS in Abidjan city, Côte D’Ivoire–West Africa: Case study of the watershed of Bonoumin-Rivièra Palmeraie
Singh Managed flood storage option for selected levees along the lower Illinois River for enhancing flood protection, agriculture, wetlands, and recreation: First report, stage and flood frequencies and the Mississippi backwater effects
Metler et al. Effect of wind waves and wind tides on the optimum control of large lakes
Mahon et al. Water resources development alternatives for the Mississippi Alluvial Plain in Eastern Arkansas
Kafle Numerical Simulation of Backwater Effects by a Downstream Dam using HEC-RAS: A case of SunKoshi-Marin Diversion Headworks, Nepal

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