基于数值模拟的房柱式采空区关键层-岩柱系统破坏机理分析

管永伟，魏义强，杨振江，李志军，武 军

(山西省交通科学研究院，山西 太原 030006)

基于数值模拟的房柱式采空区关键层-岩柱系统破坏机理分析

管永伟，魏义强，杨振江，李志军，武 军

(山西省交通科学研究院，山西 太原 030006)

为了提高浅埋房柱式采空区稳定性评价的精度和可靠性，本文以吕梁某高速公路桥梁下伏采空区为研究区，建立采空区二维模型，运用有限元数值模拟的方法，研究关键层厚度对关键层-岩柱系统应力场的影响。研究发现关键层-岩柱系统发生了应力重分布，应力集中主要位于岩柱的帮角、关键层跨中表面及岩柱正上方关键层上表面。随关键层厚度减小，岩柱应力呈先增后减趋势，关键层应力先缓慢增大后急剧增大。建立关键层-岩柱系统稳定性与关键层厚度之间的定量关系，确定关键层局部破裂、岩柱失稳、岩柱保持长期稳定时的关键层临界厚度分别为3.2 m、5.8 m和6.7 m。将数值模拟成果应用到实际工程中，发现与实际资料的吻合度较高。

房柱式采空区；关键层-岩柱系统；数值模拟；地应力；临界厚度

0 引言

山西省矿产资源丰富，资源开采后遗留的采空区诱发地质灾害问题对公路工程影响越来越突出[1]。我国相当一部分金属矿藏仍采用房柱式开采工艺，开采后，采空区留存有大量的岩柱，岩柱和岩层顶板组成结构系统，在静载作用下，其稳定性不仅取决于岩柱的整体强度，还与顶板的力学特性密切相关[2-5]。国内外对浅埋房柱式采空区岩层覆岩运动规律、残留岩柱稳定性及浅埋岩层顶板控制技术进行了一系列的研究，取得了丰富的成果。Schaubs P M等[6]在Budryk-Knothe 理论基础上引入Fourier 二维积分变换形成数值模拟方法研究采空区的稳定性。贺广零等[7]依据板壳理论和非线性动力学理论对采空区岩(煤)柱-顶板系统失稳机理进行了研究。杨敬轩、孟达、史红等[3，4，8]分别运用数学方法、力学方法、数值模拟方法对房柱式开采顶板、煤柱的承载力及煤房、煤柱尺寸等方面进行了系统的研究并对其稳定性进行评价。杨路平[9-10]等进行了下伏采空区桥隧、路面变形及应力耦合分析。

由于浅埋房柱式采空区具有浅埋深、薄基岩、上覆厚松散赋存特征，岩柱受力具有特殊性[11]。单一岩柱失效或者顶板失稳将使载荷转移到其邻近岩柱上并引起该岩柱过载，引发采空区中的岩柱群失稳多米诺效应，引起覆岩破坏和地表塌陷[4]，因此对房柱式采空区顶板-岩柱系统稳定性进行研究具有非常重要意义。本文以山西某铝土矿采空区为例，从关键层理论出发，将关键层-岩柱群系统作为研究对象，考虑二者之间的相互的作用，运用有限元数值模拟的方法对关键层-岩柱的应力分布进行分析，探索关键层破断、岩柱失稳时关键层临界值，最后，进行工程实例验证。

1 数值模拟

1.1 模型建立

在关键层理论模型的框架内，建立关键层-岩柱系统二维有限元模型(图1)。模型建立时，假设：(1)岩柱宽、高、房跨均为定值，关键层(灰岩)厚度H为变量；(2)岩柱及顶板关键层为均匀、弹性各向同性连续介质；(3)不考虑构造应力；(4)垂直荷载简化为均布荷载P施加于关键层，取P=1.29 MPa。图中模型尺寸按照某铝土矿实际参数设置，模型长52 m，岩柱宽5 m、高3 m，房跨8 m，底板为厚层灰岩厚15 m。模型的地层物理力学参数见表1。选用弹塑性本构关系，匹配摩尔-库仑准则进行计算。

图1 关键层-岩柱系统模型Fig.1 Model of the key strata-pillar system

岩层杨氏模量E/GPa泊松比μ密度ρ/(g·cm-3)内聚力C/MPa内摩擦角Φ/(°)关键层250192826423岩柱040242113150

1.2 仿真结果分析

