露天矿边坡岩体状态变化与多灾害耦合作用关系研究

2023-12-23崔铁军李莎莎

防灾减灾学报 2023年4期

崔铁军，李莎莎

（沈阳理工大学环境与化学工程学院，辽宁沈阳 110159）

0 引言

矿业生产是国民经济的支柱之一，煤矿生产更是重中之重。煤矿生产主要分为露天和井工两种方式，相较于井工方式，露天开采具有明显优势，包括开采受限制少、大型机械多、资源回收率高、生产效率好、成本低、生产条件好、建设速度快等；同时缺点是占用土地多，对地表和环境破坏严重，受气候影响直接等。具体来说，露天开采是开放式的开采方式，形成人造矿坑和边坡。这些矿坑和边坡的原始地质环境可能较为复杂且有水系发育；特别是随着露天矿深度的增加，边坡坡度和范围也将逐渐增加，由于设计和已知地质条件的缺失极易造成失稳塌方，受雨水作用形成泥石流，或是矿坑周围造成大范围沉陷和地表裂缝。无论是哪种灾害都将影响露天矿的安全生产，同时影响周边建筑及人们的生产生活。特别是在闭坑后，由于缺少生产带来的资金支持，对露天矿的维护措施将进一步弱化，这将促使更大范围的矿坑变形和周围地质条件的改变，直到水文和地质条件重新平衡。进一步而言，各种自然灾害和水文地质条件对露天矿矿坑和边坡的破坏将加强，例如水系作用、矿震、煤层自燃等都将改变岩体状态，而岩体重新平衡的过程就是矿坑和边坡的破坏过程，况且这些灾害过程可能同时出现，具有耦合作用，会加快破坏过程，加强破坏程度。多灾害耦合作用将对资源型城市和距露天矿较近的城镇造成严重的地质灾害隐患，因此岩体状态变化与多灾害耦合作用关系是亟待解决的重要科学和工程问题。

对矿业领域多灾害耦合研究的文献较少。明确提出露天矿边坡工程中各类灾害演化过程耦合研究的是王来贵教授，他在2014 年发表的露天矿边坡工程系统演化过程一文中详细介绍了灾害演化耦合关系[1]。包括多灾害耦合叠加模型的区域地震风险评估[2]，深部采动响应与灾害防控[3]，瓦斯突出多物理场参数演化[4]，多因素耦合模型的滑坡预警[5]，多因素作用的极震区泥石流过程分析[6]，多灾害耦合矿井强冲击灾后恢复治理[7]，多场耦合作用下瓦斯与煤自燃协同预防[8]，动力灾害风险精准判识及监控预警[9]，煤与瓦斯突出动力灾害多参数耦合预警[10]，多灾害耦合的贴近突出煤层安全开采[11]，煤矿动力灾害多参数耦合测定[12]，多场耦合煤矿动力灾害模拟研究[13]等。这些研究多集中于多因素对灾害的影响，而少有多灾害之间耦合作用对岩体应力应变的影响研究。而实际情况下露天矿经常出现多灾害连发、相互促进、发生发展的综合灾害演化过程。这种灾害演化过程与各灾害之间的关系、出现时机、发生程度等有着显著的相关性，因而需要确定不同条件下可能出现的灾害及其耦合情况，但目前缺乏相关研究。

对于上述问题，笔者在数学方法上提出了空间故障网络理论[14]，可描述自然灾害演化过程，确定灾害过程经历的事件、影响因素、逻辑关系和演化条件等。但实际中该方法需要大量数据，难以对矿区具体位置进行灾害演化过程分析。结合笔者已有的矿山灾害模拟研究[15-18]，论文提出了使用模拟方法研究露天矿边坡岩体状态变化与多灾害耦合作用关系来替代露天矿边坡自然灾害演化过程的研究。设置了该露天矿在未来50 年内可能出现的灾害及其耦合情况，通过模拟得到渗水、矿震和煤层自燃与岩体塑性、位移和应力之间的关系，为露天矿灾害预测防治提供了依据。

1 露天矿水文地质条件与模型建立

