基于FLAC3D 的临涣矿马头门优化方案数值模拟分析

李承涛 吴梦宇

（安徽理工大学土木建筑学院，安徽 淮南 232001）

0 引言

马头门是煤矿井筒与煤矿运输巷道的连接部分，是矿井的重要组成部分，包括了马头口上井筒段、下井筒段及马头门两侧的连接大铜室，因为其位于深部地层，故其受力情况多变且复杂[1]。马头门具有断面大、服务年限长、特别容易受到应力集中等特点[2]，因此在马头门与井壁交界部位就容易出现较严重的应力集中现象，造成马门头发生破坏[3]。马头门发生破坏的原因有很多，除了马头门位置不佳、断面大、设计不合理、支护方案效果不好、施工过程质量问题等原因外，还包括马头门所处的地层特征，由于马门头多位于软弱岩层中，加上马头门支护施工难度大，在地应力特别是在深部地层高地应力的作用下非常容易发生破坏［4-5］。所以，研究马头门支护结构的稳定性对矿井建设的安全非常有必要，能够保障矿工的生命安全和减少煤炭资源的损失，这不仅会使我国在煤炭工程领域安全生产方面拥有更加先进的技术，也可以促使马头门在其他领域得到更广阔的应用前景。

目前，对于马头门的研究成果交多，尤其是采用数值模拟的研究方法具有成本小、速度快的优势，逐渐得到应用，其中应用FLAC3D软件研究马头门的方法效果较好[6]。姚直书[7]等利用FLAC3D技术对开挖和支护全过程5个阶段进行了模拟，探究了马头门在开挖和掘进时围岩受所变化的特性，展现了开挖和掘进过程对马头门的影响；田柯[8]通过FLAC3D软件对马头门容易产生应力集中和围岩容易破坏的区域进行了定性分析，并且对应力集中系数和深度做出了定量的确定。但利用FLAC3D软件进行马门头优化设计的研究仍然不足。

基于此，以临涣矿中央风井中马头门施工支护过程为例，利用FLAC3D软件开展马门头支护方案优化设计的数值分析，以供相关技术人员参考。

1 工程模型建立

1.1 FLAC3D软件介绍

本研究采用FLAC3D软件进行数值模拟分析。FLAC3D软件是一款用于仿真计算的软件，能够对土质、岩石等多种不同性质材料的受力特性以及对具有塑性、流动特性的材料进行三维的数值分析。该软件在三维网格中通过调整多面体单元从而实现拟合实际的结构。该软件可采用线性或非线性的单元材料，当外力作用在材料上时，发生屈服流动的材料，其三维网格能够发生相应的变形和移动，这就是大变形模式。该软件的技术特点为：显式拉格朗日算法以及混合-离散分区技术，正是采用了这些技术，FLAC3D才能非常精准地对材料的塑性破坏和流动进行模拟。正因为它不需要刚度矩阵，所以它在极少的内存下就可以求解极多的单元结构。

FLAC3D软件拥有以下3种计算优势：

（1）混合离散法：当模拟的物体是塑性破坏或者塑性流动时，该方法相比采用“离散集成法”的有限元法更佳；

（2）动态运动方程：即使某个模拟的系统状态为静态，数值模拟过程中仍然使用动态的运动方程，故在模拟物体具有不稳定的状态和过程时，FLAC3D软件不存在数值上的障碍；

（3）“显式解”方案：“显式解”方案在线性本构关系与非线性本构关系的应力-应变求解上所花费的时间是相同的，与之相对的是隐式求解方案，该方案对于非线性问题处理的工作量比“显式解”方案多得多。因为它不需要存储刚度矩阵，这就意味着，在极少的内存下就可以求解极多的单元结构，并且在处理大变形和处理小变形的模拟上所消耗的时间相差无几。

1.2 工程概况

临涣矿跨度为13km，宽度约为4～5km，面积49.6624km2。该矿开采的模式为：使用立井、采用分区通风以及集中出煤的方式，矿建工程包括中央风井、-650m回风车场及回风石门、-450m回风联巷。临涣矿中央风井中马头门位于约-450m的位置，为单侧马头门，井筒净直径为6.5m，马头门及巷道净半径为2.6m。

地面与马头门之间包含3层岩层，自上至下依次为泥岩、粉砂岩及泥岩，厚度分别为20.04m、3.27m以及16.69m，岩层物理与力学参数如表1所示，马头门及井筒衬砌物理与力学参数如表2所示。

表1 岩层物理与力学参数

表2 衬砌物理与力学参数

1.3 模型建立

根据实际工程深立井、马头门及其连接巷道的平面布置、断面尺寸和水平位置，考虑边界效应的影响范围，使用Rhino+Griddle软件建立模型网格，通过将网格导入FLAC3D建立三维有限元计算模型。模型尺寸为长50m、宽40m、高50m，单元总量为1509175。在模型侧面及底面都施加位移约束，在模型顶面施加11.2MPa的竖向应力，如图1所示。

