基于深度积分法和不确定性分析的固废堆填土边坡失稳致灾风险定量评价：以深圳“12•20”光明滑坡为例

Quantitative Risk Assessment of Solid Waste Landfill Slope Instability-Induced Disasters Based on Depth-Integrated Method and Uncertainty Analysis:A Case Study of Shenzhen’s “12•20” Guangming Landslide

摘要

城市固废堆填形成的高陡边坡突发失稳常伴随大量建筑损毁和人员伤亡等灾害.为定量评估人工堆填滑坡的致灾风险，以深圳“12•20”光明滑坡为典型案例，采用基于深度积分法的Massflow软件构建动态数值模型，再现了滑坡失稳启动、高速运动及堆积的演进全过程.通过耦合滑坡物质冲击压力计算方法，定量评估了滑坡运动体对周边建筑物的冲击破坏效应，进而引入不确定性分析框架，将岩土体参数及滑动面特性视为随机变量，通过拉丁超立方抽样构建概率分析模型，揭示了多重因素影响下滑坡运动距离的统计特征，并建立了滑坡运动超越概率与建筑物冲击破坏效应间的关联模型.在此基础上，考虑承灾体易损性，对建筑物及其附近人员进行风险评价.研究结果表明：深度积分法有效捕捉了滑坡的动态演化特征，模拟所得滑裂面形态、运动距离（1 139 m）及堆积形态与实际监测数据高度吻合；滑坡运动物质对建筑物的冲击压力呈现先快速增大至峰值后逐渐衰减的规律，且峰值压力随距滑坡后缘距离增加而显著减小；通过基于参数不确定性的滑坡运动距离超越概率分析，实际堆积范围完全落于计算所得95%置信区间内，建筑物损毁区域与危险性区划结果高度匹配.本研究提出的耦合深度积分法与不确定性分析的承灾风险评价框架，为人工堆填滑坡的定量风险评估提供了新方法，对“无废城市”建设中的固废堆填安全保障具有重要理论价值与实践意义.

Abstract

Urban solid waste landfill-induced high-steep slope instability often leads to secondary disasters such as extensive building destruction and casualties. To quantitatively evaluate the disaster risk of artificial landfill landslides, this study takes Shenzhen’s “12•20” Guangming landslide as a case study, employing the depth-integrated method-based Massflow software to construct a dynamic numerical model that reproduces the entire evolution process of landslide initiation, high-speed movement, and deposition. By coupling with impact pressure calculation methods for landslide masses, the destructive impact effects on surrounding buildings were quantitatively evaluated. An uncertainty analysis framework was introduced, treating geotechnical parameters and sliding surface properties as random variables. A probabilistic analysis model was developed through Latin hypercube sampling, revealing the statistical characteristics of landslide travel distance under multifactorial influences and establishing a correlation model between landslide motion exceedance probability and building impact damage. Building vulnerability was further considered when conducting risk assessments for structures and nearby personnel. The study demonstrates that the depth-integrated method effectively captures the dynamic evolution of landslides, with simulated failure surface morphology, travel distance (1 139 m), and deposition pattern highly consistent with field monitoring data; impact pressure from landslide masses on buildings exhibits a rapid rise to peak values followed by gradual attenuation, with peak pressure decreasing significantly with distance from the landslide source; probabilistic analysis of landslide travel distance based on parameter uncertainty shows the actual deposition zone entirely falls within the 95% confidence interval, and building damage areas align with hazard zoning results. The proposed risk assessment framework integrating depth-integrated modeling and uncertainty analysis provides a novel methodology for quantitative evaluation of engineered landfill landslides, offering significant theoretical and practical implications for solid waste landfill safety management in the context of “zero-waste city” development.

Graphical abstract

关键词

风险定量评价 / 固废堆土 / 光明滑坡 / 深度积分法 / 不确定性分析 / 工程地质.

Key words

quantitative risk assessment / solid waste stockpile / Guangming landslide / depth-averaged method / probabilistic uncertainty analysis / engineering geology

0 引言

随着我国城市化进程加速，城市固废堆填土风险管控已成为保障生态安全与推进“无废城市”建设的迫切需求（Yang et al.，2016）.《中共中央国务院关于全面推进美丽中国建设的意见》指出，我国将于2027年实现60%的“无废城市”建设覆盖率，并计划于2035年达成全国范围的“无废城市”建设目标（柴立元等，2025）.过去经济高速发展时期遗留了大量城市固体废物堆积问题，据2024年中国生态环境统计年报，2023年全国一般工业固体废物产生量达42.7亿 t，处置量8.3亿 t，综合利用率不足61%（约25.8亿 t），相关固体废物堆积区生态环境保护与修复工作面临严峻挑战.除了影响生态环境，固废堆填土失稳灾害频发.2023年，重庆市渝北区双井村弃土场因固废堆填土失稳发生灾害，造成3人死亡、2人受伤及直接经济损失约957.274 2万元（重庆市渝北区应急管理局，2023）；2024年，江西省萍乡市湘东区杨家田地块截至违规堆存540万 t工程渣土，形成高悬头顶的“渣山”，持续威胁周边建筑安全（中华人民共和国生态环境部，2024）；新疆维吾尔自治区克拉玛依市建筑垃圾填埋场违规堆填固体废物达1 700万 t，进一步加剧了区域生态环境压力.这些典型案例表明，固废堆积物失稳不仅直接威胁居民生命财产安全，更对周边建筑物、生态环境产生不可避免的次生灾害链式影响.因此，亟需开展固废堆积物失稳风险评估，构建固废堆填土次生灾害链风险评估体系.

