地下空区顶板破坏特征物理实验及数值模拟

张 泉，陶治臣，王 磊，陈晓青

(1.辽宁科技大学 矿业工程学院，辽宁 鞍山 114051；2.鞍山钢铁集团矿山有限公司，辽宁 鞍山 114051)

由于地下矿床开采或其他地质原因，地下留有大量空区[1]，随着时间的推移导致空区覆岩的移动和地表的下沉。地表下沉不仅会危害人类生命和经济财产，而且也会导致地貌和地质变化，影响环保安全[2-3]。空区顶板破坏及覆岩移动规律一直是业内专家的研究热点。

近年来，随着地下空区安全问题的出现，国内外专家对地下空区顶板的安全厚度和极限跨度进行了大量研究。如苗成文等[4]通过长宽比梁板法和数值模拟方法计算分析了安家岭露天矿地下空区顶板合理安全厚度和极限跨度。李超亮[5]通过FLAC3D和RFPA软件对铁矿地下空区顶板稳定性进行模拟，得到了空区顶板安全厚度与跨度的关系。FREIDIN等[6]通过有限元数值模拟法研究了空区覆岩的应力状态和稳定性，得到了空区顶板厚度和跨度的安全值规律。王文先等[7]通过现场充分调查并采用多种不同理论分析方法计算了空区顶板安全临界厚度和采空区顶板的极限跨度。柳小波等[8]通过一种基于薄板理论的空区顶板厚度计算方法，揭示了空区顶板安全厚度与顶板几何尺寸、矿石力学性质的关系。

空区稳定性涉及因素不仅包含空区顶板厚度和跨度，其他影响因素也非常多，如采矿方法、构造应力、覆岩岩性等，是一个极其复杂的矿山问题，国内外学者通过从不同角度或采取不同方法已经对空区稳定性进行了大量的细致研究。如郭广礼等[9]根据不同岩层的特点建立了空区顶板相似物理模型，通过模型计算得到了空区稳定安全的埋深。张建等[10]根据矿区现场地质条件，在基于工程实践的基础上，建立了不同深度的空区相似物理模型进行空区顶板稳定性分析。张海波等[11]根据矿区现场地质条件，通过有限差分软件FLAC3D对空区进行了模拟研究。杨金林等[12]根据矿区现场地质条件，通过UDEC软件建立了地下空区模型，研究了不同开采条件及外界因素的影响下空区的稳定性。汪杰等[13]基于Kachanov蠕变损伤理论构建了采场顶板蠕变损伤模型，深入分析了空区顶板破坏失稳时间。刘洋等[14]通过综合考虑区域地质水文、岩石强度性等多项因素，建立了基于AHP —模糊评价失稳预测模型的空区稳定性评价指标体系。目前学者对于空区顶板稳定性分析预测大多方法单一，鲜有两种及以上分析方法综合评价空区顶板稳定性问题。

本文基于弹性薄板理论，配制空区顶板的相似模拟材料，通过构建空区模型研究空区顶板的破坏特征，与此同时，采用数值模拟手段从位移、应力等方面研究均布荷载下空区顶板移动变形规律，降低单一因素影响的偏离性，提高采空区顶板稳定性评价的合理性。

1 相似模拟实验

依据王正红[15]基于弹性力学对薄板理论的基本阐述，本文选用薄板理论对空区顶板所产生的力学反应形式进行分析。为便于计算分析，本文模型采用以下假设与规定：①将顶板视为连续均匀介质，空区各板面的强度都是一致且规则分布；②模型中的顶板厚度取其实际厚度的平均值；③忽略由于各向载荷的增加所引起的微小位移；④只考虑空区顶板受到的纵向上方的力。

用C20混凝土制作试件，试件外形为长方体形状，顶部及底部长38 cm，宽28 cm，高度为45 cm，顶板厚度为15 cm，其内部形成采空区，侧帮的厚度为5 cm，每组制作三个模型。

材料规格配比表如表1所示。试件示意如图1所示。

表1 材料规格配比表Table 1 Material specification and proportion

图1 地下空区物理模型Fig.1 Physical model of goaf

试件在辽宁科技大学矿业工程学院实验室单轴压力试验机上进行，实验加载速率设置为0.05 mm/min，采用分布均匀施压的方式，准备12个4 cm×4 cm×1 cm的铁块，放置在顶板的上部，每个铁块的间距为4 cm，并且距顶板的四个边界保持同等的距离，均匀分布在顶板的上方，通过铁块布置进行加载，如图2所示。

图2 均布加载铁块布置情况Fig.2 Layout of uniformly loaded iron blocks

2 实验结果分析

经过多次实验结果对比，发现得到的破坏规律大体一致，因此选用其中一个模型进行试验结果分析。

图3 加载试件应力-应变曲线图Fig.3 Stress-strain curve of loaded specimen

由图3应力-应变曲线可知该模型所能承受的最大应力约为6.4 MPa。加载过程可分为压密阶段、线弹性变形阶段、非线性变形阶段、残余变形阶段四个阶段。其中在位移产生约0.2 mm之后发生突变，模型产生大幅度位移。

