一种适用于混合产流区的陆面水文模式的建立方法

摘要

本发明提供一种适用于混合产流区的陆面水文模式的建立方法，包括采用双预热方式开展数值大气模式WRF模拟和/或预报,获得气象驱动场；构建WRF‑Hydro‑DQH陆面水文模式产流模块；构建WRF‑Hydro‑DQH陆面水文模式汇流模块；编译运行WRF‑Hydro‑DQH陆面水文模式；调整校正所述WRF‑Hydro‑DQH陆面水文模式的参数；进行所述数值大气模式WRF与校正后的所述WRF‑Hydro‑DQH陆面水文模式在混合产流区的水文变量耦合模拟和/或预报。本发明提出一种适用于混合产流区的陆面水文模式的建立方法，基于近两年国际关注度高的WRF‑Hydro模式无法准确模拟/预报具有区域特性的水文过程现状，针对混合产流区的产汇流机理，对模式系统进行了改进，提高了陆面水文模式在我国混合产流区上的模拟精度。

说明书

技术领域

本发明涉及水文气象的技术领域，特别是一种适用于混合产流区的陆面水文模式的建立方法。

背景技术

WRF-Hydro作为近年来逐渐被国际关注的陆面水文模式，随着模式系统研发，其可服务应用于包括陆气耦合研究、山洪灾害预测、区域水文气候影响评估、水资源季节性预报、水资源调度等多个水文气象环节。虽然如此，其在降雨产流机制描述上，在研发之处由于需要考虑到大尺度的模拟要求，缺少针对区域或局部特性上的技术支撑，因此在模拟时，相应的模拟效果因地域而差异较大。尤其在混合产流区模拟时，误差会被放大。针对超渗产流区暴雨洪水突发性强，汇流时间快等特点，并没有一种适用于混合产流区的陆面水文模式。

当前水文模型对超渗与蓄满混合产流机制的概化与实际仍有较大差异。迄今为止，识别复杂降雨产流过程下超渗产流和蓄满产流机制，并构建混合产流模型仍是水文模拟研究的一个薄弱环节。无论是概念性，还是半分布式再到全分布式的各种模型，受处理空间范围、模拟时间尺度和不同的模拟用途产生了不同的混合产流模型结构。比如：萨克拉门托模型（SAC）用变动透水面积以及饱和下渗率反映超渗和蓄满产流机制；包为民等提出的垂向混合径流模型（VMM）将超渗与蓄满过程在垂向上进行了结合；改进的广义流域负荷函数模型（GWLF）通过用函数方程取代固定的“传导系数”来更好地估计干旱地区和潮湿地区的水文状况；可变下渗与新安江组合模型（VXM）将VIC下渗方程融合到新安江模型蓄满产流结构；河北雨洪模型（HBM）采用先超后蓄的混合产流结构，在下渗中对Horton经验方程进行改进并考虑了表层土壤含水量对下渗量的影响；VIC-3L模型引入下渗与蓄满曲线，假设非蓄满面积上发生超渗产流；WRF-Hydro的简单水平衡产流方案（SIPM）中使用基于土壤水力实验的经验公式配合Richards方程去描述混合产流机制。

对于上述具有混合产流结构的模型，VMM模型考虑了局部产流问题，使用了下渗能力和蓄满能力曲线表达产流的空间分布，忽略了地表和地下两个过程之间的空间可变土壤水交互过程；SAC模型成功地将超渗、蓄满过程（基于具有不同蓄水量的两个土壤层结构）与Horton下渗方程（用于估算从上层到下层的下渗）耦合在一起，但却忽略了超渗、蓄满过程中的局部产流问题；GWLF、VXM和HBM也存在类似的结构缺陷，这一定程度上由概念性的模型结构决定。与此相比，网格化的分布式水文模型在表征水文变量的空间异质性上有自身的优势，但也受制于模型自身结构，例如分布式的VIC-3L虽然能够解决局部产流以及不同产流界面之间的交互，却没有提供完整的产流过程模块，缺乏壤中流的侧向运移过程。陆面模式与水文耦合似乎是最符合的求解方案，但容易被忽视的是陆面模式与水文模型在描述问题侧重点的不同决定了其耦合尺度并不一致。以往陆面模式一般耦合在气候模式中，对产流结果的精度要求较低，这也导致产流结构多通过单一的超渗或蓄满去实现产流描述。WRF-Hydro是少数的考虑混合产流机制的陆面水文耦合模式。目前极少数研究考虑了土壤结构对土壤水力性质的影响，局部水力实验参数在模式中的应用效果仍有待探究。抛开系统复杂度，WRF-Hydro中通过土壤水力实验参数对下渗过程的表述会影响模式对产流过程的模拟效果，并在模式不同模块调用中被持续放大。

