CN109118718B - 泥石流发生降雨i-d曲线阈值构建方法、流域泥石流预警方法 - Google Patents

Abstract

本发明公开泥石流发生降雨I‑D曲线阈值构建方法、泥石流预警方法。本发明I‑D曲线阈值构建方法首先获取用于Richard模型与Van Genuchten方程数值模拟所需的流域下垫面数据，在给定降雨强度条件下启动水文模拟过程，沿时间序列判断流域各栅格稳定性与产流状况，计算水土混合物容重，记录符合控制点的数据以拟合I‑D曲线，最终得到前期雨量条件下一组I‑D阈值。本发明同时提供相应的流域泥石流预警方法。本发明方法摒弃通过历史数据资料统计分析拟合曲线的构思，实现了从测量泥石流本身特征动态指标出发构建泥石流发生条件阈值的突破。方法能够有效解决无资料流域采用I‑D曲线阈值实施泥石流预警的问题。

Description

泥石流发生降雨I-D曲线阈值构建方法、流域泥石流预警方法

技术领域

本发明涉及降雨诱发泥石流的预警方法，特别是涉及一种借助降雨 强度-降雨历时曲线阈值(I-D曲线阈值)实施的泥石流预警方法，属于 泥石流灾害防治技术领域。

背景技术

降雨是触发泥石流的关键因子。学者们采用降雨参数(或参数之间 的组合)对应的泥石流事件来研究触发泥石流的降雨条件经验阈值。借 助已建立的阈值关系和持续的降雨监测数据，实现对泥石流预报的目的。 目前研究最为广泛的经验关系是Cain(1980)最先提出的降雨强度-降雨 持时阈值(I-D曲线阈值)，其一般的表达方式为I＝δD^β，其中I表示降 雨的平均雨强(mm/h)(y-axis)，D表示降雨的持续时间(h)(x-axis)，δ与β 为经验系数。I-D曲线阈值可以揭示触发区域内泥石流所需的最低降雨强 度与降雨持时，利用触发泥石流的下包络线判定是否有泥石流发生。泥 石流发生的基本判定条件为：高于下包络线的降雨强度-降雨持时组合都 认为可以触发泥石流。

现有技术中，建立I-D曲线阈值的普遍方式均需要依赖大量包含降雨 数据与泥石流灾害事件的历史数据，通过对这些历史数据的统计分析拟 合得到一个I-D曲线方程。这类方法主要存在五方面缺陷：其一、为采用 统计方法拟合I-D曲线必须流域内降雨数据监测数据与流域内泥石流活 动观测数据，在实际工程中，由于泥石流流域往往属于经济生产条件落 后地区，流域内没有建立起完善的气象灾害观测点，缺乏流域内降雨数 据监测数据与流域内泥石流活动观测数据，因而无法构建本流域适用的 I-D曲线阈值，也就无法利用其实施泥石流预警。其二、I-D曲线阈值是 一条区域性的经验曲线。因而理论上统计得到的泥石流点越集中，其气 候、地质因素的差异越小，所得到I-D曲线阈值的变异性越小，基于I-D 曲线阈值的泥石流预报就越准确。但实际的情况是由于拟合I-D曲线方程 所用数据的往往精度较低，因而造成所建立的过程可重复性较差。其三、 依靠经验方法建立的I-D曲线阈值没有考虑下垫面对泥石流发生所产生 的作用，包括对泥石流形成中所需的固体物源补给与降雨形成的径流量 两个必要条件充分考虑。其四、统计分析方法拟合I-D曲线方程必须依赖 于大量包含降雨数据与泥石流灾害事件的历史数据，对于无资料地区或 者资料不足的地区，就难以使用该方法实施泥石流预警。其五、现有技 术中，每套基于I-D曲线阈值的泥石流监测/预警方案都仅设置一条I-D曲 线阈值。当降雨强度-降雨持时组合位于I-D曲线阈值之下的坐标区域时， 表示无泥石流，当降雨强度-降雨持时组合位于I-D曲线阈值之上的坐标 空间时，表示有泥石流。这样的泥石流预警方案在原理上忽视了I-D曲线阈值的约束条件，在实际预警中忽视了I-D曲线阈值的预警边界。

发明内容

本发明的目的就是针对现有技术的不足，提供一种基于泥石流物理 过程的I-D曲线阈值构建方法，以及利用该方法实现的降雨诱发泥石流 的预警方法。

