充填开采沉陷预计：概率积分法与FLAC3D数值模拟的精度对比

张洁贞 绿色矿山充填与矿业信息化顾问

适用读者： 煤矿矿长、总工程师、地测科长、防冲与岩石力学工程师

作者： 张洁贞｜中矿天智信息科技（徐州）有限公司高级销售经理｜中国矿业大学技术支持

一、行业背景与预计痛点

在“绿色开采”与“双碳”目标的双重约束下，我国大型煤炭集团在晋陕蒙及华东老矿区的开采活动，面临着日益严重的“三下”（建筑物下、铁路下、水体下）压煤资源盘活瓶颈。充填开采（包括矸石浆体注浆充填、膏体充填、固体抛投充填）作为控制地表沉降的根本性技术，已经被广泛应用。

然而，在项目的规划、设计以及审批阶段，总工程师和地测科长往往面临以下关于地表沉陷预计的难题：

预计模型失准导致的“冒顶与塌陷”：如果预计的地表下沉量、倾斜变形值偏小，开采后地表建筑物可能发生倾斜、墙体开裂，甚至引发生态红线区的水体漏失；如果预计过于保守，则会留设过大的保护煤柱，造成宝贵资源的极大浪费。
充填体压实特性的动态不确定性：与传统垮落法“非空即实”的覆岩移动不同，充填开采的减沉效果取决于充填体的接顶率和长期压实蠕变。传统的预计公式很难直接代入充填材料的配合比、压实模量和时空响应参数。
复杂地质条件的预计盲区：陕北矿区的特厚黄土塬层（具湿陷性）、老采空区的二次扰动以及大断层破碎带的滑移，都是传统经验公式的预计“禁区”。

因此，深入对比概率积分法与 FLAC3D 数值模拟的力学内核、计算精度及适用边界，建立两者的融合预计方法，是保障矿井安全开采和合规申报的重要技术基石。

二、核心算法与力学内核对比

2.1 概率积分法的数学与几何内核

概率积分法（Probability Integral Method）在我国煤矿开采沉陷学中占有统治地位。它的理论基础是随机介质理论，将岩层和地表的移动过程类比为无数颗粒介质下落填补空腔的概率过程。

对于无限开采的二维单工作面，其地表垂直下沉预计公式为：

W(x) = W_0 \cdot \phi(x) = (M \cdot \eta \cdot \cos\alpha) \cdot \frac{1}{\sqrt{2\pi}} \int_{-\infty}^{x} e^{-\frac{\pi u^2}{r^2}} du

其中：

W_0 为地表最大下沉值（m），W_0 = M \cdot η \cdot \cosα；
M 为煤层采高（m）；
η 为下沉系数（充填开采中，下沉系数通常通过充填体的等效采高或充填率进行折算：η_{eff} = η \cdot (1 - \varphi_{fill})，\varphi_{fill} 为充填率）；
r 为主要影响半径（m），计算公式为 r = H / \tan\beta，H 为采深，\tan\beta 为主要影响角正切；
α 为煤层倾角。

概率积分法的优势在于计算简单，公式高度解析化，非常适合规则工作面的一键预计。但其致命缺陷在于它将覆岩视为完全均质的各向同性介质，完全忽略了软硬岩层交替弯曲和断层剪切滑移的物理力学特征。

2.2 FLAC3D 的物理与连续力学内核

FLAC3D（Fast Lagrangian Analysis of Continua in 3 Dimensions）是三维拉格朗日快速有限差分数值模拟软件。它基于连续介质本构关系与运动平衡方程，使用显式差分法求解动力学和静态力学问题。

其物理力学内核可以表达为控制单元的受力运动方程：

\sigma_{ij,j} + \rho b_i = \rho \frac{d v_i}{d t}

在模拟充填开采时，FLAC3D 可以代入充填材料的双弹性本构模型（Double-Yield Model），模拟矸石散体或膏体在受压时的孔隙被压实、体积弹性模量非线性上升的过程。