某露天矿东西长为6600 m，南北宽为2200 m，总面积达13.2 km2。其位于城市边缘，与市区非常接近，随着城市发展和露天矿开采活动的进行，该矿北帮边坡的应力应变状态变化对临近城市造成了严重影响。最显著的是地表出现的不均匀沉降，长达百米的地表裂缝及周围建筑产生的裂隙。加之矿坑北帮以北1200～1600 m有较大河流经过，且矿坑位置原为两条支流的河道，使得矿区周边水系发育较好，据测算在不进行排水情况下矿坑将在50 年左右完全注水形成湖泊。研究涉及的矿坑边坡和河流等的位置如图1所示。

图1 露天矿、市区及河流地形图Fig.1 Topographic map of open pit mine，urban area and river

表1 模型内各组土层名称

研究模拟范围为北帮E200 到E800 范围，如图1 中虚线内区域。模型南北长1100 m，东西宽600 m，底部海拔-600 m，最高处海拔70 m。该模型北侧边界距河流平均1300 m。模型中包括花岗片麻岩、玄武岩、白垩系砂岩、绿色泥岩、油母页岩、煤岩等，同时存在多条断层。根据不同岩层的位置和物理力学参数设置模型及断层接触面，模型层理结构、主要研究区域和地震作用位置如图2 所示，编号与岩层对应关系如表1所示，各岩层参数如表2所示。

图2 左图表示模型中不同区域岩层的岩体种类。中图表示模拟过程中岩体状态变化的主要区域，区域1 是大高差边坡，主要由绿色泥岩构成；区域2 是小高差边坡，主要由玄武岩、煤岩和油母页岩构成。右图表示矿震荷载施加位置，即模型的下底面和北侧边界，以及发生自燃的煤层位置和接触面。

岩体模型采用FLAC3D 的摩尔库伦本构模型建立，网格为四面体网格，边长为10～30 m。为了模拟断层设置两处接触面，分别为ZG_008与ZG_007 的接触面1，ZG_007 与ZG_006 的接触面2，接触面参数如表3所示。

表1至表3参数是经过实地测量获得的，有一些参数通过差值法得到。

图2 模拟区域构造图Fig.2 Simulated regional structure

表2 各组土体参数表

表3 接触面参数

2 多灾害耦合作用设置

考虑到在自然状态下矿坑预计将在50 年左右被周围水系注满，同时考虑到坑底距地表约为400m，将坑内水位高度变化50m 作为一个模拟阶段，那么全过程模拟分为9 个阶段进行，分别是-380m（第一次）、-380～-330（第二次）、-330～-280（第三次）、-280～-230（第四次）、-230～-180（第五次）、-180～-130（第六次）、-130～-80（第七次）、-80～-30（第八次）、-30～20（第九次）。模拟过程的9 个阶段如图3 所示。另外根据不同高度的水的渗入速度，可判断9 个阶段分别对应的持续时间，但在这里不是重点不做赘述。

图3 潜水层设置与水位变化过程图Fig.3 Water plan setting and water level change process

考虑到矿区较为常见的灾害形式，设置在渗水过程中伴随着矿震和煤层自燃现象。模拟主体顺序为自燃、渗水和矿震，自燃和矿震的发生在渗水期间是随机发生的。

煤层自燃现象是煤矿开采中经常发生且不可避免的自然灾害，因此首先考虑煤层自燃对边坡的影响。自燃的模拟考虑到实际煤燃烧过程中的物理化学性质变化，突出表现在体系缩小、强度减小及密度减小几方面，因此对煤层中煤的参数进行调整。另考虑煤自燃主要发生在煤层露头部分及一定深度内，结合已有研究成果[17-18]，设置煤自燃范围是从露头延伸入煤层260 m 左右，最深处距地表约178 m，如图2中右图煤层所示。上述9 个阶段内均发生一次自燃模拟，参考表2 对煤层参数进行适当折减，折减系数为0.94。但考虑到煤层实际位置，在被水覆盖后煤层自燃较为困难，因此在第6 次模拟后不再设置煤层自燃。