发明内容

为了解决上述的技术问题，本发明提出一种适用于混合产流区的陆面水文模式的建立方法，基于近两年国际关注度高的WRF-Hydro模式无法准确模拟/预报具有区域特性的水文过程现状，针对混合产流区的产汇流机理，对模式系统进行了改进，提高了陆面水文模式在我国混合产流区上的模拟精度。

本发明的目的是提供一种适用于混合产流区的陆面水文模式的建立方法，包括以下步骤：

步骤1：采用双预热方式开展数值大气模式WRF模拟和/或预报,获得气象驱动场；

步骤2：构建WRF-Hydro-DQH陆面水文模式产流模块；

步骤3：构建WRF-Hydro-DQH陆面水文模式汇流模块；

步骤4：编译运行WRF-Hydro-DQH陆面水文模式；

步骤5：调整校正所述WRF-Hydro-DQH陆面水文模式的参数；

步骤6：进行所述数值大气模式WRF与校正后的所述WRF-Hydro-DQH陆面水文模式在混合产流区的水文变量耦合模拟和/或预报。

优选的是，所述步骤1包括对于每场降雨事件，初期的模式初始条件由对WRF陆面预热与WRF陆面-大气预热两部分组成。

在上述任一方案中优选的是，所述步骤2包括通过对混合产流区多年的历史雨洪资料的测算，采用了考虑前期土壤湿度反映下渗能力的下渗方程替换模式原描述：

其中，

在上述任一方案中优选的是，基于地形指数的空间下渗能力方程提供流域指定分辨率的网格最大下渗能力

其中，

在上述任一方案中优选的是，所述步骤3包括汇流模块在马斯京根-康吉洪水演算中增加随流速变化改变的河道渗漏损失模块，设定WRF-Hydro中的河段假定为梯形且沿整个河段的高度

其中，

在上述任一方案中优选的是，所述步骤3还包括在低流速期间对河道截面湿周

其中，

所述河道截面湿周

在上述任一方案中优选的是，所述步骤3还包括根据所述修正后湿周

其中，

在上述任一方案中优选的是，所述WRF-Hydro陆面水文模式汇流模块包括Noah-Mp陆面模块和水文汇流模块。

在上述任一方案中优选的是，所述步骤4包括给出改进产汇流结构后模式的具体实现方法以及在混合产流区内具体需要率定的参数。

在上述任一方案中优选的是，所述步骤4包括以下子步骤：

步骤41：将扣除蒸发渗漏截留后的网格净雨输入陆面模式，当净雨强度超过下渗能力时地表发生超渗产流；

步骤42：当土柱单元饱和（水位超过规定的水位上限）后，将超出土柱的水（类似于蓄满产流量概念）与下渗产生的超渗产流量在网格上累加；

步骤43：侧向运移过程中流入河道网格的水量则进行河道汇流演算。

在上述任一方案中优选的是，在所述WRF-Hydro陆面水文模式汇流模块中，原始的下渗能力

其中，

在上述任一方案中优选的是，在所述WRF-Hydro-DQH陆面水文模式中求得的下渗能力受下渗系数约束，此时空间上单位时段下渗量

其中，

在上述任一方案中优选的是，所述步骤5包括采用并行的动态多维搜索算法对敏感参数进行率定，提出相应计算节点的分配方法。

在上述任一方案中优选的是，所述分配方法为在处理器原有N个节点基础上分配N-1个节点用于解集计算，1个节点用于初始的计算任务分配与模型输出判断。

在上述任一方案中优选的是，所述步骤6包括将WRF与改进产汇流结构后的WRF-Hydro-DQH结合开展陆气耦合研究，并采用了降尺度聚解以及升尺度聚合方式对产汇流变量进行网格细分。