为实现上述目的，本发明首先提供一种泥石流发生降雨I-D曲线阈 值构建方法，其技术方案如下：

一种泥石流发生降雨I-D曲线阈值构建方法，其特征在于：

首先，获取泥石流流域数据，所述泥石流流域数据包括流域内用于 一维非饱和Richard数学模型与Van Genuchten方程水文过程数值模拟所 需的数据；确定水土耦合后的混合物容重控制点ρ′、确定水文模拟计算 的前期雨量Ar′；

其次，在给定降雨强度I的条件下启动水文模拟过程，联合一维非 饱和Richard数学模型、Van Genuchten方程进行水文过程数值模拟，确 定用于计算流域边坡稳定性安全系数Fs的关键水文参数；

再次，计算泥石流流域每一栅格的Fs值，判断每一栅格稳定性，沿 时间序列试算当前流域内水土耦合后的混合物容重ρ，捕捉试算结果混 合物容重ρ＝水土混合物容重控制点ρ′的点，确定该点为拟合数据点，记 录拟合数据点对应的降雨历时D、降雨强度I；

最后，利用所有拟合数据点拟合得到流域在前期雨量Ar′条件下的 I-D曲线I＝δD^β，在坐标系上绘制I-D曲线阈值。

上述泥石流发生降雨I-D曲线阈值构建方法是一种基于泥石流物理 过程的I-D曲线阈值构建方法，其主要技术原理在于：一、在获得用于 水文过程数值模拟所需的数据(通常是流域下垫面数据)基础上，采用 以泥石流形成物理过程中水土耦合特征变化趋势为依据的分析手段，在 降雨-泥石流形成的时间序列上捕捉降雨强度(I)、降雨持时(D)数据 点用于拟合I-D曲线方程(即确定阈值)。具体采用的水土耦合特征变化 趋势依据是在前期研究的基础上，以泥石流体的水土混合物容重值ρ的 动态变化为依据判断泥石流形成发生进程变化。当在时间序列上当前试 算的水土混合物容重值ρ达到预设的水土混合物容重控制点ρ′时，判断 当前点是拟合数据点。二、采用ZL 201210193426.4一种基于流域水土 耦合的泥石流预报方法公开的方法实现对水土混合物容重值ρ的动态测 量。具体是，通过对流域数据(包括地形、土地利用基础资料)进行栅 格化处理，建立流域栅格单元并获取每个栅格单元的本底数据；将流域 降雨空间分布数据与栅格单元本底数据作为输入数据，利用分布式水文 模型获取降雨条件下流域内每个栅格单元的土体含水量垂向分布数据与 每个栅格单元形成的径流量数据，进而计算确定降雨条件下区域内每个 栅格单元的边坡稳定性安全系数Fs、失稳固体物质体积总量W_S、区域径 流体积总量W_w，最后计算确定水土混合物容重值ρ。三、分布式水文模 型的选择是为通过数值模拟方法完成区域水文过程模拟，以期尽量准确 地获得用于水土混合物容重值ρ的动态测量所需的关键水文参数。泥石 流灾害大多发生在雨季，因而每次降雨过程中的土体内部水文参数的变 化特征关键所在。本发明根据前期研究分析降雨作用于下垫面的一般水 文过程为：土体初始为非饱和状态，但随着降雨入渗的持续作用，土体 内部的含水量持续增加，土体内部的基质吸力(ψ)骤降，进而会导致直 面土体失稳。由此，为在水文过程模拟中充分体现土体内部的土体含水 量以及基质吸力的变化规律，本发明选择联合一维非饱和Richard数学模 型(式1-1)与VanGenuchten方程(

式2)描述 降雨下渗的物理过程。

本发明进一步提供上述泥石流发生降雨I-D曲线阈值构建方法的几 种可独立实施的优化方案，具体包括：

优化方案一、Richard数学模型初始条件

