一种地表径流模拟系统及方法

技术领域

本发明涉及仿真模拟技术领域，尤其涉及一种地表径流模拟系统及方法。

背景技术

洪涝灾害是发生频率高、危害范围广且对国民经济影响最为严重的自然灾害，洪灾具有明显的季节性、区域性和可重复性。根据相关资料可知，降雨是直接导致洪灾发生的一个重要因素，在以往对洪灾进行径流生态调研实践中可知，径流生态调控可以通过植被、土壤等介质对大气降水的截留与吸收，强化雨水就地利用与异地叠加利用的能力，在减少水土流失的同时，能够最大程度地增加土壤蓄水量，改善区内水文环境。

地表径流是指由降水或冰雪融化形成的、沿着流域的不同路径流入河流、湖泊或海洋的水流。模拟径流的成流过程能够为径流生态调控提供数据支撑，如何实现地表径流的动态模拟，以直观查看地表径流随时间变化成为亟待解决的技术问题。

发明内容

针对现有技术中存在的问题，本发明提供一种地表径流模拟系统，包括：

数据预处理模块，用于获取待模拟区域的数字高程模型、降水网格数据和土地网格数据，并对所述数字高程模型进行网格化处理得到包含多个区域子网格的区域网格，随后根据所述降水网格数据和所述土地网格数据处理得到各所述区域子网格关联的降水数据和土地类型数据；

径流模拟模块，连接所述数据预处理模块，用于构建水动力控制方程，将各所述区域子网格关联的所述降水数据和所述土地类型数据作为所述区域网格的当前时刻的初始状态数据，对所述水动力控制方程进行迭代求解得到后续各时刻各所述区域子网格的地表径流的向上位移量，并对关联有所述向上位移量的各所述区域子网络进行渲染着色得到所述待模拟区域的地表径流动态模拟结果。

优选的，所述降水网格数据包含多个降水子网格，每个所述降水子网格关联有相应的第一经纬度坐标和单位时间降雨量；所述土地网格数据包含多个土地子网格，每个所述土地子网格关联有相应的第二经纬度坐标和随时间变化的渗水变化曲线；则所述数据预处理模块包括：

网格化处理单元，用于对获取的所述待模拟区域的所述数字高程模型进行网格化处理得到包含多个所述区域子网格的所述区域网格；

各所述区域子网格关联有相应的第三经纬度坐标；

第一关联单元，连接所述网格化处理单元，用于针对每个所述区域子网格，根据各个所述降水子网格的所述第一经纬度坐标和所述区域子网格关联的所述第三经纬度坐标处理得到覆盖所述区域子网格的各所述降水子网格的覆盖面积占所述区域子网格表面积的覆盖比例，随后将对应的各所述覆盖比例作为权重，将对应的各所述降水子网格关联的单位时间降雨量进行加权求和得到所述区域子网络关联的所述降水数据；

第二关联单元，连接所述网格化处理单元，用于针对每个所述区域子网格，根据各个所述土地子网格的所述第二经纬度坐标和所述区域子网格关联的所述第三经纬度坐标处理得到覆盖所述区域子网格的各所述土地子网格的覆盖面积占所述区域子网格表面积的覆盖比例，随后将对应的各所述覆盖比例作为权重，将对应的各所述土地子网格关联的所述渗水变化曲线中对应时刻的渗水率进行加权求和得到所述区域子网络关联的所述土地类型数据。

优选的，所述网格化处理单元包括：

网格化子单元，用于对获取的所述待模拟区域的所述数字高程模型进行网格化处理得到包含多个第一子网格的数字高程网格；

处理子单元，连接所述网格化子单元，用于针对每个所述第一子网格，将所述第一子网格以及包围所述第一子网格的多个其他所述第一子网格作为一个所述区域子网格，得到所述区域网格。

优选的，所述径流模拟模块包括：

第一求解单元，用于将每个所述区域子网格作为一个微元控制体，根据所述区域网格的当前时刻的所述初始状态数据，采用有限体积法并基于水质量瞬时守恒对所述水动力控制方程进行求解，得到当前时刻的各所述区域子网格的地表水质量；