本发明提出了一种适用于混合产流区的陆面水文模式的建立方法，数值大气模式WRF可以提供未来短期内（一周）相对精确的气象要素模拟/预报，与WRF-Hydro耦合后能够有效的提高降雨径流模拟/预报的预见期，同时陆面水文模式系统WRF-Hydro可以预测研究区内未来的水热演变情况，这对于弄清历史极端暴雨洪水衍生发展至关重要，改进后的陆面水文模式系统WRF-Hydro-DQH能够更为精准的再现区域的径流时空分布，为相关水利部门提供区域面上包含径流过程在内的水文要素的局部特征，从而能够有针对性的开展防洪抢险救灾，同时更有利于精准制定区域防洪规划，节省大量经济与人力资源成本。

WRF-Hydro-DQH是指改进WRF-Hydro产汇流结构后适用于混合产流区的陆面水文模式。

WRF-Hydro是指近年来逐渐被国际关注的陆面水文模式。

MC方程为马斯京根-康吉河道洪水演算方程。

附图说明

图1为按照本发明的适用于混合产流区的陆面水文模式的建立方法的一优选实施例的流程图。

图2为按照本发明的适用于混合产流区的陆面水文模式的建立方法的混合产流区的陆面水文模式的一实施例的模拟/预报流程示意图。

图3为按照本发明的适用于混合产流区的陆面水文模式的建立方法的WRF嵌套方案的一实施例的示意图。

图4为按照本发明的适用于混合产流区的陆面水文模式的建立方法的一优选实施例的概化河道横截面示意图。

图5为按照本发明的适用于混合产流区的陆面水文模式的建立方法的一优选实施例的改进后写入的流域地形指数示意图。

图6为按照本发明的适用于混合产流区的陆面水文模式的建立方法的陆面水文模式的一实施例的参数聚解聚合示意图。

具体实施方式

下面结合附图和具体的实施例对本发明做进一步的阐述。

实施例一

如图1所示，执行步骤110，采用双预热方式开展数值大气模式WRF模拟和/或预报,获得气象驱动场；对于每场降雨事件，初期的模式初始条件由对WRF陆面预热与WRF陆面-大气预热两部分组成。

执行步骤120，构建WRF-Hydro-DQH陆面水文模式产流模块。通过对混合产流区多年的历史雨洪资料的测算，采用了考虑前期土壤湿度反映下渗能力的下渗方程替换模式原描述：

其中，

上面的公式即为考虑表层土壤含水量的下渗方程，该描述方式将流域看作整体，用下渗能力分布曲线来概括流域下垫面空间变化对超渗产流过程的影响，即认为分布曲线在整个流域空间上是固定的，并未考虑空间网格间下渗能力表达的差异。下面将构建地形指数与下渗能力的函数关系，实现模式对下渗能力的空间离散。

基于地形指数的空间下渗能力方程提供流域指定分辨率的网格最大下渗能力

其中

执行步骤130，构建WRF-Hydro-DQH陆面水文模式汇流模块。汇流模块在马斯京根-康吉洪水演算中增加随流速变化改变的河道渗漏损失模块，设定WRF-Hydro中的河段假定为梯形且沿整个河段的高度

其中，

在低流速期间对河道截面湿周

其中，

所述河道截面湿周

根据所述修正后湿周

其中，

执行步骤140，编译运行WRF-Hydro-DQH陆面水文模式，给出改进产汇流结构后模式的具体实现方法以及在混合产流区内具体需要率定的参数。其中，WRF-Hydro陆面水文模式汇流模块包括Noah-Mp陆面模块和水文汇流模块。在本步骤中，执行步骤141，将扣除蒸发渗漏截留后的网格净雨输入陆面模式，当净雨强度超过下渗能力时地表发生超渗产流。执行步骤142，当土柱单元饱和（水位超过规定的水位上限）后，将超出土柱的水（类似于蓄满产流量概念）与下渗产生的超渗产流量在网格上累加。执行步骤143，侧向运移过程中流入河道网格的水量则进行河道汇流演算。

在所述WRF-Hydro陆面水文模式汇流模块中，原始的下渗能力

其中，

在所述WRF-Hydro-DQH陆面水文模式中求得的下渗能力受下渗系数约束，此时空间上单位时段下渗量

其中，