水的渗透作用对边坡的影响模拟使用潜水面（water plan）实现。water plan 之下代表岩体被水淹没，认为岩体被水完全浸泡是饱和状态岩土体。对应于水位变化的9个阶段，设置9个water plan 来模拟9 次水位上升过程对岩体的影响。由于本构的内部机制，water plan 之下的岩体以饱和形式模拟，那么关键问题是设置适合的water plan。参考该露天矿与北侧最大水系距离，模拟区域与该河流的平均距离为1300 m，因此设置9 个water plan 的焦点在模型北侧水平距离1300 m 处。另外模型南侧边缘到河流水平距离约2400 m，按照实际测量存在3/1000 的水力梯度，因此模拟河流位置应比实际位置下降7.2 m，同时考虑到该河流流域海拔约为39～70 m，及河流下部的渗透情况，模拟位置再下降10 m，最终各water plan 的焦点海拔设为20 m。因此water plan的焦点假设位置在海拔20 m且距模型南侧2400 m，如图3所示。

在露天矿矿坑注水后极易发生矿震和小型地震。在每个模拟阶段内都设置3 次矿震，每次矿震约为三级，对应的最大运动速率约为0.075 m/s。施加的地震荷载为正弦曲线波，且为相同纵波和横波的混合形式，震动周期为0.7 s，最大速率为0.075 m/s，速率随时间的变化而变化。地震正弦波施加的位置为模型最北侧边界平面和模型底部边界平面，如图2 右图中设置的正弦地震波位置。

综上所述，模拟矿坑及边坡在50 年自然状态下，研究区域模拟过程可总结为：根据水位变化分为9个阶段进行模拟，前6阶段模拟过程为自燃、地震和渗水的耦合结果，后3 阶段为矿震和渗水的耦合结果。在水充满当前高度的模拟过程中，随机设置发生自燃和地震，自燃发生1次，矿震发生3次。随机发生通过间隔计算步实现，虽然随机灾害影响的矿坑和边坡状态变化结果不同，但总体状态变化的趋势相同，可以为类似灾害过程耦合提供分析对照。

3 模拟结果与分析

上述多灾害耦合作用模拟过程共分9 个阶段，每阶段中的自燃和矿震存在随机顺序，但总体上先执行自燃，然后水位变化，最后模拟矿震，因此前6阶段有3次模拟平衡，分别为自燃、渗水和矿震，后3阶段只有2次平衡。由于研究需要和篇幅所限，仅研究了9 阶段的最终平衡状态，即矿震结束后的平衡状态。获得岩体状态变化最具特征的量，包括塑性区变化、位移变化和主应力变化。

多灾害耦合作用下模拟区域矿坑及边坡的9阶段塑性区变化过程如图4所示。

图4 塑性区变化图Fig.4 Change of plastic zone

在模拟过程中塑性区表示了岩体变形和受力等超过弹性变形极限，岩体进入塑性的状态，在应力消失后变形依然有残余。对岩体而言塑性区则代表了受拉和受剪破坏的区域，而岩体的抗拉和抗剪强度较低，因而塑性区代表了岩体在模拟过程中可能受破坏的区域。从图4 的塑性区变化过程可知岩体状态变化的特征。边坡顶部存在两条明显的红色区域，是受拉破坏形成的塑性区。这两个位置对应了两条断层（对比图2 右图），因而在岩体运移过程中都是较薄弱的环节，只要下部岩体状态发生变化这两处断层就会出现较大位移形成塑性区，这是9阶段模拟结果的共同特点。