图1 水平马头门数值模型

对于所研究区域，中央风井井筒断面为圆形，巷道及马头门硐室断面为直墙半圆拱形。初期支护时，马头门及井筒均采用锚、喷、网的支护形式，马头门及井筒的锚杆长度为2.5m，间距、排距为0.8m×0.8m，预应力为60kN；马头门处锚索长度为6.2m，间距、排距也为0.8m×0.8m，预应力为120kN，井筒及马头门衬砌浇筑采用C50混凝土。在数值计算模型中，巷道、回风联巷及其连接处共分4个步骤进行模拟，具体描述如下：

（1）马头门上部井筒掘进开挖，设置锚杆等参数，模拟让压过程；

（2）马头门上部井筒二次衬砌，计算平衡；

（3）马头门掘进开挖，设置锚杆等参数，模拟让压过程；

（4）马头门二次衬砌，计算平衡。

2 马头门初期支护方案数值分析

模型建立完成后，对井筒进行开挖。首先对上井筒段进行开挖，开挖完毕后对其进行支护；然后对马头门进行单侧开挖，开挖完毕后进行二次衬砌。利用FLAC3D软件得到其竖向位移、竖向应力以及最大主应力，将云图导入Tecplot软件进行分析，得出位移及应力云图如图2~图4所示。由图2可知，马头门开挖后进行了锚网喷支护，距离巷道、马头门空间越近，对应的岩体位移越大,此时巷道顶板最大下沉降位移超2.5cm，底板底鼓位移已超4.5cm。根据图3和图4可知，竖向应力最大达到了38MPa，最大主应力达到了110MPa，支护失效风险极大。因此，需要对支护方案进行优化。

3 马头门优化支护方案数值分析

3.1 优化方案一数值分析

考虑到有限差分软件模拟单元需连续，导致围岩应力释放产生的变形将作用到衬砌结构上，这与实际情况存在差异。方案一为在井筒及马头门衬砌施加前，将围岩应力释放产生的位移清零，其他支护参数都保持不变，导出位移及应力云图见图5～图7。

对围岩应力释放产生的位移清零进行计算后，对比图5和图2不难发现：优化方案一处马头门处衬砌整体下沉，其中地板下沉约0.5cm，顶板下沉位移超过0.7cm，最大位移已经降低了许多，说明方案一的优化设计对马头门的稳定性有很大的提高；而通过对图3与图6，图4与图7的对比发现，竖向应力和最大主应力在优化前后应力值大小未发生较大改变，说明优化方案一并不能改变竖向应力的大小。

图2 z向位移（单位：cm）

图3 z向应力（单位：MPa）

图4 最大主应力（单位：MPa）

3.2 优化方案二数值分析

优化方案一的模拟结果证实，现有的措施基本能满足支护要求。为了安全考虑，优化方案而在优化方案一的基础上对马头门15m范围内的岩体进行注浆加固，用以改善岩性，依旧采用Tecplot软件得到的位移云图和应力云图如图8～图10所示。根据图8可以看出，相较于图5中马头门整体位移变化不大，说明注浆加固对马头门附近围岩的位移影响不大；对比图6和图9，可以清楚地看出竖向应力也未发生明显变化，只是从原来的38MPa降低为34MPa，对马头门整体影响不大；对比图7和图10能明显看到最大主应力的巨大变化，图7中最大主应力由110MPa变为图10中约40MPa，说明马头门15m范围内的岩体进行注浆加固对降低围岩最大主应力有很好的效果，综合未优化、优化方案一、优化方案二这三个支护方案来看，虽然优化方案一已经基本满足马头门支护的效果，但优化方案二的支护效果更佳、结构稳定性更好。

图5 优化方案一z向位移（单位：cm）

图6 优化方案一z向应力（单位：MPa）

图7 优化方案一最大主应力（单位：MPa）

图8 优化方案二z向位移（单位：cm）

图9 优化方案二z向应力（单位：MPa）

图10 优化方案二最大主应力（单位：MPa）

4 结束语

通过FLAC3D软件建立了临涣矿马门头数值计算模型，分析了单侧马头门采用锚网喷支护以及二次衬砌下的位移场和应力场，由于产生的位移较大，马头门支护失效风险较大，故需要进行优化方案设计。针对初期设计，设计了两种优化方案，一种为在井筒及马头门衬砌施加前，将围岩应力释放产生的位移清零，其他支护参数都保持不变；另一种为在优化方案一的基础上对马头门15m范围内的岩体进行注浆加固，并对比优化前后的位移场和应力场。通过未优化与优化方案一对比可知，马头门处位移大幅下降，竖向应力和最大主应力未发生变化；通过优化方案一与优化方案二的对比可知，马头门处位移、竖向应力未发生明显变化，但最大主应力发生明显变化，故本模拟采用优化方案二的支护效果更佳、结构稳定性更好。