执行步骤150，调整校正所述WRF-Hydro-DQH陆面水文模式的参数，采用并行的动态多维搜索算法对敏感参数进行率定，提出相应计算节点的分配方法。其中，分配方法为在处理器原有N个节点基础上分配N-1个节点用于解集计算，1个节点用于初始的计算任务分配与模型输出判断。

执行步骤160，进行所述数值大气模式WRF与校正后的所述WRF-Hydro-DQH陆面水文模式在混合产流区的水文变量耦合模拟和/或预报，将WRF与改进产汇流结构后的WRF-Hydro-DQH结合开展陆气耦合研究，并采用了降尺度聚解以及升尺度聚合方式对产汇流变量进行网格细分。

实施例二

如图2所示，一种适用于我国北方混合产流区的陆面水文模式，包括以下几个步骤：

步骤1、采用双预热方式开展数值大气模式WRF模拟/预报获得气象驱动场；如图3所示，均采用3层的WRF嵌套结构，其中预报情景下，WRF选择了GFS数据用于生成预报气象驱动场；模拟情景下，WRF选择了ERA-5数据用于生成模拟气象驱动场，并且获得的气象场是在模拟/预报开始前采用陆面预热与WRF陆面-大气预热双预热的方式获得的。

步骤2、WRF-Hydro-DQH陆面水文模式产流模块构建；对WRF-Hydro陆面水文模式的下渗模块进行了改进，主要体现在考虑了前期土壤湿度的下渗方程，并通过与地形指数概念结合实现了每个栅格范围内不同的下渗能力曲线分布。

步骤3、WRF-Hydro-DQH陆面水文模式汇流模块构建；对WRF-Hydro陆面水文模式的河道汇流方式进行了改进，主要体现在在马斯京根-康吉河道洪水演算中加入了考虑河道损失的方程，将流速与河道损失相联系，该方程是在模式原有参数基础上的补充，仅增加了一个新的河道渗漏变量。

步骤4、WRF-Hydro-DQH陆面水文模式编译运行；给出改进产汇流结构后模式的具体实现方法以及在混合产流区内具体需要率定的参数。

步骤5、WRF-Hydro-DQH陆面水文模式参数调整校正；采用了并行的动态多维搜索算法对敏感参数进行率定，给出了相应计算节点的分配方案，提高了模式运行效率以及模拟/预报准确度。

步骤6、数值大气模式WRF与改进的WRF-Hydro-DQH在我国北方混合产流区的水文变量耦合模拟/预报；将WRF与改进产汇流结构后的WRF-Hydro-DQH结合开展陆气耦合研究，并采用了降尺度聚解以及升尺度聚合方式对产汇流变量进行网格细分，降尺度比例设定一般在3-10倍。

进一步，所述步骤1中对于每场降雨事件，初期的模式初始条件由对WRF陆面预热与WRF陆面-大气预热两部分组成，其中对陆面预热的目的是为了获得相对准确的模式下边界场（尤其是土壤含水量条件），并在结束时刻生成Restart文件，对于我国北方的混合产流区而言，降水集中发生在每年的7-9月份的某几天或十几天内，其余月份的降水则占比很少，另外研究的暴雨洪水场次在时间上跨度较长，在过去研究基础上并考虑模式运行效率后，选择了15天的时间来给WRF中的陆面模式提供相对准确的初始陆面条件（陆面模式预热）；WRF陆面-大气预热时调用Restart文件中的下边界场并读取驱动数据资料中的侧边界场（温度、辐射、风场等），此部分预热时间在6-24h，对于降水前期边界和背景场均直接调用Restart文件，此部分预热的目的是协调两种边界场，确保初始WRF模式运行的稳定性，设置Restart文件的另一个目的是为了增加模式的运行效率，节约运行成本。

进一步，所述步骤2中通过对混合产流区多年的历史雨洪资料的测算，采用了考虑前期土壤湿度反映下渗能力的下渗方程：

（1）替换模式原描述：

其中，

上述方程在区域内共用一个下渗曲线，无法满足陆面水文模式更高的下渗能力空间信息要求，这里采用了基于地形指数的空间下渗能力方程提供流域指定分辨率的网格最大下渗能力

其中，