在众多的固废堆填失稳案例中，2015年深圳光明滑坡作为城市人工堆填体失稳，其独特的致灾机理、高速远程运动特征及显著的社会危害性，持续引发学界高度关注（钟兴荣，2024）.该滑坡体积达270万 m³，滑动距离超1 km，造成77人遇难及重大经济损失，其复杂演化过程为灾害评估提出了重大挑战（Yin et al.，2016）.针对此类灾害的精细化评估需求，数值模型凭借其低成本、高效率及可重复性等优势，成为滑坡灾害研究的主流工具（唐辉明等，2022）.针对深圳"12•20"光明滑坡等典型案例，现有方法可归纳为3类：有限元/差分法、无网格大变形方法以及深度积分法.有限元/差分法（如Geostudio中的SEEP/W模块）在解析降雨渗流与地下水耦合作用方面具有优势（Yin et al.，2016；Wang et al.，2020；Gao et al.，2021），但受网格畸变限制，难以准确模拟高速远程滑坡的动力学特征，且无法预测堆积范围、运动速度等关键动力学参数.而无网格大变形方法如光滑粒子流体动力学方法（Gao et al.，2019；Peng et al.，2019；Li et al.，2020）及物质点法（孙玉进和宋二祥，2018；王升等，2022）尽管能揭示基底摩擦效应、孔隙水压力响应、应力路径演化及流变机制对滑坡破坏模式与堆积形态的控制作用，但该类方法建模复杂度高，计算效率低下，需耗费大量算力资源，难以应用于大规模固废堆填场景.深度积分法通过三维计算的降维处理（简化为二维问题），在保证对滑体内部非连续界面动态演化过程精准追踪的同时，大幅度提高了计算效率.Iverson and Ouyang（2015）创新性地将库仑破坏准则与侵蚀模型耦合，并据此开发Massflow数值模拟软件.以深圳光明滑坡为典型案例，该软件成功揭示高孔隙水压力对滑体流动性的显著增强机制.然而，传统数值模型存在固有局限：其将岩土体参数及滑动面特性简化为确定性变量，忽视了滑坡演化过程中水文地质条件与工程地质条件的时空变异性（蒋水华等，2024）.事实上，滑坡运动受地下水渗流、降雨入渗、岩土体非均质性等多源随机因素的共同影响，且现场监测数据因测量误差、采样偏差及人为操作等因素存在显著不确定性（孙小平等，2021）.这种理想化假设可能导致模型预测结果与实际情况产生较大偏差，削弱风险评估的可靠性.

为突破传统模型的确定性假设局限，可靠度理论为固废堆填滑坡研究提供了新的量化框架.针对深圳光明滑坡案例，Sun et al.（2021）将孔隙水压力系数与内摩擦角视为随机变量，通过构建极限状态函数评估滑坡运动超越概率，实现了危险性区划的定量表达.Zhang et al.（2021）进一步考虑土壤性质的空间变异性，通过蒙特卡洛法抽样绘制出失效概率-伤亡人数曲线，为人员伤亡风险评估提供了可量化的决策依据.这些研究虽在评估体系构建上取得突破，但仍存在两方面局限：其一，研究重点多集中于滑坡运动学特征的单体模拟，针对滑坡高速远程运动与建筑物冲击效应的耦合机制研究相对薄弱；其二，现有建筑物冲击力计算往往显式表征结构细节（孙玉进和宋二祥，2018；王升等，2022），这种高精度建模方式虽能捕捉局部应力分布，却因计算效率低下（单次模拟需数小时），难以满足大规模场景的快速评估需求（Wang et al.，2023；Sun et al.，2025）.更重要的是，当前研究尚未建立滑坡运动距离超越概率与冲击破坏效应之间的关联模型，对承灾体易损性的考量停留在单一建筑物或人口要素层面，缺乏对固废堆填滑坡过程中建筑物易损性与人口易损性的系统整合，这在一定程度上制约了风险评估体系的完整性.

基于此，本研究利用Massflow软件构建深圳光明滑坡的动态数值模型，对滑坡全运动过程进行模拟分析，尤其对于建筑物周边运动物质的速度和堆积厚度进行了重点分析；在此基础上，将动态模型与拉丁超立方抽样进行有机结合，建立多重因素影响下的滑坡空间影响概率分析模型并进行危险性区划；进一步构建滑坡运动距离超越概率与冲击破坏效应耦合的风险评估体系，同步纳入建筑物易损性和人口易损性指标，为人工堆填滑坡的预测预警及防治提供理论支撑.

1 深圳光明滑坡基本概况

2015年12月20日，深圳市光明新区发生特大城市固废堆填体失稳滑坡灾害.该滑坡主滑方向N20°W，最大运动距离达1 100 m（垂直落差113 m）.如图1a所示，平面形态呈现“漏斗-扇状”扩散特征（最大宽度630 m，最小宽度150 m），具有典型的高速远程滑坡动力学特性（李坤等，2022）.地质剖面揭示（图1b）滑坡体呈双层结构特征：下部为花岗岩基岩构成滑床，上部为人工堆填体（厚度3~16 m，均值6~7 m），自上而下由建筑渣土层与采石场弃土层组成，堆积表面呈6.1°缓倾角（许强等，2016）.此地质模型为后续数值模拟提供了基础边界条件.

灾害发生于原红坳采石场遗址改建的余泥渣土受纳场，场地设计库容4.0×10 m³（封场标高95 m），但事故时实际堆填量已达5.83×10 m³（超容45.8%），堆体后缘标高突破设计值达160 m.场地特殊的地形条件与严重超量堆载形成耦合致灾效应：采石场遗留的基岩凹坑构成天然滑床，北侧基岩凸起形成侧向临空面，高势能堆体在重力作用下沿软弱滑动带启动（图1c）.未规范建设的导排水系统导致长期饱和渗流，进一步弱化了堆填体底部抗剪强度.基底地形约束使滑体运动过程中呈现扇形扩散，显著放大了灾害链式效应.据统计，此次灾害形成泥流覆盖面积约38×10 m²，滑出渣土量约2.7×10 m³，共造成73人死亡、4人失踪、17人受伤，并导致33栋建构筑物损毁（含24栋厂房、3栋宿舍及6栋民宅；Yin et al.，2016）.

2 基于连续深度积分法的确定性分析动态数值模型

本研究基于Ouyang et al.（2017）验证的基础上，采用流变连续介质假设构建动态模型，利用深度积分法建立不可压缩流体运动方程，同步耦合孔隙水压力效应与冲击荷载计算，为建筑物易损性分析提供了精确的动力学参数支撑.

2.1　确定性分析模型建立

2.1.1　控制方程

控制方程组由连续性方程（式（1））和动量方程（式（2~3））构成：