第二求解单元，连接所述第一求解单元，用于将当前时刻的各所述地表水质量作为已知量代入所述水动力控制方程，迭代求解下一时刻各所述区域子网格的地表水质量；

径流处理单元，分别连接所述第一求解单元和所述第二求解单元，用于针对求解得到的每个时刻的各所述区域子网格的所述地表水质量，根据地表水密度和每个所述区域子网格的表面积处理得到各所述区域子网格的地表径流的所述向上位移量，并建立各所述向上位移量与对应的所述区域子网格的关联关系；

径流模拟单元，连接所述径流处理单元，用于将关联有所述向上位移量的各所述区域子网络进行渲染着色得到所述待模拟区域的所述地表径流动态模拟结果。

优选的，所述水动力控制方程为欧拉方程。

本发明还提供一种地表径流模拟方法，应用于上述的地表径流模拟系统，所述地表径流模拟方法包括：

步骤S1，所述地表径流模拟系统获取待模拟区域的数字高程模型、降水网格数据和土地网格数据，并对所述数字高程模型进行网格化处理得到包含多个区域子网格的区域网格，随后根据所述降水网格数据和所述土地网格数据处理得到各所述区域子网格关联的降水数据和土地类型数据；

步骤S2，所述地表径流模拟系统构建水动力控制方程，将各所述区域子网格关联的所述降水数据和所述土地类型数据作为所述区域网格的当前时刻的初始状态数据，对所述水动力控制方程进行迭代求解得到后续各时刻各所述区域子网格的地表径流的向上位移量，并对关联有所述向上位移量的各所述区域子网络进行渲染着色得到所述待模拟区域的地表径流动态模拟结果。

优选的，所述降水网格数据包含多个降水子网格，每个所述降水子网格关联有相应的第一经纬度坐标和单位时间降雨量；所述土地网格数据包含多个土地子网格，每个所述土地子网格关联有相应的第二经纬度坐标和随时间变化的渗水变化曲线；则所述步骤S1包括：

步骤S11，所述地表径流模拟系统对获取的所述待模拟区域的所述数字高程模型进行网格化处理得到包含多个所述区域子网格的所述区域网格；

各所述区域子网格关联有相应的第三经纬度坐标；

步骤S12，所述地表径流模拟系统针对每个所述区域子网格，根据各个所述降水子网格的所述第一经纬度坐标和所述区域子网格关联的所述第三经纬度坐标处理得到覆盖所述区域子网格的各所述降水子网格的覆盖面积占所述区域子网格表面积的覆盖比例，随后将对应的各所述覆盖比例作为权重，将对应的各所述降水子网格关联的单位时间降雨量进行加权求和得到所述区域子网络关联的所述降水数据；

步骤S13，所述地表径流模拟系统针对每个所述区域子网格，根据各个所述土地子网格的所述第二经纬度坐标和所述区域子网格关联的所述第三经纬度坐标处理得到覆盖所述区域子网格的各所述土地子网格的覆盖面积占所述区域子网格表面积的覆盖比例，随后将对应的各所述覆盖比例作为权重，将对应的各所述土地子网格关联的所述渗水变化曲线中对应时刻的渗水率进行加权求和得到所述区域子网络关联的所述土地类型数据。

优选的，所述步骤S11包括：

步骤S111，所述地表径流模拟系统对获取的所述待模拟区域的所述数字高程模型进行网格化处理得到包含多个第一子网格的数字高程网格；

步骤S112，所述地表径流模拟系统针对每个所述第一子网格，将所述第一子网格以及包围所述第一子网格的多个其他所述第一子网格作为一个所述区域子网格，得到所述区域网格。

优选的，所述步骤S2包括：

步骤S21，所述地表径流模拟系统将每个所述区域子网格作为一个微元控制体，根据所述区域网格的当前时刻的所述初始状态数据，采用有限体积法并基于水质量瞬时守恒对所述水动力控制方程进行求解，得到当前时刻的各所述区域子网格的地表水质量；

步骤S22，所述地表径流模拟系统将当前时刻的各所述地表水质量作为已知量代入所述水动力控制方程，迭代求解下一时刻各所述区域子网格的地表水质量；