通过不断改变模型关键层厚度参数进行模拟，获得多组关键层-岩柱系统应力场特征，为了消除边界效应，取中间两跨进行分析。图2为H=3 m和H=6 m时的应力分布云图。图中，矿体开采遗留岩柱，改变了局部的区域应力场，应力发生了重分布。岩柱的帮角、与帮角接触的关键层底部、跨中表面及岩柱正上方的关键层上表面等位置均出现了应力集中现象，对比图2(a)、图2(b)发现，H=6 m时，岩柱处于应力升高区，关键层总体上为应力降低区；H=3 m时，岩柱应力明显降低，关键层应力升高，应力集中现象也表现的更为显著和复杂，此时关键层出现明显“X”型形态的力学特征，预示着其可能已经发生了断裂。

为了方便对比，将关键层下表面(岩柱上表面)位置处的应力大小绘制到其对应水平位置，图3显示，随着关键层厚减小，关键层应力逐渐增大，高应力异常逐渐向跨中位置靠拢，H=3 m，跨中下表面应力σ下=4.5 MPa；岩柱应力先增大后减小，高应力异常先向岩柱中心位置靠拢后又发散。提取关键层上表面应力数据并绘制应力分布曲线，以H=6 m和H=3 m为例(图4)。图中显示，岩柱正上方关键层顶部及跨中顶部应力较为集中，随后应力值均出现明显的跌落，最小应力值出现在岩柱帮角正上方位置，分析发现岩柱正上方关键层顶部在水平方向上应力梯度高于跨中；结合图2可知，岩柱正上方关键层顶部应力集中程度高于跨中底部(二者均为受拉区)。因此，岩柱正上方关键层顶部位置优先拉张开裂，原因在于岩柱正上方关键层受到岩柱的反力作用，能够形成更大的扭矩，更容易拉裂。

1.3 临界值确定

图5给出了岩柱正上方的关键层上表面与岩柱帮角应力随关键层厚度变化的曲线，图中显示，H=9 m～3 m，随关键层厚度减小，岩柱应力呈现先增后减趋势，关键层应力呈快速增大趋势。H=5.8 m时，岩柱出现应力最大值，σ柱max=5.3 MPa，说明此时岩柱整体已进入屈服状态，结合云图2可知，岩柱的破坏是两侧塑性区逐步向核心弹性区扩展的过程，表现为片帮缩颈；关键层应力增长速率经历了缓慢增大、急剧增大后又放缓的过程，H=3.2 m时增长速率最大，表明顶板关键层发生了局部破裂。

综上，H≤3.2 m，即σ关键层≥4.46 MPa，关键层局部破裂，岩柱应力转移降低，岩柱处于弹性状态；当3.2 m5.8 m，岩柱部分完好，关键层也处于弹性变形状态。据此得到岩柱强度[σ柱]=5.3 MPa，关键层局部拉裂强度[σ关键层]=4.46 MPa。关键层局部受拉屈服不代表整个关键层的垮落，而是随时间或者关键层厚度的继续减小，屈服区不断扩展和延伸，直到关键层形成了完整的“O-X”破裂带，因此取关键层临界破裂厚度H=3.2 m；岩柱的稳定性是由岩柱本身的强度和岩柱应力所决定的，当岩柱应力超过其强度时，将造成岩柱失稳破坏。可用稳定性安全系数(fs)来评定岩柱的稳定性(fs=[σ柱]/σ柱),fs>1.5，岩柱能保持长期稳定，此时σ柱=3.53 MPa，对应H=6.7 m。因此，关键层局部破裂、岩柱失稳时关键层临界安全应力分别为4.46 MPa和3.53 MPa，临界厚度分别为3.2 m 和6.7 m。

图2 关键层-岩柱系统应力场的分布Fig.2 Stress distribution of the key strata-pillar system

图4 关键层上表面应力分布曲线Fig.4 Plots of stress distribution in top surface of the key strata

图5 关键层-岩柱应力随关键层厚度变化曲线Fig.5 Plots of key strata-pillar distribution Vary with the thickness

2 实例分析

2.1 研究区概况

吕梁环城高速公路某大桥穿越某铝土矿采空区，2016年6月起桥梁南侧K6+441-K6+747段即0号桥台至7号桥墩之间，出现多处地表裂缝和地表沉降现象，不良地质现象初步判定为采矿形成的采空区坍塌导致。钻探资料显示，桥梁下方多处掉钻，主要集中在5号桥墩附近，图6为5号桥墩下伏采空区分布情况。

矿井2008年始采，2016年4月停采，采空区埋深51.3～62.8 m，采空区面积210 000 m2。矿区为典型的黄土高原中低山丘陵地貌。地层从上到下有：第四系黄土，石炭系泥岩、灰岩、铝土岩及奥陶系灰岩。矿井开采石炭系上统本溪组铝土矿层，矿层均厚3 m，倾角平缓，采用房柱式开采工艺，采场之间留连续矿柱，矿柱尺寸5 m×5 m，房跨8 m。地质资料显示，桥梁下伏采空区基本顶(老顶)为坚硬、层状灰岩，稳定发育，厚度6.5～7.4 m，对采场上覆直至地表的全部地层的活动起控制作用，称为关键层。