\frac{\partial (h)}{\partial t}+\frac{\partial (hu)}{\partial x}+\frac{\partial (hv)}{\partial y}=\mathrm{0},（1）

\frac{\partial (hu)}{\partial t}+\frac{\partial \left(h{u}^{\mathrm{2}}+\frac{k_{a/p}g{h}^{\mathrm{2}}}{\mathrm{2}}\right)}{\partial x}+\frac{\partial (huv)}{\partial y}=k_{a/p}gh\frac{\partial z_b}{\partial x}-\frac{{\left({\tau }_{zx}\right)}_b}{\rho }，（2）

\frac{\partial (hv)}{\partial t}+\frac{\partial (huv)}{\partial x}+\frac{\partial \left(h{v}^{\mathrm{2}}+\frac{k_{a/p}g{h}^{\mathrm{2}}}{\mathrm{2}}\right)}{\partial y}=k_{a/p}gh\frac{\partial z_b}{\partial y}-\frac{{\left({\tau }_{zy}\right)}_b}{\rho }，（3）

式（1）~（3）中，h为流动体厚度，u、v分别为x、y方向深度平均流速；g为重力加速度，\rho为流体密度；z_b为基底高程，{\tau }_{zx}、{\tau }_{zy}为基底剪切应力分量；k_{a/p}为侧向土压力系数，通过朗肯土压力理论（Christen et al.，2010）确定，如式（4）和（5）：

k_a=\frac{\mathrm{1}-\mathrm{s}\mathrm{i}\mathrm{n}\phi }{\mathrm{1}+\mathrm{s}\mathrm{i}\mathrm{n}\phi }(主动状态),（4）

k_p=\frac{\mathrm{1}+\mathrm{s}\mathrm{i}\mathrm{n}\phi }{\mathrm{1}-\mathrm{s}\mathrm{i}\mathrm{n}\phi }(被动状态),（5）

2.1.2　基底阻力模型

考虑孔隙水压力效应，采用改进的Coulomb摩擦模型（Xing et al.，2014），其公

式如式（6）所示：

\tau =\sigma \left(\mathrm{1}-r_u\right)\mathrm{t}\mathrm{a}\mathrm{n}\phi,（6）

式（6）中，r_u=\frac{u}{{\sigma }_n}为孔隙水压力系数，反映液化程度（0≤r_u≤1）；u为孔隙水压力，{\sigma }_n为基底有效正应力，\phi为动力摩擦角.

2.1.3　冲击荷载计算

采用Armanini冲击压力模型计算滑坡运动过程中单位深度对于建筑物的冲击力（Armanini，1997），将其考虑为动压力和静压力的耦合作用：

P_a=P_d+P_s，（7）

P_d=\rho v^{\mathrm{2}},（8）

P_s=\mathrm{0.5}\rho gh,（9）

式（7）~（9）中：P_{\mathrm{d}}为动压力； P_{\mathrm{s}}为平均静止压力； \rho为滑坡运动物质的密度； v为冲击面平均流速； h为流动体到达某位置的厚度； g为重力加速度.

总冲击荷载通过对于冲击力沿深度方向长度L积分得到：

F_a={\int }_{\mathrm{0}}^L\left(P_d+P_s\right)\mathrm{d}A=\rho h{v}^{\mathrm{2}}L+\mathrm{0.5}\rho g{h}^{\mathrm{2}}L .（10）

2.2　模型参数分析

首先对研究区地形进行简化，并构建二维结构化网格（网格单元尺寸设定为1 m，共计1 933个网格单元）.根据现场影像资料，将受影响的7栋建筑物按空间位置从东至西依次编号为A~G，如图2所示.模型关键参数包含4个物理力学指标：粘聚力c、内摩擦角\phi、土体密度\rho及孔隙水压力系数r_u.材料参数取值参照（Yin et al.，2016）的环剪试验成果，土体本构采用理想弹塑性模型，基底有效粘聚力c=\mathrm{0}\mathrm{k}\mathrm{P}\mathrm{a}，土体密度\rho =\mathrm{1}\mathrm{750}\mathrm{k}\mathrm{g}/\mathrm{c}{\mathrm{m}}^{\mathrm{3}}，内摩擦角\phi =\mathrm{24}°（峰值强度与残余强度参数统一取值）.孔隙水压力系数取r_u=\mathrm{0.75}（高杨，2018）以表征孔隙水压力对土体抗剪强度的折减效应.

控制方程的数值求解采用Godunov型有限体积法（Ouyang et al.，2015），确保计算精度与稳定性.模拟总时长设定为60 s，采用自适应变时间步长算法，初始时间步长设为0.05 s.计算过程中通过CFL数（courant⁃friedrichs⁃lewy数）控制时间步长，该无量纲参数用于平衡计算稳定性与收敛效率，根据文献（Ouyang et al.，2013）所知，取值为0.25.单次滑坡全过程模拟耗时27.3 s（计算机配置为Intel（R）Core（TM）i7⁃8750H CPU @ 2.20GHz 2.21 GHz）.

2.3　确定性分析结果

2.3.1　滑坡运动及堆积特征

深圳光明滑坡动态演化过程的模拟结果如图3a所示，其运动特征可划分为3个典型阶段：初始破坏阶段（0~13 s）：剪切破坏首先出现在坡趾位置，并沿潜在滑动面向深部扩展，至第13 s形成贯通性滑裂面，标志着滑坡正式启动；加速运动阶段（13~27 s）：贯通滑体在重力作用下开始整体加速，速度-时间曲线（图3b）显示该阶段速度增长速率达到峰值.至27 s时，前缘速度达18.6 m/s；稳定堆积阶段（27~60 s）：滑体后缘呈现阶梯状渐进破坏模式，主体部分经历高速远程运动后动能逐渐耗散.至第50 s时，运动距离达到1 139 m并保持相对稳定，最终堆积形态（图3c）与物质点法模拟结果（王升等，2022）及实际堆积特征均表现出良好一致性.

