CN111931440A - 一种考虑地表积水深变化的坡面多流向漫流分析方法 - Google Patents

Abstract

本发明涉及一种考虑地表积水深变化的坡面多流向漫流分析方法，包括：地表积水深的确定、坡面动态汇流路径确定、坡面多流向水量分配权重确定和坡面多流向漫流水量确定。本发明针对平原洼地区地表起伏小，积水面积大、变化快的特点，实现了对地表积水深动态变化及其与周边单元的水量流入流出关系转化的准确刻画，解决了常规数字高程模型应用失真的问题，精细模拟了坡面多流向漫流过程，有效改善了平原低洼区地表产汇流过程模拟。该方法面向平原洼地区水文模拟，相关结果可为平原洼地区水文模型构建与应用提供基础支撑。

Description

一种考虑地表积水深变化的坡面多流向漫流分析方法

技术领域

本发明涉及水利工程技术领域，尤其涉及平原低洼区地表产汇流过程模拟，具体为一种考虑地表积水深变化的坡面多流向漫流分析方法。

背景技术

平原低洼区具有地势低平、地表水流流向随机的特点，流向多不单一且不固定，易受地表积水变化影响而发生改变，长期以来对平原低洼区地表产汇流过程的准确模拟是个难题。为此，坡面多流向漫流理论被提出并不断完善，其基本思想是认为水流应向邻域中所有高程较低的方向进行分配。目前，相关理论方法均是依靠常规数字地形数据(一般为数字高程模型“DEM”)，研究坡度、地表形态、粗糙度等对地表水流多流向流量分配的影响。然而，从平原区特别是平原湿地沼泽区来说，当地洼地星罗棋布，地表起伏小，积水面积大、变化快，采用传统不变的数字高程模型(DEM)难以准确反映地形地貌变化，会对流向和流量确定产生较大的干扰，影响地表产汇流过程模拟效果。以洼地蓄滞水过程为例，随着洼地积水，水流流出单元与汇入单元间的水力梯度逐渐得到消弭，直接影响到单元出流水量。

发明内容

本发明所要解决的技术问题，在于克服现有平原低洼区坡面多流向漫流过程模拟技术存在的缺陷，构建考虑地表积水深度变化的动态数字地形，并基于此创造性提出包含坡面动态汇流路径确定和多流向漫流水量计算的坡面汇流计算方法。该发明可有效刻画平原洼地积水变化及其对地表产汇流过程的影响，提高平原低洼区地表产汇流过程模拟精度。

本发明的目的是这样实现的：

一种考虑地表积水深变化的坡面多流向漫流分析方法，所述方法包括四个部分：地表积水深确定、坡面动态汇流路径确定、坡面多流向水量分配权重确定和坡面多流向漫流水量确定，各部分的具体步骤如下：

地表积水深确定：

步骤一：获取研究区内气象站点在模拟期内的气象数据，所述气象数据包括逐日降水、气温、风速、日照、湿度，选择参证气象站，利用反距离平方加权插值算法将参证气象站点的气象数据空间展布到研究区模拟栅格单元上；

步骤二：采用Green-Ampt入渗模型模拟计算栅格单元的降雨日累计入渗量；

步骤三：利用Penman公式模拟计算栅格单元地表当日最大可蒸发量；

步骤四：利用水量平衡方程计算模拟单元日积水深，计算公式如下：

式中，H为日积水深/mm；H₀为当日初始积水深/mm；P为当日降水深/mm；E为当日最大可蒸发量/mm；F为降雨日累计入渗量/mm；

坡面动态汇流路径确定的具体步骤：

步骤1：以日为汇流模拟时段，赋予研究区内每个模拟栅格单元5个基本属性：坐标x、坐标y、地形高程z、日积水深H、是否已参与汇流计算；已参与汇流计算标记为ture，否则标记为False；

步骤2：逐汇流模拟时段内，引入当前模拟栅格单元日积水深H，将其与模拟栅格单元的地形高程z相加，以此作为坡面漫流、河道汇流和漫溢计算所需的“虚拟高程”，记作DH；

步骤3：遍历模拟栅格单元，寻找DH最大值所在的模拟栅格单元，确定流向和出流量，然后更新模拟栅格单元的DH，并将结束汇流的模拟栅格单元标记为True，即表示本时段内该栅格不再参与汇流计算；

步骤4：逐次遍历，直至所有模拟栅格单元的汇流路径和出流量全部计算完成；

坡面多流向水量分配权重确定的具体步骤：

