学术研究报告：考虑空间异质性的流域尺度滑坡水-力耦合模型
一、 作者、机构与发表信息
本文的主要作者为来自中国科学院地理科学与资源研究所资源与环境信息系统国家重点实验室的邓志超、兰恒星（通讯作者）、李郎平（通讯作者）、刘一宪和田乃曼。其中，邓志超、李郎平和刘一宪也来自中国科学院大学。兰恒星同时任职于长安大学地质工程与测绘学院。该研究发表于Journal of Hydrology期刊第667卷（2026年），文章编号134855，于2025年12月25日在线发表。
二、 研究背景与目标
该研究属于降雨诱发滑坡机理与预测模型领域。空间异质性是自然系统的基本特征，也是控制降雨诱发滑坡启动与分布的关键因素，它既体现在地表水平方向上（如土地覆盖类型变化），也体现在土壤剖面垂直方向上（如土壤结构、植被根系分布）。尽管基于水-力耦合模型（Hydro-Mechanical Coupling Model, HMCM）的物理模拟是当前最先进的滑坡分析方法，但现有模型大多仅部分考虑了空间异质性，通常只在水文模块或力学模块之一处理变异性。例如，一些集成的三维水文模型（如使用Richards方程求解器）能够较好地刻画水文过程的空间异质性，但其配套的稳定性模型往往将斜坡单元或计算网格内的力学参数均一化处理。反之，一些基于元胞自动机（Cellular Automaton, CA）理论、考虑土柱间力传递的稳定性模型能处理力学层面的水平异质性，但其依赖的水文模块往往采用简化的一维湿润锋法，无法捕捉水文层面的水平异质性，且通常忽略垂直方向上的力学参数变化。因此，缺乏一个能够同时捕捉水平和垂直异质性，从而更真实地再现滑坡启动耦合过程的综合建模框架。
基于此，本研究旨在提出并验证一个考虑空间异质性的流域尺度滑坡水-力耦合模型。其核心目标是：通过耦合一个支持三维水文参数异质性的集成水文模型（GEOTOP）与一个能够表征水平和垂直力学异质性的新型层状土柱稳定性模型（Layered Column Stability Model, LCSM），构建一个统一的空间离散化策略，实现对水文和力学过程空间异质性的综合考量，从而克服现有模型的局限，提高对降雨诱发滑坡（尤其是群发性滑坡）空间分布与破坏深度的预测精度。
三、 详细研究流程与方法
本研究包含方法构建、案例验证与结果分析三个核心流程，具体细节如下：
1. 水-力耦合模型（HMCM）的构建与耦合策略 * 水文模型（GEOTOP）：研究选用GEOTOP 3.0作为水文模块。该模型基于数字高程模型（DEM）网格进行空间离散化，地表每个二维网格作为遵循D8拓扑的径流单元，地下部分垂直划分为多层。其水文核心采用有限体积法求解三维Richards方程，模拟变饱和流，能够动态解析饱和与非饱和状态，捕捉侧向渗流和基质吸力的瞬态变化，从而有效考虑水文层面的三维空间异质性。 * 稳定性模型（LCSM）：研究团队自主开发了层状土柱稳定性模型（LCSM）。该模型基于GEOTOP的数据结构和元胞自动机理论，创建了“滑坡单元”（Landslide Cell, LC）对象。每个LC的土壤被分层，并能与GEOTOP输出的各层非饱和含水量、干土重等进行一一映射。稳定性计算采用改进的无限斜坡法，但关键创新在于：（a）考虑垂直异质性：计算每层土体的抗剪强度时，分别累加干土应力、水应力和植被应力，避免了参数均一化带来的误差。公式能够自动识别每层土壤的含水状态（非饱和、瞬态饱和、稳定饱和）并应用相应的有效应力公式，下方土层在计算自重应力时可使用浮容重。（b）考虑水平异质性：基于初步计算的每个LC安全系数（Fs）分布图，模型执行迭代的力传递。不稳定LC首先向上坡邻居请求抗滑力补充，若不足则向下坡邻居请求；同时，上坡不稳定LC会对下坡LC产生附加推力（降低其稳定性）。力的分配采用多轮加权算法，直至传递的累计力低于收敛阈值，从而模拟土柱间的力学相互作用。 * 模型耦合与工作流：两个模型采用基于DEM网格的统一空间离散化策略。GEOTOP首先运行，输出孔隙水压力和体积含水量的时间序列。在每个时间步后，调用LCSM：LCSM读取GEOTOP输出，为每个LC的各土层分配静态和动态参数；进行基于分层公式的初步稳定性计算；执行多轮力传播迭代，更新Fs直至收敛。耦合采用显式顺序策略，LCSM读取但不修改GEOTOP的状态变量，以确保数值稳定性。研究还对模型关键参数（如有效内摩擦角φ‘、坡度β、有效粘聚力c’、含水量θn等）进行了敏感性分析，确认了参数设置的合理性，并验证了模型在计算性能上适用于流域尺度模拟。
2. 案例研究：数据准备与模型校准 * 研究对象与区域：研究选取中国福建省武平县十方镇2022年5月发生的一次群发性降雨滑坡事件作为验证案例。研究区域为一个约22平方公里的流域。该区域为典型的低山丘陵区，降雨丰富，覆盖厚层第四纪残积砂质粘土，植被覆盖度高（以樟树、杉木、马尾松、毛竹等为主），滑坡多为浅层土质滑坡，滑动面基本位于残坡积土体内根系分布与无根系土层的界面附近。 * 数据输入：整合了多源数据构建高精度输入数据集，包括：12.5米分辨率的DEM、坡度、坡向、土地覆被数据；源自NOAA站点并经本地监测数据校正的逐小时降雨和温度数据；通过现场勘察和实验室试验获取的土壤水力学参数（饱和含水率θs、导水率Kv/Kh等）、力学参数（干密度γs、有效摩擦角φ‘、有效粘聚力c’）以及植被根系加固强度（τr）和根系深度（Rd）等。根系参数根据土地覆盖类型进行了沿深度的差异化赋值。 * 模型校准：首先，利用2024年研究区内三个站点采集的多层（30-150厘米）土壤含水量观测数据，对GEOTOP的关键土壤水力参数进行校准，优化后的模型在相应深度上表现出良好的模拟性能（平均R²=0.787）。然后，使用2022年5月1日至20日的区域平均表层土壤水分序列作为目标，调整模型以获取滑坡事件前真实的初始湿度场。校准后的模型成功重建了5月21日至28日土壤湿润过程前期的水文条件，为后续滑坡反演分析提供了可靠的初始场。
3. 滑坡事件反演与模型验证分析 * 滑坡空间分布与演变模式分析：运行耦合模型（记为Model I），反演了极端降雨期间及之后流域土壤水分和安全系数的时空演变。分析了不同土层（浅层0.6米、中层1.2米、深层2.1米）含水量以及最小Fs（即滑面位置）在降雨前、极端降雨期和降雨后六个时间切片上的空间分布规律。同时，选取了三种典型土地覆被类型（灌木/竹林、阔叶林、草地/农田）的土柱，深入分析了其土壤含水量和Fs沿深度和时间的变化。 * 考虑垂直异质性的对比验证：为了突出垂直异质性的重要性，设置了对比模型（Model II），该模型使用LCSM但不考虑深度相关的水力和力学参数（即参数均一化）。将Model I和Model II的模拟滑坡（基于Fs输出合并相邻不稳定单元生成滑坡多边形）与通过随机森林分类和目视解译建立的实地滑坡编录图进行对比验证。 * 验证指标：从多个维度评估模型性能：（a）确定性指标：比较模拟与观测的滑坡总面积、数量、空间分布（核密度分析）、在不同子流域、高程、坡度、土地覆被上的变化规律，以及滑坡面积累积分布函数。（b）概率性指标：将Fs场转换为破坏概率（Pf）场，绘制接收者操作特征（ROC）曲线，计算曲线下面积（AUC），并比较真阳性率（TPR）和假阳性率（FPR）。
四、 主要研究结果
1. 土壤水分与斜坡稳定性的时空演化规律 模拟结果显示，前期降雨已使流域浅层土壤湿度整体升高，灌丛和沟谷区域尤为明显。极端降雨开始后，浅层土壤迅速接近饱和。降雨结束后，森林区域因深根系植被的排水作用，浅层湿度快速下降，而灌丛区和沟谷区则保持高饱和状态。中层和深层土壤的水分响应存在延迟，且受植被根系引起的“瓶颈效应”（Bottleneck Effect）影响显著：在根系密集的土层下方，由于水力传导系数梯度大，水分易在根-土界面附近滞留，形成瞬态饱和区。 稳定性方面，Fs分布呈现西高东低的总体格局。前期降雨使部分高坡度区和灌丛区Fs下降，但无不稳定单元。极端降雨初期，稳定性变化不明显。随着降雨持续，森林区浅层土壤饱和度增加，不稳定区域逐渐扩大并向下坡方向延伸。在降雨峰值后约12小时，流域平均Fs达到最小，约11%的土柱变得不稳定，模型达到收敛。对典型土柱的深度剖面分析显示，在阔叶林和灌丛区，滑面主要出现在根系加固消失且含水量显著升高的深度（约1.2-1.7米），这与植被根系通过改变土壤水力参数梯度诱导“瓶颈效应”，从而控制滑动面形成的机制相符。
2. 模型验证结果与垂直异质性的重要性 * 滑坡空间分布：Model I模拟的滑坡总面积为7.995×10⁵ m²，数量为217个；Model II模拟的面积和数量分别为8.342×10⁵ m²和253个。实地观测值为面积6.954×10⁵ m²，数量189个。Model I在面积和数量上的高估程度（均为15%）均低于Model II（20%和34%）。核密度分析、与地理空间数据（子流域、高程、坡度、土地覆被）的关联分析以及滑坡面积累积分布函数均表明，Model I与观测值的一致性更好。在观测滑坡点位置上的Fs累积分布显示，Model I中约65%的滑坡落在Fs=1-1.5区间，而Model II中该比例低于50%。 * 概率性能评估：Model I的ROC曲线下面积（AUC）为0.7302，比Model II（0.6882）高出约6%。在代表性阈值下，Model I具有更高的真阳性率（TPR=0.266）和更低的假阳性率（FPR=0.076）。 * 破坏深度估计：对比两个模型模拟的土壤含水量沿深度变化及流域不稳定深度比例发现，Model II由于参数均一化，总入渗量更高，导致中、深层土壤含水量被高估，从而预测的破坏深度（在降雨后期可达1.7-2.1米）普遍深于Model I（主要集中在0.6-1.7米）。野外调查证实该区域滑坡厚度主要在0.5-1.5米，偶尔达2米，这与Model I的估计更为吻合。这直接证明了考虑沿深度（垂直）水-力学异质性，特别是与植被根系非均匀分布相关的异质性，能够显著提高滑坡破坏深度预测的准确性。
五、 研究结论与价值
本研究成功开发并验证了一个能够综合考虑水文和力学空间异质性的流域尺度滑坡水-力耦合模型（HMCM）。主要结论如下： 1. 模型有效性：所提出的耦合模型框架能够更真实地模拟降雨入渗过程、孔隙水压力演化以及斜坡稳定性的时空动态，成功再现了2022年福建武平群发滑坡的分布模式。 2. 垂直异质性的关键作用：考虑沿深度的水-力学异质性（尤其是植被根系引起的异质性）对于准确预测滑坡的空间分布和破坏深度至关重要。忽略这种异质性会导致模型高估不稳定区域的深度和范围。 3. 植被的“瓶颈效应”机制：研究阐明了植被根系通过改变土壤水力参数梯度，诱导水分在根-土界面附近滞留（即“瓶颈效应”），从而削弱下层无根系土层的抗剪强度并控制滑动面形成的物理机制。 4. 提供了平衡效率与精度的物理框架：该模型为模拟降雨诱发滑坡提供了一个在计算效率与预测精度之间取得平衡的物理框架。基于此框架，可以进一步评估空间异质性（特别是由植被引起的流域尺度地表异质性）对滑坡的影响。
该研究的科学价值在于推动了物理滑坡模型的发展，明确了在统一框架下同时处理多维度空间异质性的必要性和有效途径。其应用价值在于提高了对群发性降雨滑坡的预测能力，为滑坡风险评估和早期预警提供了更可靠的工具。
六、 研究亮点
七、 其他有价值内容
研究还讨论了模型的局限性及未来展望，包括：模型未显式解决地形变形和土壤应变过程，可能低估陡坡上因孔隙度突变引起的快速孔隙压力瞬变；在部分子流域出现的误报可能与数据支持和参数不确定性有关；未来可通过集成机器学习辅助校准、引入概率框架处理参数不确定性、结合近地表地球物理技术进行滑面间接验证等方法来增强模型的稳健性和适用性。模型模块化结构也使其有潜力适应更复杂的地质灾害环境（如冰川边坡失稳、地震扰动边坡等）。这些讨论为后续研究指明了方向。