进一步，所述步骤3中汇流模块在马斯京根-康吉洪水演算中增加了随流速变化改变的河道渗漏损失模块，如图4所示，WRF-Hydro中的河段假定为梯形，并且假设沿整个河段的高度恒定（即上面介绍的平均高度

其中，

河道截面湿周

由于河床地形不均匀的影响，当流量很低时，湿周会减少，因为水只会流向不规则河道的最低的那部分。为了解释梯形河道的这种趋势,在低流速期间对湿周进行修正：

其中，

其中，

进一步，所述步骤4中WRF-Hyro包含Noah-Mp陆面模块和水文汇流模块，将扣除蒸发渗漏截留后的网格净雨输入陆面模式，当净雨强度超过下渗能力时地表发生超渗产流；当土柱单元饱和（水位超过规定的水位上限）后，将超出土柱的水（类似于蓄满产流量概念）与下渗产生的超渗产流量在网格上累加，用于随后的水通量的侧向运移和“土壤再下渗”过程（即地下汇流和坡面汇流过程）；侧向运移过程中流入河道网格的水量则进行河道汇流演算，其中产汇流的水流流向通过D8方法获得。其中下渗过程发生在Noah-Mp陆面模式中，为了方便后期参数的自动校准，对模式参数的读取接口和Noah-Mp中相应物理描述模块分别进行了代码修改。在Noah-Mp陆面代码中添加所需变量及方程（主要包括了地形指数二维数据读取以及下渗模块），其中下渗过程方案的修改以及对于常规参数（即参数是一个值的情况）的写入调用主要在module_sf_noahmplsm.F下进行，module_sf_noahmpdrv.F是用来配合调用参数以及识别不同参数的驱动文件，为了方便调参，将一维数据写入到GENPARM.TBL。

在WRF-Hydro中，原始的下渗能力

其中，

上述方程得到了下渗能力分布方案，模式中求得的下渗能力受下渗系数（水力传导系数）约束。此时空间上单位时段下渗量

其中，

具体改进时，公式（5）中降水强度

对于表层土湿的处理方式有以下两个优点：（1）半干旱地区土层深厚，地下水埋藏较深，下渗锋面很难达到包气带底部，使用上层0-40cm（一二层）或0-100cm（一二三层，本研究中使用）内土壤湿度的和来表示降雨的下渗锋面在产流过程中可能到达土层深度内的土壤湿度，减少了中间循环变量，排除了原方程中下层土壤层湿度对超渗过程的影响。（2）放大前2层/3层土壤湿度的影响，未来有利于遥感以及同化对模式超渗过程的校准。

进一步，所述步骤5中为了加快计算资源获取高质量解集的收敛速度，在率定过程中采用了并行动态多维搜索算法（Parallel Dynamically Dimensioned Search，P-DDS），由于陆面水文模式包含大量参数，这里经过前期参数敏感性分析后选择了包含模式原有以及改进后新增加的几个关键参数变量，所需率定的参数已在步骤4种给出，P-DDS中加入了对于计算节点的分配，即在处理器原有N个节点基础上分配N-1个节点用于解集计算（计算节点），1个节点用于初始的计算任务分配与模型输出判断（处理节点）。

进一步，所述步骤6中采用WRF的3层嵌方式（每层嵌套比例为1:3）并取最内层区域用于为WRF-Hydro-DQH提供包括降水、辐射、湿度、温度、气压、风速等在内的模式输入，如图6所示，为了解决陆面水文分辨率不匹配问题，降尺度聚解以及升尺度聚合方式对产汇流变量进行网格细分，降尺度比例设定一般在3-10倍。

实施例三

本发明所采用的技术方案是针对北方超渗产流区暴雨洪水突发性强，汇流时间快等特点，对WRF-Hydro的产汇流特性进行了改进，形成一种适用于我国北方混合产流区的陆面水文模式，按照以下步骤实施：采用双预热方式开展数值大气模式WRF模拟/预报获得气象驱动场；WRF-Hydro-DQH陆面水文模式产流模块改建；WRF-Hydro陆面水文模式汇流模块构建；WRF-Hydro-DQH陆面水文模式编译运行；WRF-Hydro-DQH陆面水文模式参数调整校正；数值大气模式WRF与改进的WRF-Hydro-DQH在我国北方混合产流区的水文变量耦合模拟/预报。