步骤1)：比较模拟栅格单元和该单元相邻的8个栅格单元的DH大小，若模拟栅格单元DH大于相邻栅格单元，则认为水流流向该临近栅格单元，研究认为模拟栅格单元的水流应向邻域中所有DH较低的方向进行分配，由此对流向进行编码；

步骤2)：结合地表Manning糙率，确定模拟栅格单元在临近多个DH较低栅格单元上的水量分配权重；

坡面多流向漫流水量确定的具体步骤：

步骤Ⅰ：使用运动波模型计算水流流出栅格单元在本汇流时段内的流量，记为可出流量；计算公式如下：

S_f＝S₀ (20)

式中，A为过流断面面积，单位为m²；Q为断面流量，单位为m³/s，即为可出流量；q_L为单宽流入量，单位为m³/s/m；S_f为摩擦比降；S₀为地面坡降或河道坡降；R为过流断面水力半径，单位为m；n为坡面曼宁糙率系数；t为时间；x为河段长度；

步骤Ⅱ：计算水流流出栅格单元与流入栅格单元“虚拟高程”出现持平时的出流量d_c，以此作为该计算单元在本汇流时段的最大出流量；

步骤Ⅲ：取Q与d_c两者之间的较小值作为单元出流量最终模拟值，并按水量分配权重对各方向单元流入量进行修正。

进一步的优化，步骤一中利用反距离平方加权插值算法将参证气象站点的气象数据空间展布到研究区模拟栅格单元上，计算公式如下：

式中，z和M分别表示气象站的高程和气象数据值，dM/dz为气象数据随高程增加的线性变化率；m表示参证气象站点的总数；D和D_i分别为模拟栅格单元和第i个参证气象站点的气象数据值；z₀和z_i分别为模拟栅格单元和第i个参证气象站点的高程；λ_i为第i个参证气象站点所占权重；x₀和x_i分别为模拟栅格单元和第i个参证气象站点的横坐标；y₀和y_i分别为模拟栅格单元和第i个参证气象站点的纵坐标；d_i为模拟栅格单元和第i个参证气象站点的平面距离；

进一步的优化，步骤二中采用Green-Ampt入渗模型模拟计算栅格单元的降雨日累计入渗量，计算方法如下：

当入渗湿润锋到达第m层土壤时，入渗能力f为：

若湿润锋进入第m-1层土壤时地表持续积水，降雨日累计入渗量F由式(5)计算；若前一时段t_n-1地表无积水，当前时段t_n开始积水，降雨日累计入渗量F则用式(6)计算；

其中，

F_p＝A_m-1(I_p/k_m-1)-B_m-1 (10)

t_p＝t_n-1+(F_p-F_n-1)/I_p (11)

△θ_i＝θ_s,i-θ_i,0 (12)

式中，m为湿润锋最终到达的土层数；i表示湿润锋当前所在土层数；k_m为湿润锋所在m层时的土壤导水系数，单位为mm/day；k_i为第i层土壤导水系数，单位为mm/day；Δθ_i为第i层土壤含水量的变化量，cm³/cm³；Δθ_m为湿润锋所在m层的土壤含水量变化量，cm³/cm³；A_m-1和B_m-1为参数，无量纲；F_m-1为m-1层以上的土壤累积入渗量，mm；F_P为当前时段开始积水时的累积入渗量，单位为mm；t为时刻；t_n-1为第n-1时段，min；F_n-1为第n-1时段的土壤累积入渗量，mm；t_m-1为湿润锋进入m-1和m层土壤界面时的时间；t_p为开始积水时刻；L_i为第i层土厚，单位为mm；SW_m为第m层湿润锋处的毛管吸力，单位为mm；θ表示土壤饱和含水量，cm³/cm³；I_p为积水开始时的雨强，mm/min；θ_s,i为第i层土壤饱和含水量，cm³/cm³；θ_i,0为土壤初始含水量，cm³/cm³。

进一步的优化，步骤三：利用Penman公式模拟计算栅格单元地表当日最大可蒸发量，计算公式如下：

式中：E为当日最大可蒸发量，单位为mm；Δ为饱和水汽压曲线斜率，单位为kPa/℃；RN为净放射量，单位为MJ/m；G为地中热传导，单位为MJ/m²；ρ_a为空气密度，单位为kg/m³；C_P为空气定压比热，单位为J/kg/℃；δ_e为实际与饱和水蒸气压之差，单位为kPa；r_a为空气动力学阻抗，单位为s/m；λ为水体的气化潜热，MJ/kg；γ为空气湿度常数，kPa/℃；P_a为大气压，kPa。

