留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

双层异步迭代洪水演进模拟算法

刘修国 王玉着 刘旭东 高伟 张唯

刘修国, 王玉着, 刘旭东, 高伟, 张唯. 双层异步迭代洪水演进模拟算法[J]. 武汉大学学报 ● 信息科学版, 2016, 41(12): 1570-1576,1612. doi: 10.13203/j.whugis20140710
引用本文: 刘修国, 王玉着, 刘旭东, 高伟, 张唯. 双层异步迭代洪水演进模拟算法[J]. 武汉大学学报 ● 信息科学版, 2016, 41(12): 1570-1576,1612. doi: 10.13203/j.whugis20140710
LIU Xiuguo, WANG Yuzhuo, LIU Xudong, GAO Wei, ZHANG Wei. Double-layer Asynchronous Iterative Algorithm for Flood Routing Simulation[J]. Geomatics and Information Science of Wuhan University, 2016, 41(12): 1570-1576,1612. doi: 10.13203/j.whugis20140710
Citation: LIU Xiuguo, WANG Yuzhuo, LIU Xudong, GAO Wei, ZHANG Wei. Double-layer Asynchronous Iterative Algorithm for Flood Routing Simulation[J]. Geomatics and Information Science of Wuhan University, 2016, 41(12): 1570-1576,1612. doi: 10.13203/j.whugis20140710

双层异步迭代洪水演进模拟算法

doi: 10.13203/j.whugis20140710
基金项目: 

国家自然科学基金 No.41501584

详细信息
    作者简介:

    刘修国,博士,教授,主要从事GIS技术及其应用工程方向研究。liuxg318@163.com

    通讯作者: 张唯,博士。zw_paper@163.com
  • 中图分类号: P237.9;P333.2

Double-layer Asynchronous Iterative Algorithm for Flood Routing Simulation

Funds: 

The National Natural Science Foundation of China No.41501584

More Information
    Author Bio:

    LIU Xiuguo, PhD, professor, specializes in GIS and its application. E-mail: liuxg318@163.com

    Corresponding author: ZHANG Wei,PhD. E-mail:。zw_paper@163.com
图(7) / 表(1)
计量
  • 文章访问数:  977
  • HTML全文浏览量:  39
  • PDF下载量:  361
  • 被引次数: 0
出版历程
  • 收稿日期:  2015-02-10
  • 刊出日期:  2016-12-05

双层异步迭代洪水演进模拟算法

doi: 10.13203/j.whugis20140710
    基金项目:

    国家自然科学基金 No.41501584

    作者简介:

    刘修国,博士,教授,主要从事GIS技术及其应用工程方向研究。liuxg318@163.com

    通讯作者: 张唯,博士。zw_paper@163.com
  • 中图分类号: P237.9;P333.2

摘要: 均一化迭代是基于栅格数据的洪水演进模拟常采用的迭代策略,存在设定迭代步长与栅格单元实际所需时间不一致的问题,直接影响到栅格单元水位、流量的计算精度。针对该问题,提出了一种双层异步迭代策略。外层迭代控制洪水演进时刻,内层迭代则通过分析不同栅格单元洪水的流速特征,实现迭代步长及迭代次数的自适应设定。采用商用FloodArea软件、均一化迭代算法及双层异步迭代算法,分别对福建省万安流域暴雨洪水历史数据进行模拟。该算法模拟结果的平均误差比FloodArea小0.361 m,比均一化迭代算法小0.654 m,说明该算法能有效提高洪水演进模拟的整体精度。

English Abstract