步骤S23，所述地表径流模拟系统针对求解得到的每个时刻的各所述区域子网格的所述地表水质量，根据地表水密度和每个所述区域子网格的表面积处理得到各所述区域子网格的地表径流的所述向上位移量，并建立各所述向上位移量与对应的所述区域子网格的关联关系；

步骤S24，所述地表径流模拟系统将关联有所述向上位移量的各所述区域子网络进行渲染着色得到所述待模拟区域的所述地表径流动态模拟结果。

优选的，所述水动力控制方程为欧拉方程。

上述技术方案具有如下优点或有益效果：

1)结合数字高程模型并综合考虑影响地表径流的降水数据和土地类型数据，构建水动力控制方程进行迭代求解，实现能够基于土地类型数据随时间和降雨数据变化的规律进行地表径流模拟，并能够实现地表径流变化的动态可视化展示，为径流生态调控提供有力数据支撑，同时便于专业人员进行径流灾害评估；

2)水动力控制方程的构建和求解分步考虑了瞬时的质量守恒、能量守恒以及动量守恒，以及全局的能量守恒及动量守恒，且由于待模拟区域足够大而无需考虑粘滞问题，有效提升数据运算效率。

附图说明

图1为本发明的较佳的实施例中，一种地表径流模拟系统的结构示意图；

图2为本发明的较佳的实施例中，一种地表径流模拟方法的流程示意图；

图3为本发明的较佳的实施例中，区域子网格和降水子网格的示意；

图4为本发明的较佳的实施例中，步骤S1的子流程示意图；

图5为本发明的较佳的实施例中，步骤S11的子流程示意图；

图6为本发明的较佳的实施例中，步骤S2的子流程示意图。

具体实施方式

下面结合附图和具体实施例对本发明进行详细说明。本发明并不限定于该实施方式，只要符合本发明的主旨，则其他实施方式也可以属于本发明的范畴。

本发明的较佳的实施例中，基于现有技术中存在的上述问题，现提供一种地表径流模拟系统，如图1所示，包括：

数据预处理模块1，用于获取待模拟区域的数字高程模型、降水网格数据和土地网格数据，并对数字高程模型进行网格化处理得到包含多个区域子网格的区域网格，随后根据降水网格数据和土地网格数据处理得到各区域子网格关联的降水数据和土地类型数据；

径流模拟模块2，连接数据预处理模块1，用于构建水动力控制方程，将各区域子网格关联的降水数据和土地类型数据作为区域网格的当前时刻的初始状态数据，对水动力控制方程进行迭代求解得到后续各时刻各区域子网格的地表径流的向上位移量，并对关联有向上位移量的各区域子网络进行渲染着色得到待模拟区域的地表径流动态模拟结果。

具体地，本技术方案的地表径流模拟系统优选部署于云端，一方面能够降低对本地部署硬件的成本及算力要求，另一方面提升可扩展性。进一步地，上述数字高程模型、降水网格数据和土地网格数据可以通过云端对接多个第三方业务平台进行获取，但不并以此进行限定。

本实施例中，以网格作为处理对象，网格是一种有限元的仿真技术，由于计算机的计算是离散的，用一定区域的平均值代替该区域的分布，这种分布被称之为网格。具体地，选取与地表径流形成相关性较强的降水网格数据、土地网格数据以及数字高程模型包含的地形数据进行地表径流模拟。具体地，在对数字高程模型进行网格化处理后，可以理解的是，获取的降水网格数据和土地网格数据对应的网格尺寸与数字高程模型网格化处理后得到的区域网格的区域子网格的尺寸通常是不一致的，因此在进行地表径流模拟之前，首先需要将降水网格数据包含的降水数据和土地网格数据包含的土地类型数据映射到各个区域子网格中，具体包括：

降水网格数据包含多个降水子网格，每个降水子网格关联有相应的第一经纬度坐标和单位时间降雨量；土地网格数据包含多个土地子网格，每个土地子网格关联有相应的第二经纬度坐标和随时间变化的渗水变化曲线；则数据预处理模块1包括：