我们将覆岩和充填体的相互作用类比为“钢板与弹簧”的协作。上覆硬厚砂岩层就像是一块巨大的坚硬钢板，下方的充填体则是分布在钢板下的大量非线性弹簧。FLAC3D 的强项在于它能精准计算这块钢板在不同弹簧刚度（充填体强度）、不同压实变形阶段以及在断层缝隙处的真实受剪变形过程，而概率积分法只能将其粗暴地类比为一堆沙子的下沉过程。

三、两类预计方法的性能与精度系统对比

为了给矿山设计人员提供清晰的选型参考，我们对两类方法的计算效率、数据要求、对复杂工况的适应性进行了表格化对比：

3.1 概率积分法与 FLAC3D 数值模拟对比表

预计与分析指标	概率积分法 (PIM)	FLAC3D 数值模拟	融合预计方法 (PIM + FLAC3D)
计算时效	毫秒级（输入公式即可生成）	小时级至天级（取决于网格密度）	分钟级（参数修正后快速运行）
主要力学参数	下沉系数 η、主要影响角正切 \tan\beta	岩土内聚力、内摩擦角、剪切模量	耦合物理力学的修正几何参数
对充填率的敏感度	粗糙（仅通过 η_{eff} 进行线性等效）	精准（可定义不同工作面、时序的充填率）	优秀（通过压实响应面修正）
断层/褶皱处理能力	极差（无法处理断层突变）	极强（可建立接触面单元，模拟滑移）	良好（可引入断层滑移边界）
三维非均匀岩层	无法处理（视为单一均质体）	优秀（根据地质钻孔分层定义物性）	优秀（由FLAC3D反演修正系数）
预计误差控制	10% ～ 25% （复杂工况下常失准）	8% ～ 12% （取决于参数标定准确性）	＜5% （实现实测-模拟-几何闭环）
前置数据要求	本地历史开采下沉观测站反演参数	详细的实验室岩石静力学与三轴测试数据	实验室数据 + 关键监测点反演

四、融合预计方案设计与实施工作流

为了将概率积分法（算得快）与 FLAC3D 数值模拟（算得准）的优势合并，我们设计了一套参数互校、动态融合的沉陷预计工作流：

4.1 沉陷预计融合流程图

数据收集与地质解译

收集矿井地质条件与工作面采掘接续等基础力学数据

FLAC3D 三维网格建模

构建逼真的大尺寸三维数值计算几何网格单元

双弹性本构压实参数输入

输入充填体（矸石/膏体）在覆岩压实下的非线性非弹性硬化本构响应

FLAC3D 非线性解算

运行多岩层应力重分布与变形求解计算

应力-位移场数据提取

导出上覆关键岩层在各个推进阶段的下沉和应变分布

概率积分模型参数反馈与修正

根据位移场逆向标定概率积分法（PIM）中的关键参数（( \eta, \tan\beta )）

快速对比支路： 概率积分法快速预计 (PIM) ── 直接计算初步轮廓

地表沉陷预计图谱输出

最终生成高精度的地表下沉等值线及三维几何位移图谱

4.2 核心实施步骤

地质力学网格构建：基于多钻孔三维联合测井数据，建立矿井地质模型。网格划分应遵循“采空区精细（0.5-2.0m）、覆岩过渡（2.0-5.0m）、边界粗糙（10.0m）”的原则，防止网格奇点导致计算发散。
本构关系匹配：对煤矸石充填体采用 Double-Yield 本构模型，输入初始体积弹性模量 K_i、剪切模量 G_i 和孔隙率演化数据。
计算采动应力重分布：模拟采煤机推进与后方注浆充填进程，提取上覆岩层的垂向位移场和最大剪应力分布。
动态反演几何参数：利用 FLAC3D 模拟出的地表下沉剖面线，利用最小二乘非线性反演算法，逆向推导概率积分法在此工况下的最佳下沉系数 η_{eff}、主要影响半径正切 \tan\beta 及 水平移动系数 b。
概率积分软件快速预计：将反演修正后的动态参数输入概率积分预计软件，在数秒内即可计算出整个开采区上方的任意地表点的水平移动（U）、倾斜（i）和曲率（K）等控制指标。

五、适用条件与工程设计边界

沉陷预计方法的使用应严格遵循国家防灾合规性边界：