模型验证结果如表2所示，模拟的堆积形态与实际滑坡堆积特征高度一致.模拟最大堆积厚度为27.9 m（实际值21.3 m），相对误差12%；运动距离预测值1 139 m（实际值1 100 m），偏差仅3.5%，显著优于物质点法模型（1 260 m，偏差14.5%）.这进一步印证了深度积分法在捕捉滑坡大变形运动过程及堆积特征方面的独特优势，其高效处理复杂边界条件和非连续界面演化的能力为精准预测提供了重要支撑.但值得注意的是，模拟结果显示平均堆积厚度为8.4 m，较实际值12.9 m偏小约35%；堆积坡角3.99°亦小于实测值6.1°，差异达34%.经分析，这种偏差主要源于模型维度简化——将三维运动过程降维至二维平面，导致未能充分表征滑体物质的横向扩展行为，后续研究可通过引入横向扩散修正系数优化模拟精度.

2.3.2　滑坡对建筑物的影响

(1)冲击速度与堆积厚度.针对研究区建筑物进行冲击速度计算分析，结果如图4所示.从图4a显示，冲击速度呈现典型的三阶段演变特征：初始阶段速度快速攀升至峰值，随后经历短暂衰减期，继而再次增长形成次高峰，最终逐渐衰减至零值.此外，各监测点堆积厚度均经历快速增长后趋于稳定（图4b）.值得注意的是，建筑物A位置在30 s时刻出现堆积厚度异常波动，表现为短暂下降后于45 s左右再次回升.经分析，该现象可能与滑坡物质流经建筑物A与B之间陡坎时发生速度突变有关，前后段速度梯度导致运动物质产生回弹效应.

进一步对建筑物首次受冲击时间及最大冲击速度的时空分布规律展开研究.分析表明，建筑物受冲击时间随距滑坡后缘距离增加而滞后，但最大冲击速度并未呈现线性递减规律.如表2所示，建筑物B所承受的最大冲击速度反而高于建筑物A，这可能与运动物质在陡坎处产生速度放大效应及能量集中释放有关，与贺旭荣等（2024）基于物理模型试验观察到的能量耗散规律一致.

(2)冲击压力与冲击荷载.运动物质冲击压力计算结果如图5a所示.分析表明，滑坡运动产生的冲击压力均呈现典型的3阶段演变特征：初始阶段压力快速攀升至峰值，随后逐渐衰减至稳定水平.其中A位置记录的峰值冲击压力显著高于其他监测点，这与该位置受冲击时刻恰为滑坡运动速度峰值期的动力学机制密切相关.冲击压力时程曲线进一步揭示，滑坡运动物质对路径上建筑物的冲击压力以动压为主.随着滑坡运动停止，动压力逐渐消散，冲击压力最终退化为静压作用.

峰值冲击压力值随计算点与滑坡体距离的增加呈现显著衰减规律.具体而言，距滑坡后缘滑裂面最近的A点峰值冲击压力达475.8 kPa，而最远点E（距983 m）则降至53.2 kPa.根据Zanchetta et al.（2004）提出的框架结构建筑物易损性准则，当冲击压力超过100 kPa时将造成严重结构损毁.计算结果显示，A~D建筑物所受峰值冲击压力均超过此临界值，这为深圳光明滑坡中建筑物大面积损毁及淤埋现象提供了力学机制解释.

进一步分析冲击荷载时程特征（图5b），滑坡运动物质在建筑物位置产生的冲击荷载均表现为初期快速增长后趋于稳定的演变模式.表3数据显示，建筑物处最大冲击荷载约为稳定值的1.0~1.5倍，这验证了滑坡冲击过程中动荷载向静荷载转化的能量衰减规律，同时也反映了运动物质与建筑物相互作用的动力学特征.研究成果对于定量评估滑坡冲击荷载动态特性及其致灾效应具有重要理论价值.

3 滑坡超越概率驱动下的承灾体风险评价模型

3.1　风险评价模型建立

在确定性数值模型的基础上，引入可靠度理论框架（孙小平等，2021），将影响滑坡运动的关键物理力学参数（摩擦角\phi、密度\rho及孔隙水压力系数r_u）视为服从特定概率分布的随机变量，以量化参数不确定性对滑坡运动特征的影响.基于已构建的动态数值模型，采用拉丁超立方抽样方法（LHS）实现参数空间的高效抽样，并通过概率分析方法计算滑坡运动距离的超越概率（即滑坡运动距离L超过临界阈值L_{\mathrm{c}\mathrm{r}\mathrm{i}\mathrm{t}}的概率P\left(L>L_{\mathrm{c}\mathrm{r}\mathrm{i}\mathrm{t}}\right)），据此开展滑坡危险性区划.进一步，将滑坡运动学特征、承灾体易损性与应急响应能力进行耦合，构建致灾风险评价框架.技术实现路径如图6所示，为滑坡灾害风险的系统性评估提供理论支撑.

针对承灾体风险定量评估需求，构建多维风险指标体系.对于建筑物i，其风险值可定义为（Galasso et al.，2021）：

R_b\left(L_{crit,i}\right)=P\left(L>L_{crit,i}\right)·D_i·C_i,（11）

式（11）中，P\left(L>L_{crit,i}\right)为滑坡运动距离超越建筑物临界位置的超越概率；C_i 为资产价值系数，综合考虑建筑物类型（砖混／框架结构）、使用功能及区域经济发展水平3类指标确定权重；D_i 为建筑物易损性系数，通过冲击效应量化模型计算：

D_i={\left(\frac{F_{a,i}}{F_{crit,i}}\right)}^{\gamma },（12）

式（12）中，F_{a,i} 为基于概率抽样的建筑物平均冲击力，F_{crit,i} 为建筑物抗冲击阈值，\gamma 为损伤累积效应系数（取值范围\mathrm{1.5}~\mathrm{2.0}）.针对 F_{crit,i}的建模，考虑滑坡冲击荷载对于建筑物的单侧瞬时强冲击特性（冲击方向假设垂直于坡体），当建筑物长深比 \frac{L}{B}>\mathrm{3} 且临空面自由长度超过5 m时，可等效为悬臂梁模型.在定向荷载作用下，建筑物的挠曲变形可表示为：

\omega =\frac{F{L}^{\mathrm{3}}}{\mathrm{3}EI},（13）