网格化处理单元11，用于对获取的待模拟区域的数字高程模型进行网格化处理得到包含多个区域子网格的区域网格；

各区域子网格关联有相应的第三经纬度坐标；

第一关联单元12，连接网格化处理单元11，用于针对每个区域子网格，根据各个降水子网格的第一经纬度坐标和区域子网格关联的第三经纬度坐标处理得到覆盖区域子网格的各降水子网格的覆盖面积占区域子网格表面积的覆盖比例，随后将对应的各覆盖比例作为权重，将对应的各降水子网格关联的单位时间降雨量进行加权求和得到区域子网络关联的降水数据；

第二关联单元13，连接网格化处理单元11，用于针对每个区域子网格，根据各个土地子网格的第二经纬度坐标和区域子网格关联的第三经纬度坐标处理得到覆盖区域子网格的各土地子网格的覆盖面积占区域子网格表面积的覆盖比例，随后将对应的各覆盖比例作为权重，将对应的各土地子网格关联的渗水变化曲线中对应时刻的渗水率进行加权求和得到区域子网络关联的土地类型数据。

具体地，本实施例中，各区域子网格关联的第三经纬度坐标为每个区域子网格的四个顶点的经纬度坐标，同时，每个顶点还关联有高程数据，为统一量纲，优选将高程数据转换为海拔。

上述降水网格数据对应的原始数据通常为泰勒多边形形式，优选基于泰勒多边形形式的原始数据生成单位时间降雨量的热力图，进而对热力图进行网格化即可得到降水网格数据。

进一步具体地，如图2所示，其中，各虚线框表示各区域子网格，四个实线框表示降水子网格，可以看出，对于区域子网格a11、区域子网格a12、区域子网格a21和区域子网格a22，其完全被降水子网格b1覆盖，即覆盖比例为100％，降水子网格b1的单位时间降水量即为区域子网格a11、区域子网格a12、区域子网格a21和区域子网格a22的降水数据。对于区域子网格a13、区域子网格a23、区域子网格a24、区域子网格a25和区域子网格a26，其同时被多个降水子网格覆盖，需要进行加权求和。以区域子网格a26为例，其同时被降水子网格b1、降水子网格b2、降水子网格b3和降水子网格b4覆盖，其中可以根据区域子网格a26四个顶点的第三经纬度坐标、降水子网格b1的四个顶点的第一经纬度坐标处理得到降水子网格b1的覆盖面积占区域子网格a26表面积的第一占比r1，以此类推，可以分别处理得到降水子网格b2的覆盖面积占区域子网格a26表面积的第二占比r2，降水子网格b3的覆盖面积占区域子网格a26表面积的第三占比r3，降水子网格b4的覆盖面积占区域子网格a26表面积的第四占比r4，即可基于以下公式计算得到区域子网格a26的降水数据P：

P＝P1*r1+P2*r2+P3*r3+P4*r4

其中，P1为降水子网格b1的单位时间降雨量，P2为降水子网格b2的单位时间降雨量，P3为降水子网格b3的单位时间降雨量，P4为降水子网格b4的单位时间降雨量。

其他降水子网格的降水数据的处理方式以此类推，此处不再赘述。

同样地，区域子网络关联的土地类型数据的处理方式与降水数据相同，区别仅在于将参与计算的单位时间降雨量替换为各土地子网格关联的渗水变化曲线中对应时刻的渗水率。需要说明的是，土地网格数据用于表示对应区域的土地类型，如沙漠、草原、耕地等，可以理解的是，土地类型不同对应的渗水能力不同，本实施例中采用渗水率随时间变化的渗水变化曲线表征土地类型。其中渗水变化曲线的获取方式可以是现场实验采集，或者实验室渗水实验采集。

进一步地，土地类型用于量化当前土地含水量关于雨量的饱和曲线，其一阶导数为该土地的渗透函数(雨量大小和渗透率变化之间的函数关系，渗透率是单位面积单位时间内渗透水量占全部水量的百分比。雨量大小则是雨量函数在时间上的积分)，其中：

F(留存地表径流)＝G(1-(当前土地渗透率)*雨量)