进一步的，步骤2)：确定模拟栅格单元在临近多个DH较低栅格单元上的水量分配权重，计算公式如下：

式中，i和j表示两个不同栅格单元；DH_i和DH_j分别为栅格单元i和栅格单元j的虚拟高程；w_ji表示从DH较高栅格单元j到相邻DH较低单元i的分配水量权重；s_ji表示栅格单元j到栅格单元i间的水力梯度；n_ji为栅格单元j到栅格单元i间的坡面曼宁糙率系数；n_i和n_j分别为栅格单元i和栅格单元j的坡面曼宁糙率系数；x_i和x_j分别为栅格单元i和栅格单元j的横坐标；y_i和y_j分别为栅格单元i和栅格单元j的的纵坐标。

本发明产生的有益效果是：

本发明一种考虑地表积水深变化的坡面多流向漫流分析方法，在模拟坡面多流向漫流过程中考虑地表积水深的变化及其影响，可精细刻画平原洼地单元的积水变化及其与周边单元的水量流入流出关系转化，提高平原低洼区地表产汇流过程模拟精度。

本发明方法创造性提出了坡面动态汇流路径确定和坡面多流向漫流水量确定的具体流程，准确描述了平原洼地区域地表起伏小，积水面积大、变化快的特点，深入刻画了洼地对坡面径流的调蓄作用以及平原地表径流多流向漫散过程，实现了对地表积水深动态变化及其与周边单元的水量流入流出关系转化的准确刻画，解决了常规数字高程模型应用失真的问题，有效改善了平原低洼区地表产汇流过程模拟。该方法面向平原洼地区水文模拟，相关结果可为平原洼地区水文模型构建与应用提供基础支撑

附图说明

下面结合附图和具体实施例，对本发明作进一步详细说明。

图1是本发明的实施例1所述方法的流程图；

图2是本发明考虑地表积水深的动态汇流路径确定示意图；

图3是本发明不同流向对应编码示意图；

图4是本发明地表漫流多流向水量分配计算示意图。

具体实施方式

本实施例是一种考虑地表积水深变化的坡面多流向漫流分析方法，计算流程如图1所示。本实施例所述方法的计算过程如下：

第一步：计算确定研究区模拟栅格单元的地表积水深

1、获取研究区内气象站点在模拟期内的逐日降水、气温、风速、日照、湿度等气象数据，然后利用反距离平方加权插值算法将气象站点数据空间展布到研究区模拟栅格单元上。计算公式如下：

式中：z和M分别表示气象站的高程和气象数据值，dM/dz为气象数据随高程增加的线性变化率，若采用平面插值则dM/dz＝0；下标0表示待插值栅格点，i表示参证气象站点；m表示参证气象站点的总数；D和D_i分别为模拟栅格单元和第i个参证气象站点的数据值；z₀和z_i分别为模拟栅格单元和第i个参证气象站点的高程；λ_i为第i个参证气象站点所占权重；x₀和x_i分别为模拟栅格单元和第i个参证气象站点的横坐标；y₀和y_i分别为模拟栅格单元和第i个参证气象站点的纵坐标；d_i为模拟栅格单元和第i个参证气象站点的平面距离。

dM/dz为气象数据随高程增加的线性变化率，z和D分别为参证气象站点的高程和数据值；λ_i为第i个参证气象站点所占权重；d_i为第i个参证气象站点与待插值栅格点的平面投影距离。

2、采用Green-Ampt入渗模型模拟计算栅格单元的降雨日累计入渗量。当入渗湿润锋到达第m层土壤时，入渗能力f为：

若湿润锋进入第m-1层土壤时地表持续积水，降雨日累计入渗量F由式(5)计算；若前一时段t_n-1地表无积水，当前时段t_n开始积水，降雨日累计入渗量F则用式(6)计算。

其中，

F_p＝A_m-1(I_p/k_m-1)-B_m-1 (10)

t_p＝t_n-1+(F_p-F_n-1)/I_p (11)

△θ_i＝θ_s,i-θ_i,0 (12)

