目前,我国的主要能源是煤炭,据不完全统计,仅国有煤矿的生产矿井“三下”压煤量达到140多亿t[1-3]。随着国家经济的快速发展,煤炭需求量增加,大量的煤炭从地下采出[4],破坏了岩体内部原有的力学平衡状态,使上覆岩层不同程度地变形和破坏[5]。当开采面积达到一定范围即一般当采空区的长和宽(沿煤层倾向的水平投影长度)都超过平均采深的0.2~0.3倍时,地表开始移动,起始采场附近的岩层移动和变形将扩展到地表[1],此时的地表移动和变形将影响到位于开采影响范围内的房屋建筑、工程、河流、及管线,会改变它们原有的状态,甚至破坏[6-9]。因此,对采空区地表沉陷的研究具有重大的意义。
针对采空区地表沉陷及预测问题,国内外专家开展了大量的研究,并取得了显著的研究成果[10]。如典型曲线法和剖面函数法[11-12]、概率积分法[13-14]、基于连续介质力学方法以及影响函数法等[15-17]。但由于现场情况多变,这些预测方法适用于某种特定情况,导致预测结果不是很理想,给研究开采沉陷的学者及技术人员带来挑战。但数值模拟只需在计算机上进行,可以施加试验方法达不到的条件。赵卫强等[18]应用FLAC3D对村庄下充填开采地表沉陷规律进行了研究;郭玉芳等[19]利用FLAC3D数值模拟分析了厚松散层开采条件下地表移动变形规律;邓小龙等[20]通过数值模拟确定了采区开采结束后的地表沉陷边界和地表塌陷体积。不过,他们没有结合理论计算与现场实测,仅应用数值模拟预测,结果不具有说服力。为此,笔者以店坪煤矿5-210工作面为工程背景,采用理论计算、数值模拟与现场实测相结合的方法,更加准确地预测采空区地表沉陷范围,为煤炭的开采及房屋建筑的安全提供重要的理论依据。
5-210工作面位于店坪煤矿工业广场、王家庄村和店坪村的南面,麻塔村的北面,南岭上村的东面。从西向东依次穿过里沟、闫家沟、柳沟、麻塔梁。东西走向长度1 465 m,倾向长度240 m,煤层近水平、薄基岩、煤厚平均3.1 m,平均埋深200 m,属于浅埋煤层。最低处位于里沟中,标高为+1 075 m,雨季沟内可能有部分积水渗入地表,最高处位于麻塔梁上,标高为+1 223 m。地表沟壑纵横,麻塔梁、走马梁顶部有少量耕地,属丘陵梯田,其余为坡地。因地表附近有工业广场和村庄,需研究5-210工作面回采所引起的地表沉陷是否影响建筑物的稳定性。
考虑建模需要,根据空间解析几何及岩层移动相关理论,建立地表下沉移动轨迹坐标系,如图1所示。坐标系原点设定于采空区正中央,此处为地表沉量最大点,走向方向为x轴,沿工作面推进方向为正,反之为负,倾向方向为y轴,回风巷方向为正,运输巷方向为负,z轴为垂直方向下沉量,边界应满足:
x2/a2+y2/b2=1
式中:x、y为走向和倾向方向变量;a、b为椭球圆台的基本参数。
下沉移动范围条件:
0≤x2/a2+y2/b2≤1
图1 地表下沉移动轨迹坐标系
Fig.1 Surface sinking movement coordinate system
1)压实稳定区。当地表沉陷点在x2/lz/2-Lzd-Lzs2+y2/lq/2-Lqd-Lqs2≤1范围内时,处于压实稳定区,该区域类似于“椭圆盘”底部,在此范围内地表沉陷量基本一致,其中:lz为走向长度,m;lq为倾向长度,m;Lzd为走向低应力区长度,m;Lzs为走向应力升高区长度,m;Lqd为倾向低应力区长度,m;Lqs为倾向应力升高区长度,m。
2) 载荷影响区。当地表沉陷点x2/lz/2-Lzz-Lzy2+y2/lq/2-Lqz-Lqy2≥1且x2/lz/2-Lzz2+y2/lq/2-Lqy2≤1范围内时,处于载荷影响区,该区域类似于“椭圆盘”的边部,在此范围内任意一点都可将其等效转化到坐标系x、y轴上,如图2所示,转化后可得沉陷点与煤壁水平距离L。将转化后的沉陷点与煤壁水平距离L代入载荷影响区碎胀系数表达式中,得出该点所对应的碎胀系数,最终将其代入地表沉陷表达式中,得出该点沉陷量。其中:Lzz为走向长度在z轴方向上的长度,m;Lzy为走向长度在y轴方向上的长度,m;Lqz为倾向长度在z轴方向上的长度,m;Lqy为倾向长度在y轴方向上的长度,m。
图2 坐标转换示意
Fig.2 Schematic diagram of coordinate transformation
3)自然堆积区。当地表沉陷点在x2/lz/2-Lzz2+y2/lq/2-Lqz2≥1且x2/lz/22+y2/lq/22≤1范围内时,处于自然堆积区,在此范围内时,垮落带内岩体并未全部垮落,靠近煤壁一侧形成残余悬臂梁结构,支撑上覆砂土层,在上覆砂土层载荷作用下只产生微弱下沉。
岩石碎胀系数Kp是岩石的碎胀性可用岩石破碎后处于松散状态下的体积与岩石破碎前处于整体状态下的体积之比。岩石破碎以后的体积将比整体状态下增大,这种性质称为岩石的碎胀性。当顶板垮落到一定高度,即Σh=m/Kp-1时,就可以完全充填采空区,其上覆岩层的活动对工作面就没有明显的动压影响。因此,碎胀系数对工作面顶板管理有重要的意义。根据岩柱法计算应力平稳区垮落岩体应力,覆岩物理力学参数如下,砂土层的平均容重γs=18.5 kN/m3,基岩层的平均容重γj=23.5 kN/m3,砂土层厚度hs=30 m,基岩层厚度hj=17 m,砂土层内摩擦角φs=27°。
1)压实稳定区地表沉陷量。将数据代入垮落岩体所受覆岩应力式(1),计算出压实区垮落岩体受上覆基本顶及砂土层应力载荷为
(1)
计算得σzp=0.952 3 MPa。
根据文献[17]中全断面沉陷模型理论,由于煤层近水平,用走向半长代替a,则2a=1 465 m。选取矿井5-210工作面采空区为研究对象,根据工程类比法确定相关数据,将压实区垮落带内泥岩层及砂岩层所受应力载荷代入应力-碎胀系数关系式中,可得出压实稳定区范围内垮落岩体碎胀系数。
(2)
(3)
其中:ks、kn分别为砂岩、泥岩的碎胀系数;计算分别为1.012、1.008。
W=M-∑hiki-1
(4)
其中:W为地表沉陷量,m;M为煤层开采厚度,m;hi为垮落带各岩层厚度,m;ki为垮落带内各岩层碎胀系数。将压实区垮落带内各岩层碎胀系数及各岩层厚度代入浅埋采空区地表沉陷量表达式(4)中,得出压实区地表沉陷量W=2.604 m
2)应力升高区地表沉陷量。将应力升高区裂隙带内基岩以及上覆砂土层对垮落岩体应力与垮落岩体所受自身重力之和,分别代入砂岩、泥岩应力-碎胀系数拟合公式中,可得出应力升高区范围内砂岩及泥岩的碎胀系数与采空区边界距离关系式。
(5)
(6)
式中:kg、kc分别为垮落岩体、煤层地基系数,GN/m3;Ld为低应力区长度,m;E为基岩弹性模量,GPa;I为基岩梁惯性矩,m4。
应力升高区随与煤柱距离增大地表沉陷量增大,当应力升高区覆岩垮落达到基岩层后,地表沉陷量与煤柱距离关系变化率骤然减小,最终沉陷量与压实区沉陷量基本保持一致,约为2.6 m。
3)自然堆积区沉陷量。在此范围内,垮落岩体应力较小,呈松散堆积状态,碎胀系数最大,根据文献[17]中松散堆积区垮落岩体碎胀系数与煤壁位置关系公式:
(7)
式中:kz为自然堆积区岩石碎胀系数;c、d为回归系数;Δc为补偿参数。
垮落带内岩体并未全部垮落,靠近煤壁一侧形成残余悬臂梁结构,支撑上覆砂土层,在上覆砂土层载荷作用下只产生微弱下沉。
综上所述,采空区中部压实稳定区范围内地表沉陷量最大,达到2.604 m。
根据店坪煤矿5-210工作面所处煤岩层地质条件构建数值模拟计算模型,模型长800 m,宽400 m,高300 m,矿体埋深平均为200 m,矿体厚度为3.1 m,工作面长度240 m,推进长度500 m。模型采用Mohr-Coulomb计算准则,选用的岩石力学参数见表1。模型上表面为应力边界,煤层上部模拟了200 m上覆岩层,模型上部边界施加载荷为5 MPa,模拟上覆岩层自重边界,约束条件取两侧为限制水平方向位移的滑动支座,底部为限制垂直方向和水平方向位移的固定支座。模型共划分126 198个网格节点,118 400个节点等参单元。
根据上述建立的模拟模型及采场结构参数进行模拟,得到了店坪煤矿5-210工作面采空区地表沉陷的模拟结果如图3所示。地表沉陷的最大值为2.78 m,实测的地表沉陷值约为2.69 m,数值模拟结果与实测值的相对误差为3.34%。由此可知,数值模拟结果与实测结果有良好的耦合性。所以说,FLAC3D数值模拟方法对于预测地表沉陷是可行的。
表1 煤岩物理力学性质参数
Table 1 Physical and mechanicalproperties of coal and rock
岩层抗压强度/MPa抗拉强度/MPa黏聚力/MPa内摩擦角/(°)体积模量/GPa剪切模量/GPa石灰岩48.91.772.02404.002.40细粒砂岩65.61.402.73428.154.43泥岩56.70.931.07352.291.065煤14.91.101.60384.181.40
图3 地表沉陷最大值
Fig.3 Maximum surface subsidence
将FLAC3D数值模拟结果导入Tecplot 360软件,该软件具有强大的数据分析和可视化处理功能,可以使模拟或者实验结果快速而精准地以图表或者动画的方式显示,还可以分析复杂数据,排版多种布局和以专业的图像和动画交流结果,方便工程师和科研人员在最短的时间内对大量数据进行可视化分析,为设计和研究提供科学的依据。通过查阅资料获得,为了满足地表建筑物的稳定以及煤炭开采的安全,地表沉降控制基准确定为39 mm,对模型进行赋等值线并切片,如图4所示。
图4 地表沉陷切片
Fig.4 Surface subsidence section
地表沉陷39 mm的位置与采空区边界连线在煤柱一侧和水平线所成夹角大约为49°,由于平均埋深为200 m,根据几何关系可以算出地表下沉39 mm的位置到煤层采空区边界的水平距离是174 m,到采空区中心的水平距离是294 m,小于工业广场到5-210工作面的水平距离,且工业广场地表沉降值在39 mm之内,所以,店坪煤矿5-210工作面布置是合理的,不影响工业广场的使用。
根据煤矿地质条件和生产需求,5-210工作面布置在830水平南翼,为二采区5号煤层第7个回采工作面,西与830南翼运输巷相通,东至井田边界,南北均为实体煤,上部为3-302、3-304采空区,5号煤层与3号煤层层间距55~71 m。
结合工作面上方地表实际情况,布置18个监测点,数据采集使用上海华测GPS接收机(标称水平精度5±1 ppm ,垂直精度10±2 ppm ,进行RTK实时动态观测。在距工作面影响范围外1~2 km的开阔地带布设2个控制点,用于每次日常测量时架设基站、重置当地坐标及检核。截取地表等高线图,比例尺为1∶750,如图5所示。各测点坐标见表2。
图5 5-210工作面地表岩移测点布置
Fig.5 Layout of measuring points for surface rock movement in No.5-210 working face
2016年6月7日,工作面回采位置接近测试区域,开始监测,2016年12月20日,工作面回采结束;持续对地表沉降变形监测到2017年4月3日,地表测点位移不再发生沉降变形,测试结束。在监测过程中,10号测点位于农田中,未能测得数据,12号测点于2016年10月17日测试时,已经随地表裂缝沉降塌毁,后期无数据。其中,测点1最大沉降量约为0.813 m,测点2最大下沉量为1.032 m。测点3于8月9日后趋于稳定,最大下沉量1.951 m,测点4最大下沉量2.690 m。测点5-8受采动影响较小,测点5下沉量为17 mm,测点6下沉量为19 mm,测点7下沉量为23 mm,测点8下沉量为28 mm。测点9在工作面回采过测点200 m时趋于第一稳定阶段,下沉量为0.2 m,待工作面回采600 m后,趋于最终稳定,此时最大下沉量0.243 m。测点11和测点12两个测点均位于工作面边缘,地表下沉的稳定时间几乎同步,测点12的最终下沉量为0.893 6 m(后期无数据),测点11的最终下沉量为1.462 m。测点13和14的下沉量为10~15 mm。测点15、测点16、测点17和测点18的最终下沉量分别为1.287、1.132、1.367、1.415 m。
表2 地表测点坐标
Table 2 Surface coordinates
测点编号X/mY/m埋深/m14 172 644.619 19 509 241.674 1 167.154 24 172 580.503 19 509 235.623 1 170.793 34 172 522.206 19 509 231.991 1 172.029 44 172 472.163 19 509 210.428 1 178.606 54 172 883.639 19 509 101.848 1 122.005 64 172 830.057 19 509 064.610 1 131.076 74 172 773.701 19 509 056.877 1 153.525 84 172 734.645 19 509 055.194 1 163.832 94 172 671.324 19 509 065.321 1 178.640 104 172 503.282 19 509 068.611 1 196.891114 172 524.321 19 508 985.423 1 198.122 124 172 606.859 19 508 982.142 1 179.248 134 172 867.668 19 509 056.169 1 113.982 144 172 886.065 19 509 033.708 1 086.624 154 172 372.216 19 508 692.998 1 093.322 164 172 444.982 19 508 661.130 1 088.462 174 172 533.772 19 508 695.925 1 091.232 184 172 525.571 19 508 751.805 1 107.895
利用FLAC3D构建地表沉陷模型,根据比例尺1∶750,将其模拟结果呈现在5-210工作面等高线图上,如图6所示。
图6 地表沉陷范围
Fig.6 Surface subsidence area
测点1—测点4在影响范围内,其中测点4是5-210工作面下沉最大值点,测点5-8在工业广场附近,影响范围之外,测点1—测点8的下沉曲线图7和图8所示。
图7 测点1—测点4的下沉曲线
Fig.7 Subsidence curves of measuring points 1—4
图8 测点5—测点8的下沉曲线
Fig.8 Subsidence curve of measuring points 5—8
根据地表沉陷范围图和测点下沉曲线图可知,测点5—测点8在店坪煤矿工业广场上,下沉量稳定且在合理范围之内,不受5-210工作面开采的影响,工作面5-210的布置是合理的。
1)利用模拟与计算得出地表沉陷影响范围及工作面布置的合理边界,将5-210工作面布置在830水平南翼,为二采区5号煤层第7个回采工作面,西与830南翼运输巷相通,东至井田边界,南北均为实体煤,上部为3-302、3-304采空区。
2)通过建立数学模型,对店坪煤矿5-210工作面的地表沉陷进行了理论计算,计算结果为2.604 m;采用FLAC3D数值模拟方法,数值模拟结果为2.78 m,现场实测地表最大下沉量是2.69 m,理论计算、数值模拟与实测结果具有良好的耦合性。
3)基于理论计算和数值模拟相结合的方法,设计了5-521工作面合理位置和空间布局,取得了良好的实践效果,为类似地质条件的采空区地表沉陷预测提供了可靠的依据。
[1] 谢东海,冯 涛,袁 坚,等.采矿方法与地表沉陷预测[J].采矿与安全工程学报,2007,24(4):470-472.
XIE Donghai,FENG Tao,YUAN Jian,et al.Prediction of surface subsidence caused by underground mining methods[J].Journal of Mining and Safety Engineering,2007,24(4):470-472.
[2] 张吉雄,谬协兴,郭广礼,等.矸石(固体废物)直接充填采煤技术发展现状[J].采矿与安全工程学报,2009,26(4):395-401.
ZHANG Jixiong,MIU Xiexing,GUO Guangli,et al.Development status of backfilling technology using raw waste in coal mining[J].Journal of Mining and Safety Engineering,2009,26(4):395-401.
[3] 徐法奎.我国煤矿充填开采现状及发展前景[J].煤矿开采,2012,17(4):6-7.
XU Fakui.Current status of stowing mining and Its development prospect in China[J].Coal Mining,2012,17(4):6-7.
[4] 王金喜,李彦恒,孙利辉.采空区地表移动变形预计研究[J].矿业安全与环保,2013,40(5):5-11.
WANG Jinxi,LI Yanheng,SUN Lihui.Predictive study on surface movement and deformation in mined-out Area[J].Mining Safety and Environmental Protection,2013,40(5):5-11.
[5] 文志杰,景所林,宋振骐,等.采场空间结构模型及相关动力灾害控制研究[J].煤炭科学技术,2019,47(1):52-61.
WEN Zhijie,JING Suolin,SONG Zhenqi,et al.Study on coal face spatial structure model and control related dynamic disasters[J].Coal Science and Technology,2019,47(1):52-61.
[6] 庞 会,徐良骥.基于FLAC3D的矿区地表沉陷预测方法研究[J].煤炭技术,2018,37(1):204-206.
PANG Hui,XU Liangji.Research on prediction method of surface subsidence in mining area based on FLAC3D[J].Coal technology,2018,37(1):204-206.
[7] 潘红宇,赵云红,张卫东,等.基于Adaboost的改进BP神经网络地表沉陷预测[J].煤炭科学技术,2019,47(2):161-167.
PANG Hongyu,ZHAO Yunhong,ZHANG Weidong,et al.Prediction of surface subsidence with improved BP neural network based on Adaboost[J].Coal Science and Technology,2019,47(2):161-167.
[8] 徐乃忠,高 超,倪向忠,等.浅埋深特厚煤层综放开采地表裂隙发育规律研究[J].煤炭科学技术,2015,43(12):124-128.
XU Naizhong,GAO Chao,NI Xiangzhong,et al.Study on surface cracks law of fully mechanized top coal caving mining in shallow buried depth and extra thick seam[J].Coal Science and Technology,2015,43(12):124-128.
[9] 王正忠,孙祥柯,毕玉成,等.依兰矿区古近系软岩地层地表沉陷与覆岩破坏规律研究[J].煤炭科学技术,2018,46(S2):206-209.
WANG Zhengzhong,SUN Xiangke,BI Yucheng,et al.Study on surface subsidence and overburden failure law of paleogene soft rock in Yilan Mining Area[J].Coal Science and Technology,2018,46(S2):206-209.
[10] 贺桂成,丁德馨,刘 永,等.衡山石膏矿老采区地表沉陷的ANFIS预测[J].采矿与安全工程学报,2012,29(6):877-881.
HE Guicheng,DING Dexin,LIU Yong,et al.ANFIS prediction of the surface subsidence of the old goaf of the gypsum mine in Hengshan[J].Journal of Mining and Safety Engineering,2012,29(6):877-881.
[11] 崔希民,邓喀中.煤矿开采沉陷预计理论与方法研究评述[J].煤炭科学技术,2017,45(1):160-169.
CUI Ximin,DENG Kazhong.Research review of predicting theory and method for coal mining subsidence[J].Coal Science and Technology,2017,45(1):160-169.
[12] 邹友峰,邓喀中,马伟民.矿山开采沉陷工程[M].徐州:中国矿业大学出版社,2003.
[13] 陈绍杰,朱旺喜,李 军.2004-2013年开采沉陷类国家自然科学基金项目分析[J].山东科技大学学报,2014,33(6):58-62.
CHEN Shaojie,ZHU Wangxi,LI Jun.Analysis of mining subsidence projects funded by national natural science foundation from 2004 to 2003[J].Journal of Shandong University of Science and Technology,2014,33(6):58-62.
[14] 王 宁,吴 侃,刘 锦,等.基于Boltzmann函数的开采沉陷预测模型[J].煤炭学报,2013,38(8):1353-1356.
WANG Ni,WU Kan,LIU Jin,et al.Model for mining subsidence prediction based on Boltzmann function[J].Journal of China Coal Society,2013,38(8):1353-1356.
[15] 王金安,李大钟,马海涛.采空区矿柱-顶板体系流变力学模型研究[J].岩石力学与工程学报,2010,29(3):578-582.
WANG Jinan,LI Dazhong,MA Haitao.Study of rheological mechanical model of pillar-roof system in mined-out area[J].Journal of Rock Mechanics and Engineering,2010,29(3):578-582.
[16] 杨治林.煤层地下开采地表沉陷预测的边值方法[J].岩土力学,2010,31(S1):232-236.
YANG Zhilin.Prediction of surface subsidence in underground mining seam based on the boundary value method[J].Geotechnical Mechanics,2010,31(S1):232-236.
[17] 郝延棉,昊立新,戴华阳.用弹性板理论建立地表沉陷预计模型[J].岩石力学与工程学报,2006,25(S1):2958-2962.
HAO Yanmian,HAO Lixin,DAI Huayang.Establishing a ground settlement prediction model with elastic slab theory[J].Journal of Rock Mechanics and Engineering,2006,25(S1):2958-2962.
[18] 赵伟强,王子升,孙 珞.村庄下充填开采地表沉陷的数值模拟[J].北京工业职业技术学院学报,2018,17(1):29-33.
ZHAO Weiqiang,WANG Zisheng,SUN Luo.Numerical simulation of surface subsidence caused by filling mining under villages[J].Journal of Beijing Polytechnic,2018,17(1):29-33.
[19] 郭玉芳,孟凡迪,陈俊杰.厚松散层开采条件下地表沉陷数值模拟分析[J].煤炭工程,2014,46(6):103-105.
GUO Yufang,MENG Fandi,CHEN Junjie.Numerical simulation analysis on surface subsidence mining under thick alluvial[J].Coal Engineering,2014,46(6):103-105.
[20] 邓小龙,李丽慧,谭玉芳.矿山地下开采诱发地表沉陷的数值模拟[J].煤矿安全,2018,49(7):188-192.
DENG Xiaolong,LI Lihui,TAN Yufang.Numerical simulation of surface subsidence induced by underground mining[J].Safety in Coal Mines,2018,49(7):188-192.