初始条件：土体体积含水率θ确定为土体初始体积含水量θ_a，本发 明设计Richard数学模型初始条件由式1-2控制。上边界条件：利用超渗 产流机制控制入渗上边界，即当降水强度小于地表入渗能力(f_s)时，本 发明设计Richard数学模型上边界条件由式1-3控制。当降水强度超过地 表的入渗能力(f_s)后随即转化为径流流走，不考虑洼地积水的有压入渗， 此时地表饱和，式1-3控制的边界条件自此转化为θ＝θ_s(θ_s为某种土壤 类型对应的饱和体积含水率)。下边界条件：前期研究表明降雨入渗仅能 影响位于表层土以下40cm左右土层的土体含水量，无法与地下水位(埋 深数米)进行水分交换，故可以以此为约束条件。基于上述分析，本发 明设计的Richard数学模型下边界条件由式1-4控制。由此，优化方案中 Richard数学模型是以下方程组：

式中，θ—当前土体体积含水量，mm，

D(θ)＝K(θ)/(dθ/dψ)—非饱和土的扩散率，mm³/h，ψ依式2确定，

K(θ)—土体非饱和导水率mm/h，依式3、式4确定

R(t)—当前降水强度，mm/h，泥石流流域数据确定，

θ_a—土体初始体积含水量，％，依常规方法确定，

z—竖向坐标，以地表为原点，向下为正，mm，泥石流流域数据 确定，

L—下边界深度，m，泥石流流域数据确定，

m＝1/(1-n) 式4

式中，k₀—每层土在残余含水量条件下的非饱和入渗系数，泥石流 流域数据确定，

S_e—土体有效饱和度，％，泥石流流域数据确定，

m、n—曲线形状参数，n由泥石流流域确定，

k_s—表层土入渗系数，区域数据确定，

e—自然常数。

优化二方案、每个栅格单元的边坡稳定性安全系数Fs

造成浅层滑坡的主要因素包括降水入渗导致非饱和坡体内部基质吸 力的迅速降低，或者饱和土坡体内部孔隙水压力的提升两方面。具体考 虑到泥石流多发在雨季，在雨季来临之前，坡体一般都处于非饱和状态。 降雨入渗使得土体含水量的增加，进而造成土体的基质吸力降低，是导 致坡面土体失稳的主要诱因。在此，基于非饱和土的摩尔-库伦破坏准则 与无限边坡模型，本发明建立了用于判定泥石流流域内的两岸坡体稳定 性判定式(式5)，以优化技术方案。

式中，

—土体的内摩擦角，°，调查采样与测试试验确定，

λ—每个栅格的坡度，°，泥石流流域数据确定，

c—土体的粘结力，kPa，调查采样与测试试验确定，

ψ—土体基质吸力，kPa，式2确定

—与基质吸力相关的土体参数，取值内摩擦角

°，

γ_t—土体密度，g/cm³，调查采样与测试试验确定，

D_s—每个栅格失稳土层厚度，m，将栅格土层厚度分层代入式5试 算得到Fs<1时取对应土层厚度值。

γ_t—土体密度，g/cm³，调查采样与测试试验确定。

优化方案三：细分前期雨量Ar′及其对应I-D曲线阈值

根据现有研究成果确定用于水文模拟计算的前期雨量Ar′时，首先将 前期雨量Ar′设定为一个区间值，再将区间值划分为数个量级Ar′_i，量级 间相差5mm～15mm。然后将不同量级的前期雨量Ar′_i作为水文模型输 入条件，得到不同量级前期雨量Ar′_i对应的I-D曲线阈值，从而提高I-D 曲线阈值的灵敏性。

优化方案四：设置I-D曲线阈值的预警边界

现有技术存在技术原理上忽视I-D曲线阈值的约束条件，在实际预 警中忽视I-D曲线阈值的预警边界的缺陷。为此，本发明进一步将上述 泥石流发生降雨I-D曲线阈值构建方法的混合物容重控制点ρ′设置为一 个区间值中的多个点，每个点代表不同的控制容重，所有点组成泥石流 水土混合容重范围值。这样，在每一个前期雨量Ar′条件下，构建出一组 对应不同控制点ρ′的I-D曲线阈值，每个曲线阈值与相邻曲线阈值的间距 即是该曲线阈值的预警边界。在实际预警中，当构建所得最大控制点ρ′曲 线与最小控制点ρ′曲线阈值间距越大，表示泥石流发生的概率越大；反 之则泥石流发生概率低。这一优化能够丰富泥石流发生预警信息的内容。 尤其在预报降雨条件下实施泥石流预警时，丰富的预警信息具有更明显 的价值。