（1）采用3层的WRF嵌套结构，双预热方式开展数值大气模式WRF模拟/预报获得气象驱动场：对于每场降雨事件，初期的模式初始条件由对WRF陆面预热与WRF陆面-大气预热两部分组成，其中对陆面预热的目的是为了获得相对准确的模式下边界场（尤其是土壤含水量条件），选择了15天的时间来给WRF中的陆面模式提供相对准确的初始陆面条件（陆面模式预热），并在结束时刻生成Restart文件；WRF陆面-大气预热时调用Restart文件中的下边界场并读取驱动数据资料中的侧边界场（温度、辐射、风场等），此部分预热时间在6-24h，对于降水前期边界和背景场均直接调用Restart文件，此部分预热的目的是协调两种边界场，确保初始WRF模式运行的稳定性，设置Restart文件的另一个目的是为了增加模式的运行效率，节约运行成本。

（2）WRF-Hydro-DQH陆面水文模式产流模块改建：通过对混合产流区多年的历史雨洪资料的测算，采用了考虑前期土壤湿度反映下渗能力的下渗方程（1）替换模式原描述：

式中

上述方程在区域内共用一个下渗曲线，无法满足陆面水文模式更高的下渗能力空间信息要求，这里采用了基于地形指数的空间下渗能力方程提供流域指定分辨率的网格最大下渗能力

（3）WRF-Hydro-DQH陆面水文模式汇流模块改建：汇流模块在马斯京根-康吉洪水演算中增加了随流速变化改变的河道渗漏损失模块，WRF-Hydro中的河段假定为梯形，并且假设沿整个河段的高度恒定（即上面介绍的平均高度

河道截面湿周

由于河床地形不均匀的影响，当流量很低时，湿周会减少，因为水只会流向不规则河道的最低的那部分。为了解释梯形河道的这种趋势,在低流速期间对湿周进行修正：

其中，

上式中

（4）WRF-Hydro-DQH陆面水文模式编译运行：WRF-Hyro包含Noah-Mp陆面模块和水文汇流模块，其中下渗过程发生在Noah-Mp陆面模式中，为了方便后期参数的自动校准，对模式参数的读取接口和Noah-Mp中相应物理描述模块分别进行了代码修改。在Noah-Mp陆面代码中添加所需变量及方程（主要包括了地形指数二维数据读取以及下渗模块），其中下渗过程方案的修改以及对于常规参数（即参数是一个值的情况）的写入调用主要在module_sf_noahmplsm.F下进行，module_sf_noahmpdrv.F是用来配合调用参数以及识别不同参数的驱动文件，为了方便调参，将一维数据写入到GENPARM.TBL。

在WRF-Hydro中，公式（5）中降水强度

（5）WRF-Hydro-DQH陆面水文模式参数调整校正：为了加快计算资源获取高质量解集的收敛速度，在率定过程中采用了并行动态多维搜索算法（Parallel DynamicallyDimensioned Search，P-DDS），由于陆面水文模式包含大量参数，这里经过前期参数敏感性分析后选择了包含模式原有以及改进后新增加的几个关键参数变量，所需率定的参数已在步骤4种给出，P-DDS中加入了对于计算节点的分配，即在处理器原有N个节点基础上分配N-1个节点用于解集计算（计算节点），1个节点用于初始的计算任务分配与模型输出判断（处理节点）。

（6）数值大气模式WRF与改进的WRF-Hydro-DQH在我国北方混合产流区的水文变量耦合模拟/预报：采用WRF的3层嵌套方式（每层嵌套比例为1:3）并取最内层区域用于为WRF-Hydro-DQH提供包括降水、辐射、湿度、温度、气压、风速等在内的模式输入，为了解决陆面水文分辨率不匹配问题，降尺度聚解以及升尺度聚合方式对产汇流变量进行网格细分，降尺度比例设定一般在3-10倍。

为了更好地理解本发明，以上结合本发明的具体实施例做了详细描述，但并非是对本发明的限制。凡是依据本发明的技术实质对以上实施例所做的任何简单修改，均仍属于本发明技术方案的范围。本说明书中每个实施例重点说明的都是与其它实施例的不同之处，各个实施例之间相同或相似的部分相互参见即可。对于系统实施例而言，由于其与方法实施例基本对应，所以描述的比较简单，相关之处参见方法实施例的部分说明即可。