式中，m为湿润锋最终到达的土层数；i表示湿润锋当前所在土层数；k_m为湿润锋所在m层时的土壤导水系数，单位为mm/day；k_i为第i层土壤导水系数，单位为mm/day；Δθ_i为第i层土壤含水量的变化量，cm³/cm³；Δθ_m为湿润锋所在m层的土壤含水量变化量，cm³/cm³；A_m-1和B_m-1为参数，无量纲；F_m-1为m-1层以上的土壤累积入渗量，mm；F_P为当前时段开始积水时的累积入渗量，单位为mm；t为时刻；t_n-1为第n-1时段，min；F_n-1为第n-1时段的土壤累积入渗量，mm；t_m-1为湿润锋进入m-1和m层土壤界面时的时间；t_p为开始积水时刻；L_i为第i层土厚，单位为mm；SW_m为第m层湿润锋处的毛管吸力，单位为mm；θ表示土壤饱和含水量，cm³/cm³；I_p为积水开始时的雨强，mm/min；θ_s,i为第i层土壤饱和含水量，cm³/cm³；θ_i,0为土壤初始含水量，cm³/cm³。

3、利用Penman公式模拟计算栅格单元地表当日最大可蒸发量。计算公式如下：

式中，Δ为饱和水汽压曲线斜率(kPa/℃)；RN为净放射量(MJ/m²)；G为地中热传导(MJ/m²)；ρ_a为空气密度(kg/m3)；C_P为空气定压比热(J/kg/℃)；δ_e为实际与饱和水蒸气压之差(kPa)；r_a为空气动力学阻抗(s/m)；λ为水体的气化潜热(MJ/kg)；γ为空气湿度常数(kPa/℃)；P_a为大气压(kPa)。

4、利用水量平衡方程计算模拟单元日积水深H。计算公式如下：

式中，H为日积水深/mm；H₀为当日初始积水深/mm；P为当日降水深/mm；E为当日最大可蒸发量/mm；F为降雨日累计入渗量/mm。

第二步：确定坡面动态汇流路径

1、以日为汇流模拟时段，赋予研究区内每个模拟栅格单元5个基本属性：坐标x、坐标y、地形高程z、地表积水深H、是否已参与汇流计算(True或False，时段开始时初始默认为False)；

2、逐汇流模拟时段内，引入当前模拟栅格单元的地表积水深H，将其与模拟栅格单元的地形高程z相加，以此作为坡面漫流、河道汇流和漫溢计算所需的“虚拟高程”，记作DH(图2)；

3、遍历模拟栅格单元，寻找DH最大值所在的模拟栅格单元，确定流向和出流量(见“坡面多流向水量分配权重确定和坡面多流向漫流水量确定”)，然后更新模拟栅格单元的DH，并将结束汇流的模拟栅格单元标记为True，即表示本时段内该栅格不再参与汇流计算。

4、逐次遍历，直至所有模拟栅格单元的汇流路径和出流量全部计算完成。

第三步：确定坡面多流向水量分配权重确定

1、比较模拟栅格单元和该单元相邻的8个栅格单元的DH大小，若模拟栅格单元DH大于相邻栅格单元，则认为水流流向该临近栅格单元，研究认为模拟栅格单元的水流应向邻域中所有DH较低的方向进行分配，由此对流向进行编码(图3)。

2、确定模拟栅格单元在临近多个DH较低栅格单元上的水量分配权重。计算公式如下：

式中，i和j表示两个不同栅格单元；DH_i和DH_j分别为栅格单元i和栅格单元j的虚拟高程；w_ji表示从DH较高栅格单元j到相邻DH较低单元i的分配水量权重；s_ji表示栅格单元j到栅格单元i间的水力梯度；n_ji为栅格单元j到栅格单元i间的坡面曼宁糙率系数；n_i和n_j分别为栅格单元i和栅格单元j的坡面曼宁糙率系数；x_i和x_j分别为栅格单元i和栅格单元j的横坐标；y_i和y_j分别为栅格单元i和栅格单元j的的纵坐标。

第四步：确定坡面多流向漫流水量

1、使用运动波模型计算水流流出栅格单元在本汇流时段内的可出流量。计算公式如下：

S_f＝S₀ (20)

式中，A为过流断面面积，单位为m²；Q为断面流量，单位为m³/s；q_L为单宽流入量，单位为m³/s/m；S_f为摩擦比降；S₀为地面坡降或河道坡降；R为过流断面水力半径，单位为m；n为坡面曼宁糙率系数；t为时间；x为河段长度；

2、计算水流流出栅格单元与流入栅格单元“虚拟高程”出现持平时的出流量d_c，以此作为该计算单元在本汇流时段的最大出流量。计算时已知时段初流出栅格单元和流入栅格单元的积水深，以及流向不同栅格单元的水量分配权重，通过建立水量平衡方程即可得到不同流向上的流入量以及栅格的流出量。以双向流举例如下(图4)：