优化方案五：确定特定泥石流流域的水土混合物容重控制点ρ′、前 期雨量Ar′

本发明上述泥石流发生降雨I-D曲线阈值构建方法在一般情况下需 首先获取的泥石流流域数据主要是用于一维非饱和Richard数学模型与 Van Genuchten方程水文过程数值模拟所需的下垫面数据。但对于具备流 域内降雨数据监测数据与流域内泥石流活动观测数据的泥石流流域，可 以进一步采用优化条件，即在获取的流域数据中增加流域内降雨数据监 测数据、流域内泥石流活动观测数据，利用此二数据针对性分析本流域 泥石流及降雨条件的具体指标，确定本流域内典型泥石流的水土混合物 容重并将其记为水土耦合后的混合物容重控制点ρ′、确定本流域内典型 泥石流发生的前期雨量记为水文模拟计算的前期雨量Ar′，从而提高所构 建I-D曲线阈值对特定流域的敏感性。

基于上述泥石流发生降雨I-D曲线阈值构建方法，本发明同时提供 一种泥石流发生降雨预警方法，其技术方案是：

一种利用上述泥石流发生降雨I-D曲线阈值构建方法实现的流域泥 石流预警方法，其特征在于：在泥石流流域内布置雨量监测***，包括 分布在流域内的雨量传感器与控制中心，并在***中预置不同前期雨量 Ar′及对应的一组不同控制点ρ′的I-D曲线阈值；所述不同前期雨量Ar′及 对应的一组不同控制点ρ′的I-D曲线阈值由泥石流发生降雨I-D曲线阈值 构建方法确定；运算中心执行方案一和/或方案二：

方案一：根据天气预报降雨数据，当在计算前期雨量条件Ar₁下的降 雨强度-降雨持时组合位于Ar₁对应的一组I-D曲线阈值之间区域时，发出 泥石流预警信号；

方案二：当降雨发生后，***中心依据流域内实时降雨数据计算实 时前期雨量值Ar₂与实时降雨强度-降雨持时组合，并依据Ar₂选择不同的 前期雨量Ar′及对应的一组I-D曲线阈值，若实时降雨强度-降雨持时组合 位于I-D曲线阈值之间区域，发出泥石流预警信号。

与现有技术相比，本发明的有益效果是：(1)本发明提供的I-D曲 线阈值构建方法摆脱了现有技术在历史数据资料的基础上采用统计分析 方法拟合曲线的构思，采用根据泥石流形成发生原理，以泥石流体物理 特征变化过程为条件寻找数据点用于拟合I-D曲线方程的新构思，实现 了从测量泥石流本身特征动态指标出发，而非从数理统计手段出发构建 泥石流发生条件阈值的突破。由此，本发明方法无需流域内降雨数据监 测数据与流域内泥石流活动观测数据(此二数据对现有技术是必须)，只 需要泥石流流域的下垫面数据便可启动曲线构建过程，并得到适用于本 流域的I-D曲线阈值。因而，本发明方法能够有效解决无资料泥石流流 域采用I-D曲线阈值实施泥石流预警的问题。(2)本发明采用根据泥石 流形成发生原理，以泥石流体物理特征变化过程为条件、从测量泥石流 本身特征动态指标寻找数据点用于拟合I-D曲线方程的新构思体现了泥 石流发生学机理，方法更具科学原理。(3)采用Richard方程与Van Genuchten方程相结合的水文数值模拟计算，并进一步优化初始条件与设 计边界条件，使运算所得用于泥石流水土混合物容量指标动态分析的关键水文参数更符合泥石流多发季降雨条件与多发区土体条件特征。(4) 由于采用数值模拟计算，使得可以将影响泥石流形成发生的前期雨量Ar 细分为不同量级，并最终得到与其对应的I-D曲线阈值。利用该结果， 可以在流域泥石流降雨预警中针对不同前期雨量条件启用不同的I-D曲 线阈值，从而提高I-D曲线阈值的灵敏性与整体预警有效性。(5)本发 明构建I-D曲线阈值的方法是利用泥石流发生过程中水土偶合物密度为 标准，所得阈值对应着泥石流的不同密度值。由于泥石流的观测密度具 有边界值约束(1.2g/cm³～2.3g/cm³)，所以不同前期雨量Ar′条件下的泥 石流I-D曲线阈值均具有不同的预警边界范围。在泥石流预警中，增加 利用预警边界范围的限制，即在预警过程中同时监控降雨强度-降雨持时 组合点与I-D曲线阈值边界条件两项指标，可明显有利于降低泥石流的 误报率。