式（13）中，L 为建筑物特征高度，EI 为等效抗弯刚度.根据临界安全准则，当挠曲度 \omega \le {\omega }_{\mathrm{m}\mathrm{a}\mathrm{x}} 时结构保持安全，由此可计算建筑物最大可承受的冲击力（严屹然等，2023）：

F_{\mathrm{c}\mathrm{r}\mathrm{i}\mathrm{t}}=\frac{\mathrm{3}EI{\omega }_{\mathrm{m}\mathrm{a}\mathrm{x}}}{L^{\mathrm{3}}}\cdot \frac{\mathrm{1}}{{\beta }_d}，（14）

式（14）中，{\beta }_d 为冲击荷载放大系数（取值\mathrm{1.8}~\mathrm{2.5}），用于修正脉冲荷载的动力效应.

针对人口脆弱性评估，参考Zhang et al.（2021）提出时空耦合风险模型：

R_{p\left(L_{crit,i}\right)}=P\left(L>L_{crit,i}\right)·V_i·E_i,（15）

式（15）中，V_i为人口暴露度系数，通过建筑物人口密度与空间分布特征确定；E_i为应急响应效率系数，定义为：E_i=\frac{{\mathrm{10}}^{\mathrm{3}}}{t_{resp,i}·{\eta }_{evac,i}}{\mathrm{e}}^{-k\bullet t_{\mathrm{d}\mathrm{a}\mathrm{t}\mathrm{a}\mathrm{y}}}；t_{resp,i}为预警响应时间，取抽样得到的建筑物首次受冲击的平均时间；{\eta }_{evac,i}为撤离效率系数，通过交通网络分析与历史撤离数据综合评估，t_{\mathrm{d}\mathrm{a}\mathrm{t}\mathrm{a}\mathrm{y}}为决策延误时间，k为时间敏感系数.

3.2　模型参数选取

以2.2节中提供的岩土体参数为基准值，综合现场调查（高杨，2018；孙小平等，2021；Yin et al.，2016）与《地质灾害防治工程勘察规范》，对参数变异性进行界定，如表4所示.

3.3　考虑承灾体易损性的风险评价

3.3.1　滑坡致灾危险性评价