式中，D_c0、D₁₀分别表示汇流时段初流出单元与方向1流入单元的积水深，由于单元面积相等，以积水深作为单元水量；w₁、w₂分别表示两方向上的水量分配权重；d_c表示时段末水流流出栅格单元与流入栅格单元“虚拟高程”出现持平时的出流量；d₁、d₂分别表示时段末出流方向1和2上的栅格单元入流量。由于汇流时段初积水和水量分配权重已知，可求得d_c、d₁、d₂。

3、取Q与d_c两者之间的较小值作为单元出流量最终模拟值，并按水量分配权重w₁、w₂对各方向单元流入量进行修正，即各方向单元流入量等于单元出流量最终模拟值乘以各方向上水量分配权重。

最后应说明的是，以上仅用以说明本发明的技术方案而非限制，尽管参照较佳布置方案对本发明进行了详细说明，本领域的普通技术人员应当理解，可以对本发明的技术方案(比如各种公式的运用、步骤的先后顺序等)进行修改或者等同替换，而不脱离本发明技术方案的精神和范围。

Claims (5)

1.一种考虑地表积水深变化的坡面多流向漫流分析方法，其特征在于：所述方法包括四个部分：地表积水深确定、坡面动态汇流路径确定、坡面多流向水量分配权重确定和坡面多流向漫流水量确定，各部分的具体步骤如下：

地表积水深确定：

步骤一：获取研究区内气象站点在模拟期内的气象数据，所述气象数据包括逐日降水、气温、风速、日照、湿度，选择参证气象站，利用反距离平方加权插值算法将参证气象站点的气象数据空间展布到研究区模拟栅格单元上；

步骤二：采用Green-Ampt入渗模型模拟计算栅格单元的降雨日累计入渗量；

步骤三：利用Penman公式模拟计算栅格单元地表当日最大可蒸发量；

步骤四：利用水量平衡方程计算模拟单元日积水深，计算公式如下：

式中，H为日积水深/mm；H₀为当日初始积水深/mm；P为当日降水深/mm；E为当日最大可蒸发量/mm；F为降雨日累计入渗量/mm；

坡面动态汇流路径确定的具体步骤：

步骤1：以日为汇流模拟时段，赋予研究区内每个模拟栅格单元5个基本属性：坐标x、坐标y、地形高程z、日积水深H、是否已参与汇流计算；已参与汇流计算标记为ture，否则标记为False；

步骤2：逐汇流模拟时段内，引入当前模拟栅格单元日积水深H，将其与模拟栅格单元的地形高程z相加，以此作为坡面漫流、河道汇流和漫溢计算所需的“虚拟高程”，记作DH；

步骤3：遍历模拟栅格单元，寻找DH最大值所在的模拟栅格单元，确定流向和出流量，然后更新模拟栅格单元的DH，并将结束汇流的模拟栅格单元标记为True，即表示本时段内该栅格不再参与汇流计算；

步骤4：逐次遍历，直至所有模拟栅格单元的汇流路径和出流量全部计算完成；

坡面多流向水量分配权重确定的具体步骤：

步骤1)：比较模拟栅格单元和该单元相邻的8个栅格单元的DH大小，若模拟栅格单元DH大于相邻栅格单元，则认为水流流向该临近栅格单元，研究认为模拟栅格单元的水流应向邻域中所有DH较低的方向进行分配，由此对流向进行编码；

步骤2)：结合地表Manning糙率，确定模拟栅格单元在临近多个DH较低栅格单元上的水量分配权重；

坡面多流向漫流水量确定的具体步骤：

步骤Ⅰ：使用运动波模型计算水流流出栅格单元在本汇流时段内的流量，记为可出流量；计算公式如下：

S_f＝S₀ (20)

式中，A为过流断面面积，单位为m²；Q为断面流量，单位为m³/s，即为可出流量；q_L为单宽流入量，单位为m³/s/m；S_f为摩擦比降；S₀为地面坡降或河道坡降；R为过流断面水力半径，单位为m；n为坡面曼宁糙率系数；t为时间；x为河段长度；

步骤Ⅱ：计算水流流出栅格单元与流入栅格单元“虚拟高程”出现持平时的出流量d_c，以此作为该计算单元在本汇流时段的最大出流量；