附图说明

图1是小江流域的支流地形图。

图2是实施例一技术路线示意图。

图3是实施例一拟合I-D曲线阈值图。

具体实施方式

下面结合附图，对本发明的优选实施例作进一步的描述。

实施例一

如图1～图3所示，用本发明方法构建某泥石流沟道泥石流发生降 雨I-D曲线阈值。

某泥石流沟道位于云南省，属于小江流域的支流(图1)。沟内新构 造运动强烈，地震活动频繁，岩体支离破碎，再加上强烈的风化作用， 斜坡表层具有数十米厚的松散固体物质，是泥石流爆发的重要固体物源。 流域内最大降水量出现在海拔2500m～3000m地带，年降水量1200mm， 该地带正是某沟泥石流形成区。该沟成为小江流域泥石流暴发最频繁的泥石流沟。用于本实施例的数据来源于流域内观测台站积累近60年的泥 石流野外观测数据。

图2是技术路线示意图(图中“流域”是“泥石流流域”的简称)。

步骤S1、获取区域数据

确定泥石流沟道范围，获取区域数据。区域数据包括用于一维非饱 和Richard数学模型与Van Genuchten方程水文过程数值模拟所需的区域 下垫面基本数据(包括坡度图、土地利用图、土层厚度分布图、区域土 壤类型分布图、各类土壤水文参数数据、各类土壤物理力学参数数据、 植被指数分布图等)、区域DEM栅格数据。将区域下垫面基本数据配置 在各DEM栅格中。

某沟道的数字高程模型(DEM)由1:10000地形图生成，网格精度 10m。沟内主要分布3种土壤类型包括：燥红土、红黄壤和砾石土，其 中砾石土广泛分布在蒋家沟内的滑坡与崩塌处，是泥石流发育的主要固 体物质来源，每种土壤的水文学参数由国家土壤数据库获取(表1)。沟 道土地利用类型图为2000年的数据，每种土地利用的相关如参数植被指 数(NDVI)来源于MODIS数据库，分辨率250m。用于水文分析的下 垫面基础数据需要利用ArcGIS中的重采样功能，获取与DEM网格精度 相一致的数据。

表1土壤类型及用于数值模拟计算的水文力学参数

根据现有技术确定水土混合物容重控制点ρ′在1.2g/cm³～2.3g/cm³， 将其划分为：1.2g/cm³、1.5g/cm³、2.0g/cm³、2.2g/cm³4个控制点，确定 泥石流发生前期雨量Ar′＝1mm～100mm，并将其划分为5mm、10mm、 20mm、30mm、40mm、50mm、60mm、75mm、85mm、100mm等10个量级。

区域内调查采样及测试试验确定相关土力学参数。将水文力学参数、 土力学参数依据土壤类型分布分配至每个DEM栅格单元中。包括(表1、 表2)：土体内摩擦角

土体的粘结力c、与基质吸力相关的土体参数

取值内摩擦角

土体密度γ_t、表层土入渗系数f_s、土体非饱和参数α。

表2土力学参数

以下步骤S2、步骤S3、步骤S4，将步骤S1划分的10个量级前期 雨量Ar′分别作为输入条件，启动水文过程数值模拟运算(t＝0)，沿时间 序列试算流域内当前激发的水土混合物容重ρ。计算步长为1h(I变化范 围：1mm/h-I_D＝1)。

步骤S2、水文过程数值模拟

在分别给定前期降雨Ar′下，联合一维非饱和Richard数学模型(式 1-1)、VanGenuchten方程(式2)、式3、式4完成水文数值模拟，Richard 数学模型的初始条件是式1-2、上边界条件是式1-3、下边界条件是1-4。

数值模拟计算获得用于步骤S3计算的关键水文参数。

步骤S3、判断DEM栅格单元失稳状态

依式5计算每一栅格的安全系数Fs。对于Fs<1的栅格，判断为失 稳状态。

以第100行中的第125列号栅格为例，土体内摩擦角

栅格 的坡度λ＝30.2°、土体的粘结力c＝34.5kPa、土体基质吸力ψ＝2.2kPa、 与基质吸力相关的土体参数