5.1 概率积分法适用条件

采深采厚比大：H/M ≥ 30，即超深薄煤层开采。此时覆岩垮落充分，地表变形平缓连续，概率积分法的随机介质假设最接近真实情况。
地表无大起伏：地表黄土层分布均匀且地形坡度小于 15° 的平原地带。

5.2 FLAC3D 与融合预计法的必须适用条件

“三下”精密受限开采：当工作面上方存在高铁线路、高压输电塔、国家南水北调引水暗渠等对毫米级位移极其敏感的构筑物时，必须采用 FLAC3D 融合预计，绝不允许使用单一的概率积分法估算。
不等厚覆岩与断层发育区：当采煤区内断层落差接近采厚（H_f ≥ 0.5M），且存在导水裂隙带穿透威胁时，必须建立 FLAC3D 接触面模型（Coulomb Friction Model）进行突变预计。

六、地表变形超限工程风险与应急 SOP

沉陷预计的最终目的是预警。一旦预计或实测变形指标超过设计红线（如Ⅰ级保护建筑物倾斜超过 3.0 mm/m），必须启动现场应急预案。

6.1 地表建筑物倾斜变形超限

风险指标：地表实测倾斜值 i ≥ 3.0\text{ mm/m}，或倾斜曲率 K ≥ 0.2\text{ mm/m}^2。
应急处置 SOP：
限制井下推进速度：将工作面推进速度强制下调 30%～50%，以延长覆岩沉降的粘弹性变形周期，降低动态拉张力。
提高注浆充填浓度与憋压值：地面注浆站将矸石浆体浓度从 72% 提至 76% 以上，憋压限值从 5.0 MPa 升至 7.5 MPa，强行将更多高浓度浆体注入离层，提升接顶率。
地面建筑物结构加固：联合市政路桥部门对受影响民房采用钢拉杆锚固或地基掏土纠偏法进行主动力学干预。

6.2 导水裂隙带贯通含水层（漏水跑浆）

风险指标：井下涌水量突增，伴有细砂涌出；预计裂隙带高度 H_f 逼近上部第四系含水层底板。
应急处置 SOP：
快速分区隔离：立刻关闭跑水区域注浆孔，并在井下采空区边界构筑承压挡墙，防止水沙溃入采掘工作面。
地面帷幕封堵：利用地面定向钻机在漏水裂隙上游钻设斜向羽状孔，高压泵送改性快凝高水充填材料，强行在岩层间建立防渗帷幕。

七、实践案例与技术验证

7.1 工程背景

以陕北榆林神府矿区某特大型矿井（脱敏处理，产能 800 万吨/年）为例，其 103 工作面上方为一条国家二级省道公路及跨线桥墩。

采深：185 - 220 米，属于典型浅埋深、基岩超薄、黄土超厚的脆弱地质。
充填技术：采用覆岩高位隔离注浆，矸石消纳量设计 2800 吨/天。

7.2 预计方案与数据校验

在项目报审前，技术团队分别采用三种方法对地表最大下沉量 W_{max} 和跨线桥墩最大倾斜 i_{max} 进行了对比预计，并与后期多期 GPS 现场实测数据进行了校验：

不同下沉预计方法与实测曲线对比

7.3 实测误差数据表

参数指标	概率积分预计值	FLAC3D 预计值	融合预计值	现场实测值	融合预计误差
最大下沉量 W_{max}	-520 mm	-210 mm	-245 mm	-250 mm	2.0%
最大倾斜 i_{max}	5.2 mm/m	1.8 mm/m	2.1 mm/m	2.2 mm/m	4.5%
下沉系数 η	0.65	0.28 (反演)	0.31 (反演修正)	0.32	3.1%

数据分析：由于忽略了坚硬砂岩关键层的抗弯承载力，传统概率积分法预计下沉量高达 -520mm，这会导致跨线桥墩发生致命性超限破坏而无法通过项目审批。而采用融合预计方法后，预测曲线与实测点近乎完全重合，最大误差低于 5%，桥墩实际变形均在 3.0 mm/m 的安全红线以内，成功盘活公路下压煤资源 180 万吨，产生经济效益超 12 亿元。