讨论模拟过程中塑性区的特征变化，根据图2 中的区域1（大高差边坡）和区域2（小高差边坡）分别讨论，两个区域中间的坑底未受到破坏，处于原始状态。对区域1，主要存在于边坡自由面，也是决定边坡安全系数的重要区域。从图4中的9次模拟可知该部分岩体一直表现为受剪和受拉剪破坏状态。受拉剪破坏形成的塑性区（下文简称拉剪塑性区）主要集中在边坡上部，而下部是受剪产生的塑性区（简称剪塑性区）。形式上9 次模拟中剪塑性区和拉剪塑性区的总区域变化不大，且拉剪塑性区在剪塑性区之上。随着模拟次数增加，即渗水水位上升，拉剪塑性区逐渐增加，剪塑区逐渐减小，但剪塑性区始终在拉剪塑性区下方。解释这一现象，受剪区域代表了不稳定的滑坡，而从基岩上产生的滑坡必然使之上岩体与基岩之间产生剪切变形，形成剪塑性区。而拉剪塑性区更接近于自由面，这部分岩体自由度更大约束更小，在模拟过程中拉剪塑性逐渐增加。提供这部分生成条件的主要是矿震震动；另一方面由于水位的上升，静水压力也提供了形成上述现象的条件。区域2 主要是由玄武岩、煤岩和油母页岩组成，该区域岩体状态较为多样，包括剪塑性区、拉剪塑性区和拉塑性区。随着模拟次数增加，这些塑性区都在逐渐增加，其特点是剪塑性区仍然在该区域的底部，其上是拉剪塑性区，最上方是拉塑性区。分析形成原因主要是受到煤层自燃和水位影响，受地震影响较小。因为这部分边坡高度低、体积小，并未形成如区域1 的大高差边坡，对震动正弦波作用的响应不大，形成的剪塑性区不大。煤层自燃导致煤体体积缩小，其物理力学性质弱化，导致上覆岩体沉降；沉降岩体产生明显位移形成拉剪塑性区；之上区域由于拉剪塑性区的支撑作用并未产生剪塑性区，而受到侧向岩体的拉力维持稳定，因此区域2 顶部形成了拉塑性区。另一方面，由于水位升高，第七至九次模拟都未考虑煤的自燃，期间区域2 岩体状态变化不大。

多灾害耦合作用下模拟区域矿坑及边坡的9阶段位移变化过程如图5所示。

图5 位移变化图Fig.5 Change of displacement

位移更为直观的表示了边坡变形状态特征。如图5 中9 阶段模拟的位移变化可分为3 个阶段，前三次为第一阶段，特点是water plan 在岩体之下。这时水渗入岩体间隙形成饱和状态，water plan 之上岩体并未受到水的作用仍处于原始状态，下部浸水岩体的骨架系统未变化。区域1 的变形更为明显，但变形量较小，对边坡的安全性影响较小，不构成滑坡破坏条件。随着水位升高，第四至第七次模拟结果的位移量较上一阶段增加了一个数量级。大变形区域集中于被水浸没的部分，且主要变形区域由区域1 转为区域2。这是由于区域2 下部经历了煤层自燃，自燃后煤体受水侵蚀强度进一步降低，体积进一步缩小，造成上覆岩体沉降，沉降将进一步促使周围岩体受拉剪作用形成更大范围的沉降。但由于煤层自燃的深度存在极限并未影响到区域1，因此随后的模拟中区域1并未发生大变形。由于水位上升，煤层自燃难易进行，第八、九次模拟结果差别不大。经过更细致的模拟结果分析，当水位较低时矿震对边坡影响较大；而水位较高浸没边坡后矿震对边坡的影响较小。

多灾害耦合作用下模拟区域矿坑及边坡的9阶段主应力变化过程如图6所示。

图6 主应力变化图Fig.6 Change of principal stress

最大主应力代表了岩体在某一方向上的最大受力情况，在岩体中最大应力位置存储了较大的弹性势能有释放倾向，可能受拉、压或剪应力。同样最大主应力的变化也可分为三个阶段，分别为第一至第四次模拟，这时最大主应力的数量级不变，数值差别不大。因为这时最大主应力来源于岩体本身，水位并未完全浸没底部岩体，模型整体上压应力占主导地位，拉压应力变化较小；模型边坡不同程度的承受拉应力，这是由于震动造成的受拉状态，因为此时水位和煤层自燃都不作用于自由面岩体。第二阶段第五至第七次模拟，这时随着水位升高，由于水渗透进入岩体产生静水压力（浮力），因而该阶段模型中的压应力变化不大，而拉应力出现增加和集中。这同样是水位变化的作用，而随着水位升高，矿震对边坡影响的程度降低；煤层自燃在该阶段也停止，自燃造成了煤层弱化，其上部岩体由于煤层体积流失和强度下降有下沉的趋势，导致这部分岩体受周围拉应力增加。第三阶段第八和九次模拟中，由于水位已到达高位，矿震对最大主应力的影响变得更小；煤层自燃停止对主应力的影响停止，因而这两次模拟结果的变化较小。