H(当前土地渗透率)＝V(土地饱和度-雨量)的一阶导数。

本发明的较佳的实施例中，网格化处理单元11包括：

网格化子单元111，用于对获取的待模拟区域的数字高程模型进行网格化处理得到包含多个第一子网格的数字高程网格；

处理子单元112，连接网格化子单元111，用于针对每个第一子网格，将第一子网格以及包围第一子网格的多个其他第一子网格作为一个区域子网格，得到区域网格。

具体地，考虑到降水网格数据和土地网格数据对应的网格尺寸通常较大，而上述第一子网格的网格尺寸相对较小，而在进行地表径流模拟时通常对精细度的要求不高，基于此，本实施例中，针对每个子网格，将包围第一子网格的多个其他第一子网格作为一个区域子网格参与计算，一方面降低运算量，另一方面能够提升数据平滑型，消除锯齿。优选的，可以以九宫格的形式选取包围一个区域子网络的多个其他第一子网格，即一个区域子网格可以包含9个第一子网格。

在处理得到各区域子网格关联的降水数据和土地类型数据之后，即可以该降水数据和土地类型数据表征的初始状态开始进行径流模拟，具体地，径流模拟模块2包括：

第一求解单元21，用于将每个区域子网格作为一个微元控制体，根据区域网格的当前时刻的初始状态数据，采用有限体积法并基于水质量瞬时守恒对水动力控制方程进行求解，得到当前时刻的各区域子网格的地表水质量；

第二求解单元22，连接第一求解单元21，用于将当前时刻的各地表水质量作为已知量代入水动力控制方程，迭代求解下一时刻各区域子网格的地表水质量；

径流处理单元23，分别连接第一求解单元21和第二求解单元22，用于针对求解得到的每个时刻的各区域子网格的地表水质量，根据地表水密度和每个区域子网格的表面积处理得到各区域子网格的地表径流的向上位移量，并建立各向上位移量与对应的区域子网格的关联关系；

径流模拟单元24，连接径流处理单元23，用于将关联有向上位移量的各区域子网络进行渲染着色得到待模拟区域的地表径流动态模拟结果。

具体地，本实施例中，考虑到在进行径流模拟时选取的待模拟区域的面积足够大，可以忽略流体粘滞的影响，即水动力控制方程为欧拉方程，具体表达式如下：

这是一个矢量表达式，其中，ρ表示水密度，v表示流速，Q表示重力势能分量。

上述

由于土地存在渗水现象，因此在全局状态下，水质量不守恒，仅能量守恒及动量守恒，则上述水动力控制方程存在两个未知数，分别为水质量和时间。虽然全局状态下水质量不守恒，但在Δt足够小时，即瞬时状态下，可以认为水质量是守恒的，换言之，不用考虑瞬态项，则上述水动力控制方程转换为质量守恒式，表达式如下：

此时，仅具有水质量一个未知量需要求解，基于此，在瞬时的质量守恒、能量守恒以及动量守恒下，可以将每个区域子网格作为一个微元控制体，并认为与该微元控制体存在对流的分别为位于该区域子网格上下左右的四个区域子网格，进而采用有限体积法对上述质量守恒式进行求解即可得到当前时刻的各区域子网格的地表水质量。随后考虑瞬态项，此时认为地表水质量为常数，即可迭代求解下一时刻各区域子网格的地表水质量。在下一刻地表水质量求解完成后，可以将下一刻地表水质量作为初始状态数据，再求解再下一时刻的地表水质量，依次迭代即可。

进一步地，在求解得到上述地表水质量后，水密度已知，即可得到每个区域子网格的地表水体积，进而将地表水体积除以该区域子网格的表面积即可得到该区域子网格的地表径流的向上位移量，随后建立该向上位移量与区域子网格之间的关联关系，从而渲染着色即可得到可视化的地表径流动态模拟结果。上述渲染着色包含但不限于采用UE5引擎实现。