刘修国, 王玉着, 刘旭东, 高伟, 张唯. 双层异步迭代洪水演进模拟算法[J]. 武汉大学学报 ● 信息科学版, 2016, 41(12): 1570-1576,1612. doi: 10.13203/j.whugis20140710
引用本文: 刘修国, 王玉着, 刘旭东, 高伟, 张唯. 双层异步迭代洪水演进模拟算法[J]. 武汉大学学报 ● 信息科学版, 2016, 41(12): 1570-1576,1612. doi: 10.13203/j.whugis20140710
LIU Xiuguo, WANG Yuzhuo, LIU Xudong, GAO Wei, ZHANG Wei. Double-layer Asynchronous Iterative Algorithm for Flood Routing Simulation[J]. Geomatics and Information Science of Wuhan University, 2016, 41(12): 1570-1576,1612. doi: 10.13203/j.whugis20140710
Citation: LIU Xiuguo, WANG Yuzhuo, LIU Xudong, GAO Wei, ZHANG Wei. Double-layer Asynchronous Iterative Algorithm for Flood Routing Simulation[J]. Geomatics and Information Science of Wuhan University, 2016, 41(12): 1570-1576,1612. doi: 10.13203/j.whugis20140710
  • 近年来,随着全球气候的变化,洪涝灾害频频发生,对人类的生产生活造成了极大的危害[1, 2]。结合实际地形信息,模拟洪水的演进过程,进而动态预测灾害的影响范围,是当前资源环境和气象灾害领域的研究热点之一[3]。洪水演进过程的模拟主要有水文模型法、数值模拟法以及水动力学建模法等[4-6]。前两类方法考虑的因素较为全面,通过复杂精细的数学模型来重建洪水演进的全过程,但其中涉及的水文专业信息较多,数据的获取和建模均有一定难度。水动力学建模主要考虑地形和水流因素,模拟水流在不同位置的流动状态,实现简单。研究者设计实现了大量优秀的水动力学建模算法。如Bates等[7]针对河道内与河漫滩洪水演进的不同特征,分别利用一维运动波及二维扩散波对模型进行简化。Hagen等[8]结合阿富汗的地形特征提出了一种简化的水动力学洪水演进模型。Gemmer[9]对基于GIS的水动力模型—FloodArea模型模拟过程的水动力参数进行修改,使其更好地模拟洪水淹没风险。孙海等[10]基于DEM提出“环形”算法,模拟溃坝洪水演进过程等。由此可见,以栅格地形为基础数据,采用圣维南方程模拟洪水动力,通过不同栅格像元的迭代实现水流过程的模拟,是水动力学洪水演进最为常见的方法之一。但现有的研究主要集中在针对地形特征的简化和圣维南方程的解算方法上,对于演进过程中不同栅格单元水流变化的迭代过程研究则较少。

    基于栅格单元的迭代策略可分为均一式和异步式两种。前者假设每个栅格单元在相同的时间间隔内完成迭代过程,是水文分析常见的简化方式。但不同的栅格单元往往具有不同的坡降、水深与流速,每个栅格单元相对应的迭代步长应有所不同。均一化迭代忽略了这一现象,使得模拟结果存在较大误差,复杂地形区域的误差尤为明显。后者充分考虑了每个栅格单元的自身特征,努力寻求迭代次数和间隔时间的平衡,因此更符合洪水演进的真实情况。但目前尚未见到具体迭代策略的实现方式。

    为真实模拟洪水演进过程,本文提出了基于双层异步迭代策略的洪水演进模拟算法。利用福建省万安流域2010年6月18日的暴雨洪水历史数据进行实验,分别用气象行业较为通用的FloodArea商业软件、均一化迭代算法以及本文提出的双层异步迭代算法对不同时间段的洪水淹没过程进行模拟。

    • 洪水数据采用与DEM同分辨率的栅格数据来表示,栅格单元的值表示该栅格单元内水深的平均值。以栅格单元为单位,利用曼宁公式[11]计算中心栅格单元流向邻域栅格单元的流速:

      (1)

      式中,V为栅格断面流速(m/s);Rh为水力半径(m),用栅格单元所对应的水深值来表示;I为最大水面比降,n为曼宁系数。

      计算中心栅格流向采用D8算法[12],通过计算中心栅格单元与周围相邻8个栅格单元的水面比降,取最大者所在的栅格为水流向的栅格。水面比降计算如式(2):

      (2)

      式中,Ii为水面比降;h、hi分别为中心栅格与邻域栅格的自由液面高度(m),即栅格高程值与水深值之和;di为栅格距离(m),定义为:

      (3)

      确定流向以及流速后,即可计算流向邻近栅格的水量:

      (4)

      式中,Q为流向邻近栅格水量(m3);V为栅格流速(m/s);Δt为迭代步长(s)A为栅格截面面积(m2)。由于栅格之间的过流断面是矩形,断面面积可通过式(5)计算:

      (5)

      式中,d表示栅格距离(m);hflow表示两个栅格间的断面水深(m),可通过Gemmer[9]提出的方法计算。

      根据水量Q及栅格单元面积S,通过式(6)最终可得到结果淹没水深H

      (6)
    • 均一化迭代法[9],是指在每一次迭代计算中,赋予所有栅格相同的迭代步长,计算其流量、水深等信息。具体迭代过程为:① 利用DEM数据获得水文特征参数信息,如水面比降、水流方向等;② 根据洪水数据初始化淹没水深矩阵;③ 设置迭代步长Δt;④ 根据式(4)、式(6)迭代计算所有栅格的断面流量Q进而计算栅格水深HZ等信息;⑤ 一次迭代计算完毕,将Δt添加到累计时长Tsum中;⑥ 当Tsum小于总模拟时长T时,重复步骤③~⑤;否则迭代计算结束。此时的栅格水深数据HZ即为淹没模拟时长T的洪水演进淹没结果,栅格值即淹没水深,均一化迭代过程如图 1所示。

      图  1  均一化迭代过程

      Figure 1.  Uniform Iteration Process

      DEM中每个栅格单元的最大水面比降不同,其对应的水流速度及流向邻域单元所需时间各异。而均一化迭代过程中所有栅格的迭代步长相同,实际所需时间不一致,影响栅格单元水位、流量的计算精度。统计2010年6月18日3时福建省将乐县万安流域DEM中各中心栅格单元流向邻域单元所需的时间,耗时分布如图 2所示。

      图  2  流向邻近栅格单元所需时间统计图

      Figure 2.  Statistical Chart of Required Time for Flowing to Adjacent Grid

      图 2可知,中心栅格单元流向邻域单元的所需时间主要集中在5~129 s之间,少数分布在129~253 s之间。在较大跨度的时间区域内,均一化迭代法降低了洪水演进模拟的精度。为解决这一问题,本文提出了双层异步迭代算法。通过外层迭代控制洪水演进时刻,内层异步迭代实现迭代步长及迭代次数的自适应设定,动态模拟洪水的演进过程。

    • 双层异步迭代洪水演进模拟算法包含外层迭代与内层异步迭代。外层迭代继承了均一化迭代算法的思想,整幅栅格数据设定相同的迭代步长,当所有栅格单元计算完毕后进行下一次迭代计算,直到累计迭代步长大于等于洪水演进的模拟时长为止。外层迭代的主要目的是通过给定所有栅格数据一个相同的迭代步长来控制洪水演进的模拟时刻,本节称此处的迭代步长为全局迭代步长Δt

      内层异步迭代主要针对均一化迭代中对每个栅格单元计算断面流量与水深这一过程进行分析,通过分析中心栅格单元流向邻域栅格单元所需的时间Tflow与栅格迭代步长ΔTcell的关系,进行多次断面流量Q及水深H的迭代计算。

      ΔTcell是指在内层迭代中每个栅格单元自身的迭代步长,由两部分构成。一是来自外层迭代的全局迭代步长Δt;二是来自前一次迭代计算后该栅格单元余下的滞留迭代步长ΔTpause。ΔTpause是指当水流向邻域单元所需的时间大于给定的迭代步长时,水在该单元产生滞留,同时迭代步长累计到该栅格单元中,并参与下次迭代计算,称此时的迭代步长为滞留迭代步长。栅格单元不同,其ΔTpause也往往不同,因此ΔTcell各异。

      内层异步迭代过程的异步特性主要体现在两个方面,一是不同的栅格单元往往具有不同的ΔTcell;二是由于ΔTcell及流速的不同,导致各个栅格单元的迭代计算次数亦不同。异步迭代过程的目标是建立迭代次数的自适应选择机制,使得大部分栅格单元的迭代步长与流向下一个栅格单元所需的时间能保持一致。

      异步迭代过程是由一个栅格单元开始逐步向多个栅格单元计算发展的多次迭代过程。在内层迭代过程中,采用栅格单元缓冲池队列存储中心单元及其在ΔTcell内水量可流出的邻域单元。中心单元的ΔTcell由全局迭代步长Δt及该单元前一次迭代计算的ΔTpause构成,而邻域单元的ΔTcell由流向邻域单元的实际时间Tflow及自身ΔTpause构成。

      根据ΔTcellTflow的关系,内层迭代步长可概括为以下3种情况。

      1) ΔTcell<0.5Tflow时,即内层迭代步长低于实际水流时间的一半,则该单元水量滞留,ΔTcell累加到滞留迭代步长ΔTpause中;

      2) 0.5Tflow≤ΔTcellTflow时,即内层迭代步长大于等于实际水流时间一半且小于等于实际水流时间,则调整ΔTpause为0,Tflow为ΔTcell。根据式(4)计算Tflow步长内该中心单元流向邻域单元的水量Q,进而利用式(6)计算水量Q对应的水深H,同步更新中心单元及邻域单元淹没水深数据HZ

      3) ΔTcellTflow时,即内层迭代步长大于实际水流时间。此时直接将ΔTpause设为0且根据式(4)计算Tflow步长内该中心单元流向邻域单元的水量Q,进而利用式(6)计算水量Q对应的水深H,同步更新中心单元及邻域单元水量。

      为了兼顾模型精度与运算速度,要合理选取全局迭代步长Δt。双层异步迭代法中,Δt的选取准则是使尽量多的栅格单元在一次迭代计算出现异步情况2)。由于栅格所需时间基本符合正态分布形态,Δt可以通过式(7)计算:

      (7)

      式中,Pi表示Tflowi出现的概率;Tflowi∈(μ-2σ,μ+2σ);μ表示Tflowi的期望;σ表示方差。

    • 实现计算模拟总时长为T的洪水演进过程,算法流程如下。

      1) 输入原始高程数据DEM(DZ)、水深数据WZ及总模拟时长T,初始化累计时长Tsum为0,根据数据WZ初始化淹没水深数据HZ

      2) 设置全局迭代步长Δt

      3) 执行内层迭代计算。依次选取栅格单元i,清空当前缓冲池。设置内层迭代步长ΔTcell为全局迭代步长Δt,将中心单元i添加到缓冲池中;

      4) 若ΔTcell小于等于0,返回步骤3);若ΔTcell大于0,对缓冲池中的栅格单元执行淹没操作,如步骤5)~7);

      5) 从缓冲池中获得要处理的栅格单元。若缓冲池为空,则调整内层迭代步长ΔTcell,返回步骤4);否则根据D8算法计算流向Dir;

      6) 若Dir小于0,返回步骤5);若Dir大于等于0,通过式(1)计算栅格断面流速V及流向邻域时间Tflow。比较Tflow与内层ΔTcell的关系,明确该单元属于哪种异步迭代情况;

      7) 若该栅格单元对应的Tflow与ΔTcell满足情况2)或3),则根据式(4)计算其流量Q、断面面积A,根据A、VTflow,计算流向中心及邻域单元的水量,同步更新水深数据HZ;若满足情况1),该单元水量滞留,则直接返回步骤5);

      8) 每次内层迭代结束后,将Δt累加到累计步长Tsum中。若TsumT,则返回步骤2)继续下一次内层迭代;反之,则异步迭代结束,此时得到的水深数据HZ即为模拟总时长为T的洪水淹没结果数据。

      双层异步迭代算法流程如图 3所示,其中异步迭代过程如图 3右侧所示。

      图  3  双层异步迭代算法流程图

      Figure 3.  Flowchart of Double-layer Asynchronous Iteration Algorithm

      水量大小与地形复杂度决定迭代算法的时间消耗,水量越大,地形越复杂,其时间消耗越多。均一化迭代算法模型简单,迭代步长固定,减少了迭代次数,但其模拟结果不够精确。双层异步迭代算法增加了一层对全局迭代步长的循环计算,其耗时相较均一化迭代略有增加,但洪水演进模拟结果更加精确。

      在实际洪水应急决策中,快速计算出淹没结果的关键在于设置双层异步迭代算法中的全局迭代步长Δt。在外层控制迭代中,Δt越小,计算所得的结果越精确,但其耗时越长。Δt取值应该位于tmintmax之间,减少内层迭代的次数,提高运行效率。其中,tmax表示最大概率时间Tflow(p),tmintmax/2。根据某时刻DEM中各中心栅格单元Tflow可获取其耗时数据集合,最大概率时间Tflow(p)则为耗时数据集合中出现次数最多的栅格单元的Tflow

    • 2010年6月18日福建省将乐县万安镇遭遇特大暴雨袭击,暴雨造成特大洪水。农田、道路、电力、通讯等基础设施毁坏极其严重,农民遭受重大损失。选择1:5万的万安流域DEM以及6月18日的降雨过程数据作为模拟洪水演进过程样本。图 4为万安流域DEM,图 5为4个气象站所监测的降雨过程线,图 5中的代号为气象站点的编号。

      图  4  万安流域DEM

      Figure 4.  DEM Data of Wanan Valley

      图  5  万安2010年6月18日降雨过程线

      Figure 5.  Rainfall Line of Wanan on June 18th,2010

      目前业内大都采用商业软件FloodArea进行洪水演进模拟,效果较好[12]。实验拟采用均一化迭代法、FloodArea软件及本文双层异步迭代法分别模拟降雨洪水演进过程,并通过分析实际灾情考察点的淹没水深与洪水淹没结果,验证双层异步迭代算法的模拟效果。采用实际灾情平均误差表示模拟精度,即误差越小精度越高。平均误差按式(8)进行计算:

      (8)

      式中,hi表示第i个灾情考察点的模拟水深;hreali表示第i个灾情考察点的实际水深;n表示灾情考察点的个数。

    • 采用上述3种算法模拟从18日3时~16时共计13 h的洪水演进过程。提取7个典型受灾点在模拟中计算得到的最大淹没水深,将其与实际淹没最大水深对比,验证算法模拟精度,结果如图 6所示。受灾点的选取主要遵循流域分区控制原则与地形控制原则[13],优先选择受洪水威胁的中小河流并按小流域分区综合布点;并考虑山区降雨受地形抬升作用的影响,选取的受灾点应能反映时段降雨量等值线的转折变化。

      图  6  万安流域调查受灾点与模拟淹没水深比较

      Figure 6.  Investigative Disaster Point and Simulated Inundation Water Depth

      利用式(8)计算3个算法的平均误差,见表 1。由表 1可知,双层异步迭代法与灾情点实际水深平均误差最小,其模拟精度最高。但由于1:5万的DEM数据对中小流域地形表现力有限,且未考虑蒸发渗透等水文模型以及人工建筑对洪水走向改变的影响等因素,因此3个算法的模拟误差均较大。

      表 1  三种模拟方法的平均误差/m

      Table 1.  Average Error of Three Models Results/m

      均一化迭代法 FloodArea 异步迭代法
      平均误差 1.014 0.721 0.360

      在洪水淹没形态分布方面,分别选取降雨2 h、6 h、10 h、13 h共4个时次的模拟结果分析比较,如图 7所示。

      图  7  降雨2 h、6 h、10 h和13 h淹没模拟对比

      Figure 7.  Inundation Simulation Contrast After Raining in 2,6,10 and 13 Hours

      通过分析3个算法淹没水深的空间分布可以看出,降雨2 h后初始水量较小,地面汇水少,模拟结果淹没状况较为一致。随着降雨过程的演进,地面水量增大,双层异步迭代算法地形敏感度较强,水流向地势低洼地方汇集。降雨6 h后,双层异步迭代法淹没结果在0.5~1 m的淹没范围较大,降雨10 h、13 h后,双层异步迭代算法与FloodArea软件的淹没结果较为一致,而均一化迭代算法由于每次迭代过程中每栅格单元最多流向一个邻域,地形坡降较小时会出现水流滞留,产生较小的“河网分叉”[14]。可见,双层异步迭代算法的模拟结果在淹没形态上与其他两个算法基本一致,但与实际灾情点淹没水深更接近,准确性更高。

    • 本文提出的双层异步迭代算法较好地解决了均一化迭代算法中迭代步长与单个栅格单元水流演进时间不一致的问题。通过外层迭代控制洪水演进的时间尺度,内层迭代计算水量信息,实现迭代次数的自适应选择。实验表明,双层异步迭代算法可以有效地模拟自然地形下的洪水演进过程,能实时获取洪水的淹没范围和水深分布信息,为防汛指挥、风险图制作、风险评估等提供有力的依据。本文仅对迭代策略进行改进,未考虑蒸发渗透、下垫面类型等因素对洪水演进过程的影响,引入这些因素以完善模型是下一步待开展的研究工作。

参考文献 (14)

目录

    /

    返回文章
    返回