《洪水演进数值模拟计算中网格流出率的修正方法.pdf》由会员分享,可在线阅读,更多相关《洪水演进数值模拟计算中网格流出率的修正方法.pdf(13页珍藏版)》请在专利查询网上搜索。
1、(10)申请公布号 CN 103559403 A (43)申请公布日 2014.02.05 CN 103559403 A (21)申请号 201310556360.5 (22)申请日 2013.11.11 G06F 19/00(2011.01) (71)申请人 张驰 地址 213000 江苏省常州市新北区天山花园 2 幢乙单元 1103 室 (72)发明人 张驰 (74)专利代理机构 南京同泽专利事务所 ( 特殊 普通合伙 ) 32245 代理人 蒋全强 (54) 发明名称 洪水演进数值模拟计算中网格流出率的修正 方法 (57) 摘要 本发明涉及一种洪水演进数值模拟计算中网 格流出率的修正方法。
2、, 包括 : 通过修正系数对负 水深网格进行修正, 获得该网格的实际水深, 克服 了在以往的移动边界计算中, 经常发生水面标高 低于地面标高的情况, 也就是计算后的水深变为 负值, 即出现负水深网格, 从而导致计算域内的质 量不能保证守恒, 稳定性变差, 甚至计算发散而得 不到结果, 直接影响到洪水演进数值计算无法顺 利进行的技术问题, 通过引入所述修正系数 , 使当出现水面标高低于地面标高时, 洪水演进数 值模拟计算得已顺利进行。 (51)Int.Cl. 权利要求书 2 页 说明书 8 页 附图 2 页 (19)中华人民共和国国家知识产权局 (12)发明专利申请 权利要求书2页 说明书8页 。
3、附图2页 (10)申请公布号 CN 103559403 A CN 103559403 A 1/2 页 2 1. 一种洪水演进数值模拟计算中网格流出率的修正方法, 包括 : 通过洪水高发地区的地形图建立一计算域内的该洪水高发地区的数字高程模型, 并 对所述计算域进行二维规则网格划分, 以建立二维规则网格模型及进行网格参数设定, 建 立所述二维规则网格模型及网格参数设定的方法包括 : A : 建立所述二维规则网格模型, 即所述计算域按照一定空间步长 (X、 Y) 进行分割 后得到二维规则网格模型, 定义 (i, j) 为该二维规则网格模型中一网格, 且该网格为正方形 网格, 即 X=Y ; B :。
4、 网格参数设定, 即设定所述二维规则网格模型中任一网格的水深为 H, 且与该网格 按 X 轴同方向的单宽流量 M, 与 Y 轴同方向的单宽流量 N ; 同时设定所述洪水演进数值计算 的时间步 n, 以及该时间步 n 对应的时间步长 t ; 通过单宽流量公式获得所述二维规则网格模型中各网格在 n+1 时间步的所述单宽 流量 M 和单宽流量 N ; 其方法, 包括 : 预设初始条件, 即, 在 n 时间步, 且设定位于所述计算域的边界上的流量输入起始界的 各边界网格的初值参数, 该初值参数包括 : 初始水深Hn, 所述单宽流量M的对应流速矢量Un, 以及所述单宽流量 N 的对应流速矢量 Vn; 在。
5、所述单宽流量公式中位于所述计算域外的网格 的初值参数与该网格相邻的一边界网格的初值参数相同, 且位于所述计算域内的网格的初 值参数的相应取值为 0 ; 以及, 所述单宽流量 M 和单宽流量 N 的初始值的计算公式 : Mn=UnHn (1) ; Nn=VnHn (2) ; 所述单宽流量公式 : 根据所述各网格在n+1时间步的单宽流量M和单宽流量N, 建立所述各网格的水深计 算式, 即 其中, n+2 时间步作为水深计算的起始时间步, H 表示所述网格水深, 表示在所述起始时间步时网格 (i, j) 的区域内的水深 ;表示在 n 时间步的 所述网格的区域内的水深, 且该时间步的水深为预设值 ;表。
6、示在 n+1 时间步的所述 网格在X轴方向上流体的流入单宽流量 ;表示相邻网格在X轴方向上的流入单宽流 权 利 要 求 书 CN 103559403 A 2 2/2 页 3 量, 即, 所述网格在 X 轴方向上的流出单宽流量 ;表示在 n+1 时间步的所述网格在 Y 轴方向上流体的流入单宽流量 ;表示相邻网格在 Y 轴方向上的流入单宽流量, 即, 所 述网格在 Y 轴方向上的流出单宽流量 ; 若步骤计算出一网格的水深 H 小于 0, 即该网格为负水深网格, 则引入一修正系数 重新对该负水深网格的水量进行修正计算, 直到该负水深网格修正为实际水深后, 再转 入步骤通过所述单宽流量公式获得所述二维。
7、规则网格模型中各网格在下一时间步的单 宽流量 M 和单宽流量 N ; 其中, 当一网格出现负水深时, 根据负水深网格在当前时间步n的单宽流量M和单宽流 量 N 计算出该负水深网格的流入总量和流出总量 所述流入总量公式 : 所述流出总量公式 : 所述计算修正系数 的公式 : 并根据公式 (1) 、(2) 、(3) 、(4) 、(6) 、(7) 修正当前负水深网格在 n+1 时间步的流入量 和流出量得出所述负水深网格在 n+2 时间步的实际水深 : 权 利 要 求 书 CN 103559403 A 3 1/8 页 4 洪水演进数值模拟计算中网格流出率的修正方法 技术领域 0001 本发明涉及一种洪。
8、水演进数值模拟计算中网格流出率的修正方法。 背景技术 0002 洪水灾害是指山丘地区在强降雨影响下, 短时间内形成具有较大洪峰流量的洪 水。我国地处东亚季风区, 山区和丘陵地区占国土地面积的三分之二。其中, 洪水灾害的优 先预防面积达 97 万, 影响人口 1.3 亿。近年洪水灾害造成的死亡人数占全国洪涝灾害死亡 人数的比例超过 70%, 成为造成人员伤亡的主要灾种。随着社会经济的发展, 洪水灾害的防 治工作越来越被重视。 0003 早先, 国际通用的洪水灾害的预测技术是对沟道、 沟口进行实地采样, 根据有可能 的灾害种类和等级确定危险指数。最具代表性的是 Aulitzky 提出的荒溪分类及危。
9、险区制 图指数法, 通过收集 9 种指标 51 个具体因子划分出不同等级的危险区。随着地理信息系 统、 数字高程模型、 遥感和卫星遥测等现代科学技术高速发展, 基于平面浅水波方程的模拟 方法被广泛应用于对流域内洪水、 洪水及尤其是山洪等灾害现象的预测和定量分析中。该 方法不受模型实验相似性理论的限制, 可快速、 精确的揭示灾害发生的原因及过程, 从而大 大提高洪水等突发洪水的预见期。最常用的平面二维模拟数值方法包括有限元法、 有限体 积法和有限差分法。针对洪水模拟, Jin 等提出了使用二阶迎风格式离散动量方程的非线 性对流项和二阶 leap-frog 格式离散线性对流项模拟 Nakdong 。
10、流域的泛滥过程。Roger 等 提出一种 MacCormack+TVD 格式的边界拟合数值模型, 通过检测每个时间步长下的水深是 否到达干枯临界值判定动边界范围。丹麦 DHI 水资源与环境研究院采用隐式交替方向算法 开发了水动力学模拟软件 Mike21。 0004 在面向洪水尤其是山洪灾害二维数值模拟中, 由于流域地形陡峭, 水流量急速变 化, 在较大计算时间步长下的动边界处理过程中 , 对流出量简单的归零处理往往导致负水 深, 造成模拟过程中的质量与动量不守恒, 最终导致计算数值不稳定甚至计算发散而得不 到结果。近年来, Satofuka 和 Mizuyama(2005) 采用一维数值模型计。
11、算了明渠洪水对山区 河床、 水坝造成的影响。而 Nakatani 等 (2008) 采用二维数值模拟模型, 基于洪水淹没深度 和沉降的变化对其形成的冲积扇进行了模拟计算。 0005 由于洪水发生地的地形具有陡峭且凹凸不平的特点, 因此洪水的边界范围、 流速 和水深都急剧变化。 在以往的移动边界计算中, 经常发生水面标高低于地面标高的情况, 也 就是计算后的水深变为负值。 导致计算域内的质量不能保证守恒, 稳定性变差, 甚至计算发 散而得不到结果, 从而直接影响到洪水演进数值模拟计算无法顺利进行。 发明内容 0006 本发明所要解决的技术问题是提供一种洪水演进数值模拟计算中网格流出率的 修正方法。
12、, 该修正方法解决了当出现水面标高低于地面标高的情况, 洪水演进数值计算无 法顺利进行的技术问题。 说 明 书 CN 103559403 A 4 2/8 页 5 0007 为了解决上述问题, 本发明提供了一种洪水演进数值模拟计算中网格流出率的修 正方法, 包括 : 0008 通过洪水高发地区的地形图建立一计算域内的该洪水高发地区的数字高程模 型, 并对所述计算域进行二维规则网格划分, 以建立二维规则网格模型及进行网格参数设 定, 建立所述二维规则网格模型及网格参数设定的方法包括 : A : 建立所述二维规则网格模 型, 即所述计算域按照一定空间步长 (X、 Y) 进行分割后得到二维规则网格模型。
13、, 定义 (i, j) 为该二维规则网格模型中一网格, 且该网格为正方形网格 ; B : 网格参数设定, 即设定 该二维规则网格模型中任一网格的水深为 H, 且与该网格按 X 轴同方向的单宽流量 M, 与 Y 轴同方向的单宽流量 N ; 同时设定所述洪水演进数值计算的时间步 n, 以及该时间步 n 对应 的时间步长 t。 0009 其中, 所述数字高程模型, 即 DEM, 用一组有序数值阵列形式表示地面高程的一种 实体地面模型, 定义 (i, j) 为该二维规则网格模型中一网格, 所述网格为正方形网格, 也称 为栅格 DEM ; 所述 i 和 j 的取值 i=1,2,3,, j=1,2,3,,。
14、 利用 i 和 j 来限定网格在所 述二维规则网格模型中的具体位置, 即网格 (1,1) , 网格 (1,2) 类似的表示形式。通过遥感 影像建立一计算域内所述洪水高发地区的数字高程模型的技术方案在现有技术中已经公 开, 这里不再重复。 0010 所述单宽流量 : 单位宽度上河流或输水管的输水流量, 这里的单位宽度即网格。 0011 所述空间步长也可以简称步长, 针对 DEM 空间的分辨力, 也就是网格的精度, 一般 为30M、 90M两种, 就是用30*30或90*90的DEM网格来表示地形, 当然也可以根据计算需要, 另外设置相应步长, 其中, X、 Y 分别表示该二维规则网格模型中一个网。
15、格的长、 宽, 这里 X=Y。 0012 所述时间步 n 为整个计算过程中的时刻间隔, 即计算步, 所述时间步 n 的取值 n=1,2,3,; 所述时间步长 t : 相邻两计算步的时间间隔, 对应的时间步长, 一般可以去 0.1s 或 0.01s, 也可以根据计算设置任意时间。 0013 通过单宽流量公式获得所述二维规则网格模型中各网格在 n+1 时间步的所述 单宽流量 M 和单宽流量 N ; 其方法, 包括 : 0014 预设初始条件, 即, 在 n 时间步, 且设定位于所述计算域的边界上的流量输入起始 界的各边界网格的初值参数, 该初值参数包括 : 初始水深 Hn, 所述单宽流量 M 的对。
16、应流速矢 量Un, 以及所述单宽流量N的对应流速矢量Vn; 在所述单宽流量公式中位于所述计算域外的 网格的初值参数与该网格相邻的一边界网格的初值参数相同, 且位于所述计算域内的网格 的初值参数的相应取值为 0。 0015 以及, 所述单宽流量 M 和单宽流量 N 的初始值的计算公式 : 0016 Mn=UnHn (1) ; 0017 Nn=VnHn (2) ; 0018 所述单宽流量公式 : 0019 说 明 书 CN 103559403 A 5 3/8 页 6 0020 0021 根据所述各网格在n+1时间步的单宽流量M和单宽流量N, 建立所述各网格的水 深计算式, 即 0022 0023 。
17、其 中, n+2 时 间 步 作 为 水 深 计 算 的 起 始 时 间 步, H 表 示 所 述 网 格 水 深, 表示在所述起始时间步时网格 (i, j) 的区域内的水深 ;表 示在n时间步的所述网格的区域内的水深, 且该时间步的水深为预设值 ;表示在n+1 时间步的所述网格在 X 轴方向上流体的流入单宽流量 ;表示相邻网格在 X 轴方向 上的流入单宽流量, 即, 所述网格在 X 轴方向上的流出单宽流量 ;表示在 n+1 时间步 的所述网格在 Y 轴方向上流体的流入单宽流量 ;表示相邻网格在 Y 轴方向上的流入 单宽流量, 即, 所述网格在 Y 轴方向上的流出单宽流量。 0024 若步骤计。
18、算出一网格的水深 H 小于 0, 即该网格为负水深网格, 则引入一修正 系数 重新对该负水深网格的水量进行修正计算, 直到该负水深网格修正为实际水深后, 再转入步骤通过所述单宽流量公式获得所述二维规则网格模型中各网格在下一时间步 的单宽流量 M 和单宽流量 N。 0025 其中, 当一网格出现负水深时, 根据负水深网格在当前时间步n的单宽流量M和单 宽流量 N 计算出该负水深网格的流入总量和流出总量 0026 所述流入总量公式 : 0027 0028 所述流出总量公式 : 0029 0030 所述计算修正系数 的公式 : 0031 说 明 书 CN 103559403 A 6 4/8 页 7 。
19、0032 并根据公式 (1) 、(2) 、(3) 、(4) 、(6) 、(7) 修正当前负水深网格在 n+1 时间步的流 入量和流出量得出所述负水深网格在 n+2 时间步的实际水 深 : 0033 0034 本发明相对于现有技术具有积极的效果 :(1) 本发明通过洪水高发地区的地形图 建立一计算域内的该洪水高发地区的数字高程模型, 即利用地形图建立数字高程模型, 由 于地形图很容易获得, 费用较低, 地形图上的等高线含有丰富的高程信息, 利用这些高程信 息生产 DEM 数据, 可以极大的发挥地形图的作用, 通过地形图建立数字高程模型的方式与 利用遥感影像建立数字高程的方式相比, 具有成本低, 。
20、而且对硬件要求低等特点 ;(2) 本发 明克服了在以往的移动边界计算中, 经常发生水面标高低于地面标高的情况, 也就是计算 后的水深变为负值, 即出现负水深网格, 从而导致计算域内的质量不能保证守恒, 稳定性变 差, 甚至计算发散而得不到结果, 直接影响到洪水演进数值计算无法顺利进行的技术问题, 通过引入所述修正系数 , 使当出现水面标高低于地面标高 (负水深网格) 时, 洪水演进数 值计算得已顺利进行。 附图说明 0035 为了清楚说明本发明的创新原理及其相比于现有产品的技术优势, 下面借助于附 图通过应用所述原理的非限制性实例说明一个可能的实施例。在图中 : 0036 图 1 为本发明的一。
21、种洪水演进数值模拟计算中网格流出率的修正方法的流程图 ; 0037 图 2 为本发明中公式 (6) 和 (7) 的水流方向示意图。 具体实施方式 0038 参照附图, 对本发明的实施方式和实施例进行详细说明。 0039 图 1 示出本发明的洪水演进数值模拟计算中网格流出率的修正方法的流程图, 现 结合图 1 所述流程图对该修正方法的各步骤作具体说明。 0040 在步骤 S001 中, 通过洪水高发地区的地形图建立一计算域内的该洪水高发地区 的数字高程模型, 并对所述计算域进行二维规则网格划分, 以建立二维规则网格模型及进 行网格参数设定。 0041 其中, 通过所述地形图建立相应地区的数字高程。
22、模型的方法如下 : 0042 DEM 数据的生成包括如下几个过程 : 说 明 书 CN 103559403 A 7 5/8 页 8 0043 (1) 数字化栅格底图的生成 0044 该过程包括如下几个步骤 : 0045 a、 地形图的扫描 : 即对购买来的地形图通过扫描仪扫描, 以获得数字化栅格底 图 ; 0046 b、 对栅格底图进行纠正, 以消除纸张变形所带来的误差。 0047 (2) 高程信息的获取 0048 高程信息的获取, 就是在栅格底图上进行数字化, 以获取生成 DEM 所必须的高程 信息, 它包括如下几个步骤 : 0049 a、 对等高线进行矢量化 ; 0050 b、 对等高线标。
23、赋高程值 ; 0051 c、 对离散高程点进行矢量化 ; 0052 d、 对高程点标赋高程 ; 0053 e、 对这些矢量化产品进行编辑、 检查、 拼接以生成拓扑关系完整的矢量图。 0054 (3) DEM 生成 0055 a、 将生成的矢量图在 ARC/NFO 软件中利用不规则三角网 (TN) 进行内插, 以使整个 区域都含有高程值 ; 0056 b、 将 TN 数据进行采样, 转换为 GRD 数据 (lat-tice 格式) , 对数据进行进一步检 查 ; 0057 c、 将生成好的 GRD 数据转成 DEM 数据 (USGS 格式) , 即获得所述洪水高发地区的数 字高程模型。 0058。
24、 其中, 建立所述二维规则网格模型及参数设定的方法包括 : 0059 A : 建立所述二维规则网格模型, 即所述计算域按照一定空间步长 (X、 Y) 进行 分割后得到二维规则网格模型, 定义 (i, j) 为该二维规则网格模型中一网格, 且该网格为正 方形网格, 即 X=Y。 0060 B : 网格参数设定, 即设定该二维规则网格模型中任一网格的水深为 H, 且与该网 格按 X 轴同方向的单宽流量 M, 与 Y 轴同方向的单宽流量 N ; 同时设定所述洪水演进数值计 算的时间步 n, 以及该时间步 n 对应的时间步长 t ; 其中, 单宽流量 M 和单宽流量 N 均为矢 量, 也可以用于表示相。
25、应流量的方向。 0061 其中, 所述数字高程模型, 即 DEM, 用一组有序数值阵列形式表示地面高程的一种 实体地面模型, 定义 (i, j) 为该二维规则网格模型中一网格, 所述网格为正方形网格, 也称 为栅格 DEM ; 所述 i 和 j 的取值 i=1,2,3,, j=1,2,3,, 利用 i 和 j 来限定网格在所 述二维规则网格模型中的具体位置, 即网格 (1,1) , 网格 (1,2) 类似的表示形式。通过遥感 影像建立一计算域内所述洪水高发地区的数字高程模型的技术方案在现有技术中已经公 开, 这里不再重复。 0062 所述单宽流量 : 单位宽度上河流或输水管的输水流量, 这里的。
26、单位宽度即网格。 0063 所述空间步长也可以简称步长, 针对 DEM 空间的分辨力, 也就是网格的精度, 一般 为30M、 90M两种, 就是用30*30或90*90的DEM网格来表示地形, 当然也可以根据计算需要, 另外设置相应步长, 其中, X、 Y 分别表示该二维规则网格模型中一个网格的长、 宽, 这里 X=Y ; 也可以采用不规则网格, 即长方形网格 (矩形网格) 的方式实现本实施例, 过程与 说 明 书 CN 103559403 A 8 6/8 页 9 正方形相同, 这里不再重复。 0064 所述时间步 n 为整个计算过程中的时刻间隔, 即计算步, 所述时间步 n 的取值 n=1,。
27、2,3,; 所述时间步长 t : 相邻两计算步的时间间隔, 对应的时间步长, 一般可以去 0.1s 或 0.01s, 也可以根据计算设置任意时间。 0065 在步骤S002中, 通过单宽流量公式获得所述二维规则网格模型中各网格在n+1时 间步的所述单宽流量 M 和单宽流量 N ; 其方法, 包括 : 0066 预设初始条件, 即, 在 n 时间步, 且设定位于所述计算域的边界上的流量输入起始 界的各边界网格的初值参数, 该初值参数包括 : 初始水深 Hn, 所述单宽流量 M 的对应流速矢 量Un, 以及所述单宽流量N的对应流速矢量Vn; 在所述单宽流量公式中位于所述计算域外的 网格的初值参数与。
28、该网格相邻的一边界网格的初值参数相同, 且位于所述计算域内的网格 的初值参数的相应取值为 0。 0067 以及, 所述单宽流量 M 和单宽流量 N 的初始值的计算公式 : 0068 Mn=UnHn (1) ; 0069 Nn=VnHn (2) ; 0070 所述单宽流量公式 : 0071 0072 0073 其中, 所述初值参数可以通过多参数水文监测仪来获得, 获得方式属于现有技术, 这里不再详细叙述。 0074 在步骤S003中, 根据所述各网格在n+1时间步的单宽流量M和单宽流量N, 建立所 述各网格的水深计算式, 即 0075 0076 其 中, n+2 时 间 步 作 为 水 深 计 。
29、算 的 起 始 时 间 步, H 表 示 所 述 网 格 水 深, 表示在所述起始时间步时网格 (i, j) 的区域内的水深 ; 表示在 n 时间步的所述网格的区域内的水深, 且该时间步的水深为预设值 ;表示在 n+1 时间步的所述网格在 X 轴方向上流体的流入单宽流量 ;表示相邻网格在 X 轴方 说 明 书 CN 103559403 A 9 7/8 页 10 向上的流入单宽流量, 即, 所述网格在 X 轴方向上的流出单宽流量 ;表示在 n+1 时间 步的所述网格在Y轴方向上流体的流入单宽流量 ;表示相邻网格在Y轴方向上的流 入单宽流量, 即, 所述网格在 Y 轴方向上的流出单宽流量。 007。
30、7 在步骤 S004 中, 判断 n+2 时间步是否有负水深网格存在。 0078 所述网格流出率修正具体步骤包括 : 若步骤 S003 中计算出一网格的水深 H 小于 0, 即该网格为负水深网格, 则引入一修正系数 重新对该负水深网格的水量进行修正计 算, 直到该负水深网格修正为实际水深后, 再转入步骤 S002 通过所述单宽流量公式获得所 述二维网格模型中各网格在下一时间步的单宽流量 M 和单宽流量 N, 该下一时间步表示一 种循环关系, 即开始下一个计算周期。 0079 在步骤 S005 中, 结合图 2, 其中, 当一网格出现负水深时, 根据负水深网格在当前 时间步 n 的单宽流量 M 。
31、和单宽流量 N 计算出该负水深网格的流入总量和流出总 量 0080 所述流入总量公式 : 0081 所述流出总量公式 : 0082 在步骤 S006 中, 用于所述计算修正系数 , 该计算修正系数 的公式 : 0083 在步骤 S007 中, 根据公式 (1) 、(2) 、(3) 、(4) 、(6) 、(7) 修正当前负水深网格在 n+1 时间步的流入量和流出量即, 按照公式 (1) 、(2) 、(3) 、(4) 获 得的在 n+1 时间步的单宽流量 M 和单宽流量 N 分别代入公式 (6) 、(7) 。 0084 在步骤 S008 中, 根据上述步骤 S007 和 S006 中的相应计算, 。
32、得出所述负水深网格 在 n+2 时间步的实际水深 : 0085 说 明 书 CN 103559403 A 10 8/8 页 11 0086 在步骤 S009 中, 输出网格的实际水深。 0087 在步骤 S010 中, 转入步骤 S002 通过所述单宽流量公式获得所述二维网格模型中 各网格在下一时间步的单宽流量 M 和单宽流量 N。 0088 显然, 上述实施例仅仅是为清楚地说明本发明所作的举例, 而并非是对本发明的 实施方式的限定。对于所属领域的普通技术人员来说, 在上述说明的基础上还可以做出其 它不同形式的变化或变动。这里无需也无法对所有的实施方式予以穷举。而这些属于本发 明的精神所引伸出的显而易见的变化或变动仍处于本发明的保护范围之中。 说 明 书 CN 103559403 A 11 1/2 页 12 图 1 说 明 书 附 图 CN 103559403 A 12 2/2 页 13 图 2 说 明 书 附 图 CN 103559403 A 13 。