从图4中可以看出，顶板破坏所呈现的裂隙发育现象为在顶板中心区域出发形成一条竖直裂隙直达顶板长边两侧中部，在到达两侧时呈放射状向左右两侧发育，裂隙发育直达侧帮部位，最终应力达到峰值时模型彻底破坏。由此可知，空区顶板受纵向加载的作用时，最易破坏的部位为中心部位，应在其后的采矿工作中，给予一定的支护手段以达到工业安全保障。

图4 空区顶板破坏情况Fig.4 Roof damage of empty area

3 地下空区顶板数值模型构建

由于开展研究的地下空区顶板、围岩都是地质材料，在力学角度来讲，顶板可以看作各向同性的弹塑性连续介质，在FLAC3D数值模拟计算中需选用莫尔库伦破坏准则来进行模拟[16]，其力学模型为：

(1)

ft=σ3-σt

(2)

式中：σ1为最大主应力，MPa；σ3为最小主应力，MPa；σt为抗拉强度，MPa；C为内聚力，MPa；φ为内摩擦角，(°)。

其中，在地下空区顶板处于拉应力作用下，当拉应力大于模型材料的抗拉强度时，即：ft大于0时，顶板将会产生受拉破坏。

为了更准确地得到顶板各部位的变化情况，在顶板中部切面设置监测点，探究顶板受到纵向载荷时，各个区域的应力场变化情况，通过模拟计算得出等值云图从而对顶板所产生的位移、应力应变等进行分析。基于莫尔—库仑强度理论的地下空区数值模型如图5所示，设置物理力学参数如表2所示。

图5 地下空区初始数值模型Fig.5 Initial numerical model of goaf

表2 数值计算岩体物理力学参数选取Table 2 Selection of physical and mechanical parameters of rock mass for numerical calculation

4 模拟结果分析

4.1 主应力分布规律

图6为地下空区顶板在加载过程中分为四个阶段垂直走向的最大主应力分布情况。通过应力云图中主应力变化情况可以看出，顶板上出现了拉应力集中区域，主要集中在顶板的中央部位。随着载荷的持续增加，在顶部周围形成一圈较为明显的压应力，压应力值逐渐增大；而在顶板的中部，最大拉应力由环状慢慢向中心收缩，最终集中在中心部位。由此可以推导出顶板产生拉应力最大的部位为顶板的中心部位。

图6 顶板最大主应力分布情况Fig.6 Distribution of maximum principal stress of roof

4.2 正应力分布规律

随着顶板接受载荷的递增，每个阶段所呈现的正应力情况有所改变，空区顶板的水平正应力分布等值云图如图7所示。顶板上水平正应力分布通过X方向应力云图进行分析，随着载荷的增加，顶板上最大拉应力逐渐向顶板中心收缩，且数值逐渐增大；在Y方向上的侧帮上方出现两段压应力区域，同时随着载荷增加而数值增大，但应力最大的区域始终是在顶板的中心区域，可推导顶板中心部位是最易产生破坏的，破坏的形式是受拉破坏。

图7 水平正应力分布情况Fig.7 Horizontal normal stress distribution

4.3 剪应力分布规律

图8为四个阶段加载后地下空区顶板的最大剪应力云图，剪应力分布存在以下规律：从顶板上方施加载荷开始，剪应力也会出现应力集中的现象；剪应力主要在顶板上形成一个“X”型分布，而不是像正应力集中于顶板的中部。随着载荷的增加，剪应力值也在不断增加，从整体来看，应力变化的幅度较小。

图8 顶板剪应力分布情况Fig.8 Distribution of roof shear stress

4.4 空区顶板变形与位移特征

在空区顶板内部布置一系列的监测点，伴随着载荷的增加监测顶板各个部位的下沉量与时间的关系，监测点的坐标为(x，14，42)，x属于{0～19}，在x轴上从中心部位开始每隔3 cm布置一个监测点，由中心出发分别命名为A、B、C、D、E，通过模拟计算自动输出每个监测点在垂直方向的位移。随着载荷的增加，顶板各个位置都出现了下沉量，顶板中央部位的下沉量相比更为明显，且越接近顶板中心，产生的位移量差值越小；由于顶板出现位移情况为左右对称，顶板两端的位移量相比于中央部位较小，因此接近侧帮的部位呈现近乎水平的状态。表3为空区顶板的垂直位移情况。

表3 顶板监测点随时间的位移量Table 3 Displacement of roof monitoring points with time

5 结论

1)物理实验结果表明，空区模型在加载作用下，顶板出现了受拉破坏，顶板中央区域最容易破坏，因此，地下空区处理应该在中央区域加强支护或加固，充填采矿法中央区域矿体的回采应当保留矿柱，矿柱不宜二次回采。

2)数值模拟实验结果表明，随着荷载增加，必然会导致空区顶板侧帮应力出现变化，顶板中央区域存在拉应力集中区，这是造成顶板出现受拉破坏的根源所在，顶板裂隙发育完全时将发生破裂甚至冒落。在载荷增加过程中，剪应力也会缓慢集中，在剪应力集中区域发生剪切破坏。

3)监测结果表明，随荷载增加，顶板各部位出现的下沉量逐渐增大，顶板中央部位的下沉量更为明显，且越接近顶板中心，三者产生的位移量差值越小。监测结果与物理试验和数值模拟实验结果较为吻合。