步骤Ⅲ：取Q与d_c两者之间的较小值作为单元出流量最终模拟值，并按水量分配权重对各方向单元流入量进行修正。

2.根据权利要求1所述的考虑地表积水深变化的坡面多流向漫流分析方法，其特征在于：步骤一中利用反距离平方加权插值算法将参证气象站点的气象数据空间展布到研究区模拟栅格单元上，计算公式如下：

式中，z和M分别表示气象站的高程和气象数据值，dM/dz为气象数据随高程增加的线性变化率；m表示参证气象站点的总数；D和D_i分别为模拟栅格单元和第i个参证气象站点的气象数据值；z₀和z_i分别为模拟栅格单元和第i个参证气象站点的高程；λ_i为第i个参证气象站点所占权重；x₀和x_i分别为模拟栅格单元和第i个参证气象站点的横坐标；y₀和y_i分别为模拟栅格单元和第i个参证气象站点的纵坐标；d_i为模拟栅格单元和第i个参证气象站点的平面距离。

3.根据权利要求1所述的考虑地表积水深变化的坡面多流向漫流分析方法，其特征在于：步骤二中采用Green-Ampt入渗模型模拟计算栅格单元的降雨日累计入渗量，计算方法如下：

当入渗湿润锋到达第m层土壤时，入渗能力f为：

若湿润锋进入第m-1层土壤时地表持续积水，降雨日累计入渗量F由式(5)计算；若前一时段t_n-1地表无积水，当前时段t_n开始积水，降雨日累计入渗量F则用式(6)计算；

其中，

F_p＝A_m-1(I_p/k_m-1)-B_m-1 (10)

t_p＝t_n-1+(F_p-F_n-1)/I_p (11)

△θ_i＝θ_s,i-θ_i,0 (12)

式中，m为湿润锋最终到达的土层数；i表示湿润锋当前所在土层数；k_m为湿润锋所在m层时的土壤导水系数，单位为mm/day；k_i为第i层土壤导水系数，单位为mm/day；Δθ_i为第i层土壤含水量的变化量，cm³/cm³；Δθ_m为湿润锋所在m层的土壤含水量变化量，cm³/cm³；A_m-1和B_m-1为参数，无量纲；F_m-1为m-1层以上的土壤累积入渗量，mm；F_P为当前时段开始积水时的累积入渗量，单位为mm；t为时刻；t_n-1为第n-1时段，min；F_n-1为第n-1时段的土壤累积入渗量，mm；t_m-1为湿润锋进入m-1和m层土壤界面时的时间；t_p为开始积水时刻；L_i为第i层土厚，单位为mm；SW_m为第m层湿润锋处的毛管吸力，单位为mm；θ表示土壤饱和含水量，cm³/cm³；I_p为积水开始时的雨强，mm/min；θ_s,i为第i层土壤饱和含水量，cm³/cm³；θ_i,0为土壤初始含水量，cm³/cm³。

4.根据权利要求1所述的考虑地表积水深变化的坡面多流向漫流分析方法，其特征在于：

步骤三：利用Penman公式模拟计算栅格单元地表当日最大可蒸发量，计算公式如下：

式中：E为当日最大可蒸发量，单位为mm；Δ为饱和水汽压曲线斜率，单位为kPa/℃；RN为净放射量，单位为MJ/m；G为地中热传导，单位为MJ/m²；ρ_a为空气密度，单位为kg/m³；C_P为空气定压比热，单位为J/kg/℃；δ_e为实际与饱和水蒸气压之差，单位为kPa；r_a为空气动力学阻抗，单位为s/m；λ为水体的气化潜热，MJ/kg；γ为空气湿度常数，kPa/℃；P_a为大气压，kPa。

5.根据权利要求1所述的考虑地表积水深变化的坡面多流向漫流分析方法，其特征在于：

步骤2)：确定模拟栅格单元在临近多个DH较低栅格单元上的水量分配权重，计算公式如下：

式中，i和j表示两个不同栅格单元；DH_i和DH_j分别为栅格单元i和栅格单元j的虚拟高程；w_ji表示从DH较高栅格单元j到相邻DH较低单元i的分配水量权重；s_ji表示栅格单元j到栅格单元i间的水力梯度；n_ji为栅格单元j到栅格单元i间的坡面曼宁糙率系数；n_i和n_j分别为栅格单元i和栅格单元j的坡面曼宁糙率系数；x_i和x_j分别为栅格单元i和栅格单元j的横坐标；y_i和y_j分别为栅格单元i和栅格单元j的的纵坐标。

Images