更进一步地，为实现地表径流动态模拟结果的平滑性，在进行渲染着色之前，还需要考虑各个区域子网格的之间的交点处的地表径流的向上位移量，优选为该交点连接的各个区域子网格的地表径流的向上位移量的平均值。以四个区域子网格连接的中心点为例，其对应的向上位移量可以是四个区域子网格的地表径流的向上位移量的平均值。

本发明还提供一种地表径流模拟方法，应用于上述的地表径流模拟系统，如图3所示，地表径流模拟方法包括：

步骤S1，地表径流模拟系统获取待模拟区域的数字高程模型、降水网格数据和土地网格数据，并对数字高程模型进行网格化处理得到包含多个区域子网格的区域网格，随后根据降水网格数据和土地网格数据处理得到各区域子网格关联的降水数据和土地类型数据；

步骤S2，地表径流模拟系统构建水动力控制方程，将各区域子网格关联的降水数据和土地类型数据作为区域网格的当前时刻的初始状态数据，对水动力控制方程进行迭代求解得到后续各时刻各区域子网格的地表径流的向上位移量，并对关联有向上位移量的各区域子网络进行渲染着色得到待模拟区域的地表径流动态模拟结果。

本发明的较佳的实施例中，降水网格数据包含多个降水子网格，每个降水子网格关联有相应的第一经纬度坐标和单位时间降雨量；土地网格数据包含多个土地子网格，每个土地子网格关联有相应的第二经纬度坐标和随时间变化的渗水变化曲线；如图4所示，则步骤S1包括：

步骤S11，地表径流模拟系统对获取的待模拟区域的数字高程模型进行网格化处理得到包含多个区域子网格的区域网格；

各区域子网格关联有相应的第三经纬度坐标；

步骤S12，地表径流模拟系统针对每个区域子网格，根据各个降水子网格的第一经纬度坐标和区域子网格关联的第三经纬度坐标处理得到覆盖区域子网格的各降水子网格的覆盖面积占区域子网格表面积的覆盖比例，随后将对应的各覆盖比例作为权重，将对应的各降水子网格关联的单位时间降雨量进行加权求和得到区域子网络关联的降水数据；

步骤S13，地表径流模拟系统针对每个区域子网格，根据各个土地子网格的第二经纬度坐标和区域子网格关联的第三经纬度坐标处理得到覆盖区域子网格的各土地子网格的覆盖面积占区域子网格表面积的覆盖比例，随后将对应的各覆盖比例作为权重，将对应的各土地子网格关联的渗水变化曲线中对应时刻的渗水率进行加权求和得到区域子网络关联的土地类型数据。

本发明的较佳的实施例中，如图5所示，步骤S11包括：

步骤S111，地表径流模拟系统对获取的待模拟区域的数字高程模型进行网格化处理得到包含多个第一子网格的数字高程网格；

步骤S112，地表径流模拟系统针对每个第一子网格，将第一子网格以及包围第一子网格的多个其他第一子网格作为一个区域子网格，得到区域网格。

本发明的较佳的实施例中，如图6所示，步骤S2包括：

步骤S21，地表径流模拟系统将每个区域子网格作为一个微元控制体，根据区域网格的当前时刻的初始状态数据，采用有限体积法并基于水质量瞬时守恒对水动力控制方程进行求解，得到当前时刻的各区域子网格的地表水质量；

步骤S22，地表径流模拟系统将当前时刻的各地表水质量作为已知量代入水动力控制方程，迭代求解下一时刻各区域子网格的地表水质量；

步骤S23，地表径流模拟系统针对求解得到的每个时刻的各区域子网格的地表水质量，根据地表水密度和每个区域子网格的表面积处理得到各区域子网格的地表径流的向上位移量，并建立各向上位移量与对应的区域子网格的关联关系；

步骤S24，地表径流模拟系统将关联有向上位移量的各区域子网络进行渲染着色得到待模拟区域的地表径流动态模拟结果。

本发明的较佳的实施例中，水动力控制方程为欧拉方程。

以上所述仅为本发明较佳的实施例，并非因此限制本发明的实施方式及保护范围，对于本领域技术人员而言，应当能够意识到凡运用本说明书及图示内容所作出的等同替换和显而易见的变化所得到的方案，均应当包含在本发明的保护范围内。