采用LHS方法生成10 000个随机样本，并将各个样本代入确定性分析模型中逐一计算滑坡运动距离.通过收敛性分析发现，当抽样次数达到1 000次时，统计结果的均值和方差趋于稳定.对10 000次抽样结果进行统计分析，绘制概率密度分布图，其中正态分布检验显示经验累积分布函数（CDF）曲线与标准正态分布CDF存在显著偏离，Kolmogorov⁃Smirnov检验结果为h=1（p<0.05），明确拒绝数据服从正态分布的假设（如图7a所示）.进一步采用广义极值分布（GEV）进行拟合分析，参数估计结果为：形状参数k=\mathrm{0.085}(k>\mathrm{0}表明数据服从Frechet型极值分布，具有厚尾特征）、尺度参数\sigma =\mathrm{150.83}（反映数据离散程度较高）、位置参数\mu =\mathrm{772.67}（略低于样本均值\mathrm{900.175}\mathrm{m}，可能受极端值影响）.拟合优度检验表明，数据总体服从GEV分布（如图7b所示），且该分布能更好表征滑坡距离的极端值特征.统计结果显示：滑坡运动距离的平均值为900.175 m，最大值1 400 m，最小值509 m，95%置信区间为［895.72 m，902.74 m］.尽管实际观测值1 100 m超出了均值置信区间，但完全落在基于广义极值分布（GEV）构建的95%预测区间（509~1 400 m）内.由于数据服从具有厚尾特征的Frechet型极值分布（形状参数k=0.085>0），极端值出现概率显著高于正态分布假设，因此该偏差属于统计规律允许的正常现象，进一步验证了计算结果的合理性.

进一步，基于抽样结果，以空间任意点与滑坡体后缘之间的距离作为动态阈值，构建多个极限状态函数.利用蒙特卡洛模拟方法，计算了滑坡运动距离超越不同阈值的概率，并绘制了运动距离-超越概率响应曲线（如图8a所示）.分析表明，滑坡运动距离的超越概率随阈值距离的增加呈现非线性衰减规律，这定量揭示了空间位置危险性与滑坡体距离的负相关关系.参考Lacasse and Nadim（2011）的灾害概率分类框架，以滑坡运动距离超越概率50%、10%、1%和0.1%为临界阈值，将深圳光明滑坡的运动危险性划分为极高、高、中、低、极低5个等级危险区（如图8a所示）.

区划结果显示，滑坡实际堆积范围完全覆盖高至极高危险区.分析发现，建筑物A~B位于极高危险区，建筑物C~G位于高危险区，表明基于超越概率的危险性区划与实际情况高度吻合，可为滑坡风险评估与管理提供科学依据.各建筑物位置处的超越概率及其随抽样次数的演化曲线如图8b所示，通过与传统蒙特卡洛方法（MCS）对比，发现LHS方法展现出更优的收敛特性，这一结论与孙小平等（2021）在大堡子滑坡中对于LHS和MCS方法的对比结论基本一致，验证了所提危险性评估方法的有效性.

3.3.2　考虑承灾体易损性的风险评价

基于滑坡运动距离的超越概率分析，结合承灾体易损性与区域经济发展水平，对研究区建筑物及其周边人口开展风险评价.根据《深圳市建筑设计规则（2015）》及事故调查报告，对建筑物类型与结构特征进行如下假设（莫婉媚和方伟华，2016）：建筑物A、B为框架结构工业厂房（层高5 m），C、D为砖混结构仓储设施（层高6 m），E、F、G为低层框架结构宿舍楼（层高3 m）.资产价值参考曾鹏等（2022）对于海口市台风次生灾害损失评估方式，评价系数依据《深圳市建筑造价经济指标》及光明区2015年商品房均价（2.8 万元/m²）综合测算，引入结构类型修正系数（框架1.2，砖混0.8）与功能修正系数（工业0.9，民用1.1）.通过改进的悬臂梁模型计算建筑物抗冲击阈值，其中临界挠度比取1/150，冲击荷载放大系数{\beta }_d取2.0，参数计算说明见表5.

风险计算结果显示（表6），砖混结构建筑物（A、B、E）因抗冲击阈值较低（0.25~0.36 MN/m），在相同冲击力作用下易损性系数D_i显著高于框架结构（C、D、F、G）.此外，框架结构仓储设施C、D虽位于高危险区（917~992 m），但通过变形协调机制使实际冲击力（0.7~0.8 MN/m）低于抗冲击阈值（0.47 MN/m），最终风险值（0.56~0.63）显著低于砖混建筑A、B.

在人口易损性评估中，构建时空耦合风险模型，具体参数取值见表7.结果显示（如表8），A、B区域人口伤亡风险分别达到18.78和11.77，远超其他区域.这主要是由于这些区域的人口密度较高，且预警响应时间相对较短，导致应急响应效率系数偏低.尽管E区域的预警时间较长，但由于其人口密度极高，仍面临较高的风险.

为更直观展示风险评价结果，通过风险分布直方图（图9）进一步解析风险构成：如图9a所示，砖混建筑A、B、E因冲击力集中（F_{a,i}>\mathrm{2.0}\mathrm{M}\mathrm{N}/\mathrm{m}）导致D_i>\mathrm{1.3}，而框架结构通过变形协调机制显著降低冲击力（），其D_i值控制在0.3以下.人口风险方面，A、B区域因预警响应时间较短（15~19 min）导致应急响应效率系数E_i偏低，而E区域尽管预警时间较长（39 min），但高密度人口暴露（20人/m²）仍使其风险值达到4.52.风险空间分布呈现出明显的“双峰”特征（图9d）：第一峰值出现在A、B区域（R_p=\mathrm{18.78}-\mathrm{11.77}），由高易损性工业建筑与短预警时间耦合导致砖混结构自身抗冲击能力较弱（F_{crit}=\mathrm{0.25}-\mathrm{0.36}\mathrm{M}\mathrm{N}/\mathrm{m}），叠加短预警时间（t_{resp}=\mathrm{13.9}-\mathrm{17.0}\mathrm{s}）导致应急措施无法及时展开，从而加剧人口伤亡风险；第二峰值出现在E区域（R_p=\mathrm{4.52}），主要受夜间高密度人口暴露（\mathrm{2.0}\mathrm{人}/{\mathrm{m}}^{\mathrm{2}}）影响——尽管该区域预警时间较长（t_{resp}=\mathrm{32.8}\mathrm{s}），但高暴露系数（V_i=\mathrm{0.80}）直接放大风险值，体现人口分布对风险的非线性放大效应.该分布模式与建筑物易损性（图9a）及人口暴露系数（表6）高度耦合，进一步验证了风险框架中“参数不确定性-暴露度-应急能力”的耦合作用机制.

4 讨论

本研究提出的耦合深度积分法与不确定性分析的承灾体易损性风险评价框架，在深圳光明滑坡案例中展现了良好的适用性.从方法论层面看，该框架通过结合动态数值模型与概率分析，实现了对滑坡运动学特征与致灾效应的量化表征，其技术路径可推广至不同地质条件与堆体规模的固废堆填场风险评估中.具体而言，在地质条件适应性方面，基于连续深度积分的数值模型可显式刻画基底摩擦特性、孔隙水压力效应及地形约束，使其适用于软弱基岩、松散堆积层等典型地质环境；针对堆体规模差异，通过调整网格分辨率与参数抽样范围，可实现对从小型渣土场到大型填埋场的全尺度覆盖.实际使用本文所提框架过程中可遵循以下步骤：首先开展现场地质调查与堆体参数反演，其次构建考虑地形数字化精度的二维或三维数值模型，再次通过拉丁超立方抽样量化参数不确定性，最后耦合承灾体易损性指标完成风险区划.

相较于现有风险评价框架（Galasso et al.，2021；孙小平等，2021；Sun et al.，2021，2025；Zhang et al.，2021），本研究框架的创新性体现在3个维度：其一，将参数不确定性与运动超越概率显式关联，突破了确定性模型的理想化假设；其二，通过冲击效应量化模型实现了建筑物易损性的动力学解析，而非仅依赖经验判据；其三，现有研究多聚焦单一承灾体类型，如Galasso et al.（2021）、Sun et al.（2021a，2025）仅考虑建筑脆弱性，Zhang et al.（2021）仅关注人员伤亡，而本研究创新性耦合建筑物易损性（D_i）与人口暴露系数（V_i），实现了“结构损毁-人员伤亡”的链式风险量化.该框架在评价指标设置上既规避了传统方法耦合有限元结构计算的复杂性（Sun et al.，2021，2025），又通过建筑物与承灾体的物理意义建模确保了科学性.

需指出的是，模型在复杂地形数字化精度与实测数据不确定性处理方面仍存在局限.例如，二维深度积分法受限于维度简化的影响，对三维空间中滑体侧向扩展行为的刻画仍显不足.此外，参数概率分布假设（如对数正态分布）与实际地质条件变异性的匹配度需进一步验证.未来研究可聚焦以下方向：开发三维深度积分模型以更精准捕捉滑坡空间演化特征；引入多源遥感数据与实时监测信息，通过数据同化技术动态修正模型输入参数；构建考虑多物理场耦合（如渗流-应力-变形）的综合性分析框架，以提升对断层破碎带或非饱和渗流主导的复杂致灾机制的解析能力；在数据匮乏地区，可依托类比工程与专家经验构建替代参数分布，但需通过敏感性分析评估输入不确定性对结果的影响.

总体而言，研究成果可为未来“无废城市”建设中的类似堆填场选址优化提供重要的理论参考.后续研究可进一步纳入建筑结构安全与抗灾能力、地质地形稳定性、社区应急管理及公众安全意识等多元因素，多维度协同优化实现城市固废处置设施的安全布局.

5 结论

本研究基于Massflow软件构建了深圳光明滑坡动态数值模型，再现了滑坡从失稳启动、高速运动至堆积演化的完整动力学过程.从冲击压力、冲击荷载、运动物质流速及堆积高度等多维度，系统揭示了滑坡运动物质对建筑物的致灾作用机制.并构建了空间概率模型实现了危险性分区，定义了建筑物承灾风险指标，明确了深圳光明滑坡承灾风险.主要结论如下.

（1）基于深度积分方法，在确定性分析框架下模拟得到滑坡运动距离为1 139 m，与实测值1 100 m高度吻合，有效验证了该数值模型在复杂滑坡运动模拟中的适用性与可靠性；

（2）滑坡运动物质对建筑物的冲击压力与冲击荷载呈现显著的"陡增-缓降"特征，即先快速达到峰值后逐渐衰减至静力状态，峰值冲击压力与荷载随计算点与滑坡体距离的增大呈非线性衰减规律；

（3）创新性提出基于参数不确定性的滑坡运动距离超越概率和建筑物易损性作为双重风险评价指标，构建了耦合冲击效应的承灾体易损性风险评价体系.滑坡实际堆积范围完全覆盖高-极高危险区，其中建筑物A~B位于极高危险区，建筑物C~G位于高危险区，该成果与深圳光明滑坡实际致灾情况高度吻合.并通过对建筑物易损性和人口易损性的评估，进一步明确滑坡致灾风险，表明所提框架在固废堆填场选址优化、风险防控策略制定及应急决策中具有重要科学价值.

参考文献

[1] Armanini,A.,1997.On the Dynamic Impact of Debris Flow.In:Armanini,A.,Michiue,M.,eds.,Recent Developments on Debris Flows.Springer,Berlin,Heidelberg,208-226.

[2] Chai,L.Y.,Ke,Y.,Wang,Y.Y.,et al.,2025.“Earth Macro⁃Circulation” of Bulk Hard⁃to⁃Dispose Industrial Solid Waste for Its Ecological Return.Strategic Study of Chinese Academy of Engineerng,27(3):98-105 (in Chinese with English abstract).

[3] Christen,M.,Kowalski,J.,Bartelt,P.,2010.RAMMS:Numerical Simulation of Dense Snow Avalanches in Three⁃Dimensional Terrain.Cold Regions Science and Technology,63(1-2):1-14.https://doi.org/10.1016/j.coldregions.2010.04.005

[4] Galasso,C.,Pregnolato,M.,Parisi,F.,2021.A Model Taxonomy for Flood Fragility and Vulnerability Assessment of Buildings.International Journal of Disaster Risk Reduction,53:101985.https://doi.org/10.1016/j.ijdrr.2020.101985

[5] Gao,Y.,2018.The Study on Long Run⁃out Mechanism of Artificial Landfill Landslide:Taking Shenzhen Guangming New District Landslide as an Example (Dissertaion).China University of Geosciences,Beijing(in Chinese with English abstract).

[6] Gao,Y.,Yin,Y.P.,Li,B.,2021.Failure Process Simulation Analysis of the Shenzhen “12•20” CDW Landfill Landslide:A Case Study.Arabian Journal of Geosciences,14(12):1094.https://doi.org/10.1007/s12517⁃021⁃07429⁃0

[7] Gao,Y.,Yin,Y.P.,Li,B.,et al.,2019.Post⁃Failure Behavior Analysis of the Shenzhen “12•20” CDW Landfill Landslide.Waste Management,83:171-183.https://doi.org/10.1016/j.wasman.2018.11.015

[8] He,X.R.,Yin,Y.P.,Zhao,L.M.,et al.,2024.Disintegration and Fragmentation Effect of High Position Rock Landslide Debris Flow Based on Large Scale Physical Model Test.Earth Science,49(7):2650-2661 (in Chinese with English abstract)

[9] Iverson,R.M.,Ouyang,C.J.,2015.Entrainment of Bed Material by Earth⁃Surface Mass Flows:Review and Reformulation of Depth⁃Integrated Theory.Reviews of Geophysics,53(1):27-58.https://doi.org/10.1002/2013rg000447

[10] Jiang,S.H.,Xiong,W.,Zhu,G.Y.,et al.,2024.Probabilitic Analysis of Reservoir Landslides Considering the Spatial Variation of Seepage Parameters under the Conditions of Rainstorm and Sudden Drop of Water Level.Earth Science,49(5):1679-1691 (in Chinese with English abstract)

[11] Lacasse,S.,Nadim,F.,2011.Learning to Live with Geohazards:From Research to Practice.GeoRisk 2011.American Society of Civil Engineers,Atlanta,Georgia,USA.64-116..https://doi.org/10.1061/41183(418)4

[12] Li,K.,Cheng,Q.G.,Lin,Q.W.,et al.,2022.State of the Art on Rock Avalanche Dynamics from Granular Flow Mechanics.Earth Science,47(3):893-912 (in Chinese with English abstract).

[13] Li,S.,Peng,C.,Wu,W.,et al.,2020.Role of Baffle Shape on Debris Flow Impact in Step⁃Pool Channel:An SPH Study.Landslides,17(9):2099-2111.https://doi.org/10.1007/s10346⁃020⁃01410⁃w

[14] Ministry of Ecology and Environment of the People’s Republic of China,2024.China Statistical Yearbook on Ecology and Environment (2024).Ministry of Ecology and Environment of the People’s Republic of China,Beijing(in Chinese).

[15] Mo,W.M.,Fang,W.H.,2016.Empirical Vulnerability Functions of Building Contents to Flood Based on Post⁃Typhoon (Fitow,201323) Questionnaire Survey in Yuyao,Zhejiang.Tropical Geography,36(4):633-641,657(in Chinese with English abstract).

[16] Ouyang,C.J.,He,S.M.,Xu,Q.,2015.MacCormack⁃TVD Finite Difference Solution for Dam Break Hydraulics over Erodible Sediment Beds.Journal of Hydraulic Engineering,141(5):06014026.https://doi.org/10.1061/(asce)hy.1943⁃7900.0000986

[17] Ouyang,C.J.,He,S.M.,Xu,Q.,et al.,2013.A MacCormack⁃TVD Finite Difference Method to Simulate the Mass Flow in Mountainous Terrain with Variable Computational Domain.Computers & Geosciences,52:1-10.https://doi.org/10.1016/j.cageo.2012.08.024

[18] Ouyang,C.J.,Zhou,K.Q.,Xu,Q.,et al.,2017.Dynamic Analysis and Numerical Modeling of the 2015 Catastrophic Landslide of the Construction Waste Landfill at Guangming,Shenzhen,China.Landslides,14(2):705-718.https://doi.org/10.1007/s10346⁃016⁃0764⁃9

[19] Peng,C.,Wang,S.,Wu,W.,et al.,2019.LOQUAT:An Open-Source GPU⁃Accelerated SPH Solver for Geotechnical Modeling.Acta Geotechnica,14(5):1269-1287.https://doi.org/10.1007/s11440⁃019⁃00839⁃1

[20] Sun,X.P.,Zeng,P.,Li,T.B.,et al.,2021.Run⁃out Distance Exceedance Probability Evaluation and Hazard Zoning of an Individual Landslide.Landslides,18(4):1295-1308.https://doi.org/10.1007/s10346⁃020⁃01545⁃w

[21] Sun,X.P.,Zeng,P.,Li,T.B.,et al.,2025.Probabilistic Vulnerability Evaluation of Buildings under Landslide Runout Impacts Considering the 3D Dynamic Interaction Processes.Engineering Geology,354:108164.https://doi.org/10.1016/j.enggeo.2025.108164

[22] Sun,X.P.,Zeng,P.,Zhang,T.L.,et al.,2021.Assessment of Exceedance Probability of Landslide Run-out Distance and Hazard Zoning.Geological Bulletin of China,40(9):1560-1569 (in Chinese with English abstract).

[23] Sun,Y.J.,Song,E.X.,2018.Dynamic Simulation of “12·20” Shenzhen Landslide. Chinese Journal of Geotechnical Engineering, 40(3): 441-448 (in Chinese with English abstract).

[24] Tang,H.M.,Li,C.D.,Gong,W.P.,et al.,2022.Fundamental Attribute and Research Approach of Landslide Evolution.Earth Science,47(12):4596-4608 (in Chinese with English abstract).

[25] Wang,D.,Wang,B.,Yuan,W.H.,et al.,2023.Investigation of Rainfall Intensity on the Slope Failure Process Using GPU⁃Accelerated Coupled MPM.Computers and Geotechnics,163:105718.https://doi.org/10.1016/j.compgeo.2023.105718

[26] Wang,K.,Zhang,S.J.,Wei,F.Q.,et al.,2020.A Case Study of the Rapid and Long Runout Landslide at Hong’ao Waste Disposal Site in Shenzhen,China.KSCE Journal of Civil Engineering,24(3):727-739.https://doi.org/10.1007/s12205⁃020⁃1399⁃x

[27] Wang,S.,Zeng,P.,Li,T.B.,et al.,2022.Initiation,Movement and Impact Simulation of Soil Landslide with Material Point Method.Journal of Engineering Geology,30(4):1362-1370 (in Chinese with English abstract).

[28] Xing,A.G.,Wang,G.,Yin,Y.P.,et al.,2014.Dynamic Analysis and Field Investigation of a Fluidized Landslide in Guanling,Guizhou,China.Engineering Geology,181:1-14.https://doi.org/10.1016/j.enggeo.2014.07.022

[29] Xu,Q.,Peng,D.L.,Li,W.L.,et al.,2016.Study on Formation Mechanism of Diffuse Failure Landslide.Journal of Southwest Jiaotong University,51(5):995-1004 (in Chinese with English abstract).

[30] Yan,Y.R.,Liu,Z.Y.,Feng,J.,et al.,2023.Research on Transmission Tower Damage Assessment Caused by Earthquake and Landslide.Technology for Earthquake Disaster Prevention,18(1):107-117 (in Chinese with English abstract).

[31] Yang,H.,Huang,X.J.,Thompson,J.R.,et al.,2016.The Crushing Weight of Urban Waste.Science,351(6274):674.https://doi.org/10.1126/science.351.6274.674⁃a

[32] Yin,Y.P.,Li,B.,Wang,W.P.,et al.,2016.Mechanism of the December 2015 Catastrophic Landslide at the Shenzhen Landfill and Controlling Geotechnical Risks of Urbanization.Engineering,2(2):230-249.https://doi.org/10.1016/J.ENG.2016.02.005

[33] Yubei District Emergency Management Bureau,2023.Mud Collapse Occurs at a Waste Soil Site in Yubei District,Chongqing.2023－12－05.http://www.ybq.gov.cn/bm/qyjj/bmdt/202312/t20231215_12712416.html(in Chinese).

[34] Zanchetta,G.,Sulpizio,R.,Pareschi,M.T.,et al.,2004.Characteristics of May 5-6,1998 Volcaniclastic Debris Flows in the Sarno Area (Campania,Southern Italy):Relationships to Structural Damage and Hazard Zonation.Journal of Volcanology and Geothermal Research,133(1-4):377-393.https://doi.org/10.1016/S0377⁃0273(03)00409⁃8

[35] Zeng,P.,Su,Z.H.,Fang,W.H.,et al.,2022.Typhoon Flooding Loss Assessment in Haikou City Based on High Precision Building Type Data.Journal of Catastrophology,37(4):155-165 (in Chinese with English abstract).

[36] Zhang,S.,Liu,Y.,Bate,B.,et al.,2021.Quantitative Human Risk Analysis of 2015 Shenzhen Dump Failure Considering Influence of Urbanization.Journal of Mountain Science,18(6):1439-1457.https://doi.org/10.1007/s11629⁃020⁃6260⁃7

[37] Zhang,S.,Wang,S.R.,Lei,M.L.,et al.,2024.Landslide Risk Level Assessment of a Landfill:A Case Study of a Landfill in Hangzhou.Environmental Sanitation Engineering,32(1):1-8 (in Chinese with English abstract)

[38] Zhong,X.R.,2024.Study on Sudden Departure Mechanism of Low⁃Potential Energy Landslides on Account of Massing Energy by Spoon⁃Shape Terrain—Taking the Landslide of Hong’ao Village Construction Solid Waste Landfill in Guangming New District of Shenzhen,Guangdong,China.Chinese Journal of Rock Mechanics and Engineering,43(10):2485-2496 (in Chinese with English abstract).

基金资助

国家自然科学基金(42577213)

湖南省地质灾害监测预警与应急救援工程技术研究中心开放课题(HECGD202501)

湖南省教育厅重点项目(23A0015)