图6 桥梁下伏采空区剖面图Fig.6 Profile of goaf under bridges

2.2 模拟结果验证

图6显示，5号墩基桩桩端下卧泥岩层距关键层约为20 m，可认为均布载荷q0施加于关键层[12]。上覆地层容重γ=2×104N/m2，厚度h=55 m，则作用在顶板上的总的均布载荷q=q0+γh=1.29 MPa。基于第1节数值模拟结果可知，岩柱保持长期稳定和整体失稳破坏时，关键层临界厚度分别为6.7 m和5.8 m，5号墩下伏关键层最小厚度为6.5 m，介于二者之间，因此桥墩不能保持长期稳定。将5号桥墩连续120 d变形监测结果绘制成曲线图(图7)。图中显示，5号墩累计变形量达23.2 mm，虽然一直处于沉降状态，但是桥梁不会产生瞬间失稳；桥梁水平位移量为11.05 mm，方向呈SSW向，变化速率逐渐平缓，说明桥梁南侧大面积采空区塌陷区对5号桥墩的影响逐渐趋于稳定。综上，采空区已严重威胁大桥的安全，须尽早处治。塌陷区内56个勘察钻孔，8个出现铝土岩矿柱，可能是该范围内关键层厚度小于其临界破裂厚度，导致顶板关键层优先失稳，引起地表塌陷，而岩柱却保持着完好状态。另外据井下观测资料显示，部分岩柱存在严重片帮现象。因此数值模拟结果与实际情况吻合。

图7 5号桥墩变形监测结果的统计Fig.7 Plots of deformation monitoring for the No. 5 pier

3 结论

建立房柱式采空区二维模型，运用有限元数值模拟的方法，研究关键层厚度对关键层-岩柱系统应力分布的影响，通过对数值模拟结果与实例资料分析，获得如下几点结论：

(1)随关键层厚度减小，岩柱应力呈先增后减趋势，关键层应力先缓慢增大后急剧增大。

(2)建立关键层-岩柱系统稳定性与关键层厚度之间的定量关系；关键层局部破裂、岩柱失稳时关键层临界应力分别为4.46 MPa和5.3 MPa，临界厚度分别为3.2 m和5.8 m；岩柱保持长期稳定的临界应力和临界安全厚度分别为3.53 MPa和6.7 m。

(3)岩柱正上方关键层上表面和跨中位置下表面均为拉应力集中区，但前者应力集中程度高于后者，优先拉张开裂；岩柱四个帮角应力集中，其失稳是由帮角塑性区向中心弹性区扩展，逐渐片帮缩颈导致。

(4)模拟结果与实际资料相吻合，关键层-岩柱系统应力有限元数值模拟方法可以作为研究浅埋房柱式采空区稳定性的可靠手段。

[1] 童立元，刘松玉，邱钰，等．高速公路下伏采空区问题国内外研究现状及进展[J]．岩石力学与工程学报，2004，23(7)：1198-1202．

TONG Liyuan，LIU Songyu，QIU Yu，et al．Current research state of problems associated with mined-out regions under expressway and future development[J]．Chinese Journal of Rock Mechanics and Engineering，2004,23(7)：1198-1202．

[2] 孙琦，张向东，杜东宁，等．浅埋采空区对路基稳定性影响的数值模拟[J]．中国地质灾害与防治学报，2015，26(2)：127-131．

SUN Qi，ZHANG Xiangdong，DU Dongning，et al．Numerical simulation of the impact of shallow seam goaf on highway safety[J]．The Chinese Journal of Geological Hazard and Control，2015，26(2)：127-131．

[3] 孟达，王家臣，王进学．房柱式开采上覆岩层破坏与垮落机理[J]． 煤炭学报，2007，32(6)：577-580．

MENG Da，WANG Jiachen，WANG Jinxue．Mechanism on the failure and caving of roof strata in pillar and house mining[J]．Journal of China Coal Society，2007,32(6)：577-580．

[4] 杨敬轩，刘长友，杨宇，等．浅埋近距离煤层房柱采空区下顶板承载及房柱尺寸[J]．中国矿业大学学报，2013,42(2)：161-168．

YANG Jingxuan，LIU Changyou，YANG Yu，et al．Study of the bearing mechanism of the coal roof and the dimension selection of the room and pillar in the shallow and close distance coal seam[J]．China University of Mining and Technology,2013,42(2)：161-168．