取值内摩擦角

土体密度γ_t＝1.8g/cm³、 栅格失稳土层厚度D_s＝1.8m。依式5计算有边坡稳定性安全系数Fs＝ 0.996。栅格失稳土层厚度D_s的计算是：将该栅格的土层厚度2m均分为 10层，每一层厚度即为0.2m。从第一层起逐层试算(即D_s依次取值0.2m、 0.4m、……2.0m代入式5试算)每一层的Fs，当第i＝9层出现Fs＜1， 确定第i层土以上为失稳土层厚度＝i×0.2m＝1.8m。

步骤S4、搜寻曲线拟合数据点

步骤S41、计算累计失稳土体总量V_s

依式6计算失稳栅格的累计失稳固体物质体积总量W_S。

式6中，W_S—区域内自t＝0起的累计失稳土体总量(m³)，T—计算 的总时间长度，自t＝0起至当前时刻(h)，N_s—区域内失稳栅格数量， 统计自t＝0起至当前时刻，A_i—每个失稳栅格的面积(m²)，D_s—每个栅 格失稳土层厚度(m)。

步骤S42、计算区域径流体积总量W_w

依水文计算数据所得每个栅格单元表层土的蓄水能力判断失稳栅格 产流状态，具体是当栅格单元表层土含水量θ临近表层土饱和含水量θ_s (θ→θ_s)时，判断栅格能够产流；

自t＝0起，依式7计算区域径流体积总量W_w。

式7中，W_w—流域径流体积总量(m³)，T—计算的总时间长度，自 t＝0起至当前时刻(h)，N_w—流域内能够产流的栅格数量，统计自t＝0 起至当前时刻，D_w—每个栅格产流深度(m)，由水文模型计算确定。

步骤S43、确定曲线拟合数据点

沿时间轴依式9计算水土混合物容重ρ(参考ZL 201210193426.4相 关内容)。当计算所得ρ＝混合物容重控制点ρ′时，确定该点为拟合数据 点，记录对应的降雨历时D、降雨强度I。

步骤S5、拟合I-D曲线方程

由所有拟合数据点拟合得到每一前期雨量Ar′量级条件下不同控制 点ρ′的一组I-D曲线方程，并在坐标系上绘制I-D曲线阈值。图3分别是 前期雨量为5mm(a)、10mm(b)、20mm(c)、30mm(d)、40mm(e)、 50mm(f)、60mm(g)、75mm(h)、85mm(i)、100mm(j)条件下的 I-D曲线阈值。以20mm(图3(c))条件下的I-D曲线阈值为例，在用 于流域泥石流预警时，该I-D曲线阈值适用的预警条件是：流域实时前 期雨量Ar₂达到20mm左右时(Ar₂→Ar′)、预警边界在图中◆点线型曲 线阈值(容重1.2g/cm³)与●点线型曲线阈值(容重2.2g/cm³)之间。当监测或预报降水的持时和强度组合介于这两条边界线之间则表明泥石流 可能会发生，需要发布泥石流预警信息。

Claims (9)

1.泥石流发生降雨I-D曲线阈值构建方法，其特征在于：

首先，获取泥石流流域数据，所述泥石流流域数据包括流域内用于一维非饱和Richard数学模型与Van Genuchten方程水文过程数值模拟所需的数据；确定水土耦合后的混合物容重控制点ρ′、确定水文模拟计算的前期雨量Ar′；

其次，在给定降雨强度I的条件下启动水文模拟过程，联合一维非饱和Richard数学模型、Van Genuchten方程进行水文过程数值模拟，确定用于计算流域边坡稳定性安全系数Fs的关键水文参数；

所述一维非饱和Richard数学模型采用如式1的方程组1：

式中，θ—当前土体体积含水量，％，

D(θ)＝K(θ)/(dθ/dψ)—非饱和土的扩散率，mm³/h，ψ依Van Genuchten方程确定，

K(θ)—土体非饱和导水率mm/h，依式3、式4确定

R(t)—当前降水强度，mm/h，泥石流流域数据确定，

θ_a—土体初始体积含水量，％，依常规方法确定，

z—竖向坐标，以地表为原点，向下为正，mm，泥石流流域数据确定，

L—下边界深度，m，泥石流流域数据确定，

m＝1/(1-n) (式4)