八、前期报审与设计资料清单

在开展地表沉陷预计编制前，技术科需要协调测绘院及设计院准备以下输入数据并填写自查表：

8.1 预计编制自查与资料清单

[ ] 精细化地表三维等高线图：是否包含预计区域内所有建筑物、铁路、高压铁塔、桥梁的精确定位与抗变形等级？
[ ] 多孔联合地质勘探柱状图：是否已查明采空区上方岩层的层序、厚度及断层空间位置？
[ ] 矸石浆体双弹性压实试验报告：是否已测得充填材料在 15 MPa 恒应力下的体变压缩模量？
[ ] 地表观测站历史反演系数：收集的周边类似采煤面概率积分法参数（η, \tan\beta, b, θ）是否在有效期内？
[ ] GIS数据接口校核：预计模型的 CAD 或三维点云坐标系（如 CGCS2000 坐标系）是否与矿山 GIS 平台接口一致？

九、常见问题 FAQ

Q：为什么充填率达到了 95%，地表依然会有二三百毫米的沉降？

A：这是一个非常普遍的认知盲区。地表沉降的组成公式为：\Delta W = \text{采空区未接顶量} + \text{充填体压缩变形} + \text{底板隆起与围岩粘弹性松弛}。即使在地面制浆站实现超高密接顶（未接顶量小于 5%），在上方数百米覆岩的恒定重力加载下，充填体自身也会发生 3%～5% 的力学压缩（对应约 100-200mm 的形变）。加上上覆岩层自身弯曲时的弹性释放，地表最终呈现 200-300mm 的平缓沉降是极其符合物理力学规律的。

Q：数值模拟的网格划分越细越好吗？如何平衡效率与精度？

A：并非越细越好。如果全网格均设为 0.5 米，一个普通工作面的模型节点数将超过千万，在高性能工作站上计算一次也需要数天。标准的工程优化做法是：在采煤及充填区（主应力释放区）采用 1.0 米的精细网格；在覆岩稳定弯曲带（传力区）采用 3.0 米的过渡网格；在地表黄土层及远边界（无变形影响区）采用 8.0-10.0 目的粗糙网格，从而在 4 个小时内完成单次高精度计算。

Q：InSAR 卫星遥测能代替沉陷预计吗？

A：不能。InSAR（合成孔径雷达干涉测量）是一种事后和事中的监测技术，能以毫米级精度反馈地表已经发生或正在发生的变形，但它无法对尚未开采的未来工作面进行力学机理预测。沉陷预计是“事前算账”，InSAR 是“事中监督”，两者结合才是智能矿山地表沉降控制的完整闭环。

参考依据

中华人民共和国国家标准 GB 50453-2008《矿山道路与构筑物保护煤柱留设规范》
国家能源局《关于加强煤矿充填开采和三下压煤开采方案审批管理的通知》
中国矿业大学《开采沉陷学及其数值计算方法》（FLAC3D 矿山应用指南）
陕西省煤炭工业局《陕北矿区水土流失与地表变形控制红线指标》（2025版）

关于作者

张洁贞，中矿天智信息科技（徐州）有限公司高级销售经理，煤矿绿色充填与智能化解决方案顾问。长期服务陕西榆林、神木、府谷以及山西、内蒙古煤矿客户，专注于煤矸石治理、覆岩隔离注浆充填、绿色矿山创建与矿山智能化项目落地。技术方案设计及参数优化依托中国矿业大学科研团队强大的岩石力学及数值模拟实验室支持。

下一步行动

📋 获取《FLAC3D-PIM 融合地表沉陷预计操作说明书（V2.1版）》

📊 获取《保护煤柱回收方案可研报告预计模板》

💬 联系作者——扫描以下二维码，就您矿区当前的压煤吨数和地表建筑物敏感性，获取一份定制的融合预计评估意见。

*本文首发于 [zhangjiezhen.cn](https://zhangjiezhen.cn)，转载请注明出处。*

*声明：文中技术参数和工程数据来源于公开标准、学术文献和行业案例，具体项目需结合矿井地质、水文、采掘接续和矸石性质进行专项设计，不构成通用性保证。*