[5] 罗周全，刘晓明，吴亚斌，等．基于Surpac和Phase～2耦合的采空区稳定性模拟分析[J]．辽宁工程技术大学学报(自然科学版)，2008，27(4)：485-488．

LUO Zhouquan，LIU Xiaoming，WU Yabin，et al．Study on cavity stability numerical simulation based on coupling of Surpac and Phasez[J]．Journal of Liaoning Technical University(Natural Science)，2008，27(4)：485-488．

[6] Schaubs P M，Ord A，German G H．Scenario testing of fluid-flow and deformation during mineralization：from simple to complex geometries FLAC and numerical modeling in aeromechanics[M]．Rotterdam：A．A．Balkema，2003：63-70．

[7] 贺广零．基于厚板理论的水平煤层顶板临界厚度的分析[J]．地下空间与工程学报，2009，5(4)：659-663．

HE Guangling．Determination of critical thickness of stiff roof in coal mine based on thick plate theory[J]．Chinese Journal of Underground Space and Engineering，2009，5(4)：659-663．

[8] 史红，姜福兴．采场上覆大厚度坚硬岩层破断规律的力学分析[J]．岩石力学与工程学报，2004，23(18)：3066-3069．

SHI Hong，JIANG Fuxing．Mechanical analysis of rupture regularity of hard and massive overlying strata of longwall face[J]．Chinese Journal of Rock Mechanics and Engineering，2004，23(18)：3066-3069．

[9] 杨路平，潘建平，郭弯弯，等．高速公路下伏采空区稳定性评价及治理[J]．重庆交通大学学报(自然科学版)，2013，32(6)：1195-1198．

YANG Luping，PAN Jianping，GUO Wanwan, et al．Stability evaluation and treatment of goaf under expressway [J]．Journal of Chongqing Jiaotong University (Natural Science)，2013，32(6)：1195-1198.

[10] 杜东宁，张向东，张淑坤，等．下伏采空区对地表高速公路变形影响研究[J]．广西大学学报(自然科学版)，2015，40(4)：984-988．

DU Dongning，ZHANG Xiangdong，ZHAND Shukun，et al．Influence of goaf on highway deformation[J]．Journal of Guangxi University (Natural Science)，2015，40(4)：984-988．

[11] 屠世浩，窦凤金，万志军，等．浅埋房柱式采空区下近距离煤层综采顶板控制技术[J]．煤炭学报，2011，36(3)：366-370．

TU Shihao，DOU Fengjin，WAN Zhijun，et al．Strata control technology of the fully mechanized face in shallow coal seam close to the above room-and-pillar gob[J]．Journal of China Coal Society，2011，36(3)：366-370．

[12] 王玉涛，毛旭阁．高速公路桥梁下伏采空区稳定性分析方法探讨[J]．公路，2013(10)：93-96．

WANG Yutao，MAO Xuge．Analysis of goaf stability under highway bridges [J]．Highway，2013(10)：93-96．

Failuremechanismofthekeystrata-pillarsysteminroomandpillargoafbasedonanumeriealsimulation

GUAN Yongwei，WEI Yiqiang，YANG Zhenjiang，LI Zhijun，WU Jun

(ShanxiTransportResearchInstitute,Taiyuan，Shanxi030006,China)

In order to improve the precision and reliability of stability evaluation in shallow buried room and pillar goaf，a 2D geological model which based on geological data of the goaf under expressway with bridges in Lvlang was established and was applied to study the effect of key strata thickness on in situ stress of the key strata-pillar system by using finite element. The results show that the key strata-pillar system causes the in situ stress to be distributed. the stress in angle of pillars , mid-span and the surface of the strata-pillar directly above pillars are highly concentrated, With the thickness of key strata decreases, pillar stress firstly increases and then decreases, and the key strata increases slowly at first and then sharply. Quantitative relations between key strata thickness and the key strata-pillar system are established. The critical key strata thickness of local failure in key strata，the instability and stability with a long-term of pillar are determined to be 3.2 m，5.8 m and 6.7 m. The conclusions above are applied to the practical project and the results are highly consistent with actual data.

room and pillar goaf; key strata-pillar system; numerical simulation; in situ stress; critical key strata thickness

10.16031/j.cnki.issn.1003-8035.2017.04.10

TD801

A

1003-8035(2017)04-0059-05

2017-01-24;

2017-03-12

山西省交通运输厅科研计划项目(07-02)

管永伟(1988-)，男，河南濮阳人，硕士，主要从事岩土工程、岩土力学数值计算研究。E-mail：guanyongcumt@163.com