式中，k₀—每层土在残余含水量条件下的非饱和入渗系数，泥石流流域数据确定，

S_e—土体有效饱和度，％，泥石流流域数据确定，

m、n—曲线形状参数，n由泥石流流域确定，

k_s—表层土入渗系数，区域数据确定，

e—自然常数；

所述Van Genuchten方程如式2：

式中，S_e、θ—含义同前，

θ_s—某土体饱和含水量，％，依常规方法确定，

ψ—土体基质吸力，kPa，依常规方法确定，

θ_r—土体残余含水量，％，依常规方法确定，

α—土体非饱和参数，调查采样与测试试验确定；

再次，计算泥石流流域每一栅格的Fs值，判断每一栅格稳定性，沿时间序列试算当前流域内水土耦合后的混合物容重ρ，捕捉试算结果混合物容重ρ＝水土混合物容重控制点ρ′的点，确定该点为拟合数据点，记录拟合数据点对应的降雨历时D、降雨强度I；

最后，利用所有拟合数据点拟合得到流域在前期雨量Ar′条件下的I-D曲线I＝δD^β，在坐标系上绘制I-D曲线阈值。

2.根据权利要求1所述的方法，其特征在于：流域边坡稳定性安全系数Fs依式5计算确定：

式中，

—土体的内摩擦角，°，调查采样与测试试验确定，

λ—每个栅格的坡度，°，泥石流流域数据确定，

c—土体的粘结力，kPa，调查采样与测试试验确定，

ψ—土体基质吸力，kPa，Van Genuchten方程确定

—与基质吸力相关的土体参数，取值内摩擦角

γ_t—土体密度，g/cm³，调查采样与测试试验确定，

D_s—每个栅格失稳土层厚度，m，将栅格土层厚度分层代入式5试算得到Fs<1时取对应土层厚度值，

α—土体非饱和参数，调查采样与测试试验确定。

3.根据权利要求1或2所述的方法，其特征在于：所述混合物容重控制点ρ′是一个区间值，再在区间值内选取多个点，每个点代表不同的控制容重。

4.根据权利要求3所述的方法，其特征在于：所述前期雨量Ar′是一个区间值，再将区间值划分为数个量级Ar′_i，各量级间相差5mm～15mm。

5.根据权利要求4所述的方法，其特征在于：所述前期雨量Ar′是1mm～100mm，并划分为5mm、10mm、20mm、30mm、40mm、50mm、60mm、75mm、85mm、100mm共10个量级。

6.根据权利要求4所述的方法，其特征在于：所述混合物容重控制点ρ′范围是1.2g/cm³～2.3g/cm³。

7.根据权利要求6所述的方法，其特征在于：所述混合物容重控制点ρ′是1.2g/cm³或1.5g/cm³或2.0g/cm³或2.2g/cm³。

8.根据权利要求4所述的泥石流发生降雨I-D曲线阈值构建方法，其特征在于：所述用于一维非饱和Richard数学模型与Van Genuchten方程水文过程数值模拟所需的数据包括坡度图、土地利用图、土层厚度分布图、区域土壤类型分布图、各类土壤水文参数数据、各类土壤物理力学参数数据、植被指数分布图。

9.利用权利要求4～8任一所述的泥石流发生降雨I-D曲线阈值构建方法实现的流域泥石流预警方法，其特征在于：在泥石流流域内布置雨量监测***，包括分布在流域内的雨量传感器与控制中心，并在***中预置不同前期雨量Ar′及对应的一组不同控制点ρ′的I-D曲线阈值；所述不同前期雨量Ar′及对应的一组不同控制点ρ′的I-D曲线阈值由泥石流发生降雨I-D曲线阈值构建方法确定；运算中心执行方案一和/或方案二：

方案一：根据天气预报降雨信息，当在计算前期雨量条件Ar₁下的降雨强度-降雨持时组合位于Ar₁对应的一组I-D曲线阈值之间区域时，发出泥石流预警信号；

方案二：当降雨发生后，***中心依据流域内实时降雨数据计算实时前期雨量值Ar₂与实时降雨强度-降雨持时组合，并依据Ar₂选择不同的前期雨量Ar′及对应的一组I-D曲线阈值，若实时降雨强度-降雨持时组合位于I-D曲线阈值之间区域，发出泥石流预警信号。

Images