颗粒流
-
2023年3月17日发(作者：自荐书模板)基于颗粒流理论的煤与瓦斯突出数值模拟研究
王锐;修毓;王刚;武猛猛;杜文州
【摘要】针对煤与瓦斯突出机理的复杂性,基于颗粒流离散元理论,利用
PFC3D(particleflowcode)软件对煤体进行340组单轴压缩及巴西劈裂试验,建立
了煤体宏观参数与颗粒流数值模拟所需细观参数之间的关系.基于寺河矿煤样单轴
压缩和抗拉物理实验的实测宏观参数,推算出模拟煤与瓦斯突出所需的细观参数.在
PFC3D流固耦合理论的基础上,利用PFC模拟软件进行了煤与瓦斯突出的数值模拟,
再现了煤体的堆积过程.研究表明:在突出准备阶段,煤壁深部瓦斯压力较高·突出孔洞
附近瓦斯压力梯度较高;发动及发展阶段,平衡状态被破坏,形成口小腔大的突出孔洞,
破坏范围由巷道高度的中间位置向顶底板处扩展;在突出过程中,工作面前方煤体靠
近煤壁的暴露面处的最大主应力向量逐渐偏离垂直方向.且应力向量连线形成一个
凸向煤壁深部的圆弧.
【期刊名称】《山东科技大学学报（自然科学版）》
【年(卷),期】2016(035)004
【总页数】10页(P52-61)
【关键词】颗粒流;煤与瓦斯突出;PFC3D;宏-细观参数
【作者】王锐;修毓;王刚;武猛猛;杜文州
【作者单位】山东科技大学矿业与安全工程学院,山东青岛266590;青岛市黄岛区
城市建设局,山东青岛266555;山东科技大学矿业与安全工程学院,山东青岛
266590;山东科技大学矿山灾害预防控制省部共建国家重点实验室培育基地,山东青
岛266590;山东科技大学矿业与安全工程学院,山东青岛266590;山东科技大学矿
业与安全工程学院,山东青岛266590
【正文语种】中文
【中图分类】TD713
煤与瓦斯突出是地应力和瓦斯压力共同作用使煤体发生破坏变形的一种复杂动力现
象。我国是煤与瓦斯突出灾害最严重的国家之一，在过去的150年里，超过1/3
的突出事故发生在我国[1]，严重威胁我国矿井的安全生产。因此研究煤与瓦斯突
出机理，预防突出事故发生是实现煤矿行业可持续发展的重中之重。目前，煤与瓦
斯突出机理尚不明确，世界各国的专家学者分别从现场调查、相似材料实验、历史
资料统计分析、数值模拟等角度进行分析，而数值模拟是研究煤与瓦斯突出机理的
主要途径，也是最为简便的方法之一。
安丰华等[2]曾利用COMSOLMultiphysics软件对煤与瓦斯突出进行数值模拟，
在研究孔隙度和渗透率受采掘活动影响的基础上，分析吸附瓦斯在突出灾变过程中
对煤体性能的影响；为了解决煤与瓦斯突出机制尚不明确、影响因素复杂的问题，
王刚等[3]对瓦斯突出的三维结构模型进行重构，引入围岩弹性潜能，构建了突出
能量条件模型和突出强度预测模型，并结合突出案例对模型进行了验证，对突出影
响因素的灵敏度进行了分析。而后王刚等[4]又将FLAC3D(fastlagrangian
analysisofcontinua)和COMET3两个商用软件相耦合：利用FLAC3D进行应力
分析，COMET3模拟煤的瓦斯吸附与解吸、扩散，研制并完善了含量法预测瓦斯
突出的计算机软件程序；黄达等[5]利用PFC颗粒流软件探讨了初始单轴静态压缩
的细观损伤程度对单轴动态压缩下单裂隙岩样力学性质的影响规律，对其能量耗散
机制进行了总结；崔峰等[6]借助PFC软件，构建了等效于煤体耦合致裂后强度的
离散元模型，以此为基础研究了耦合致裂后离散态颗粒间的铰接关系及其对垮放程
度的影响，注水压力和煤体垮放能力之间的映射关系，并结合工程实践进行了验证；
王莉等[7]利用PFC软件构建了细观尺度上巴西劈裂试验声发射的细观模拟方法，
并将试验和计算结果进行对比分析，验证了该方法的合理性；为了解揭煤区域内发
生煤与瓦斯突出可能性，高瑞元等[8]运用PFC软件对煤层进行数值模拟，并对模
拟结果和瓦斯压力进行加权计算处理，使模拟结果更真实地反映现场实际情况。胡
莹莹等[9]曾在实验室内模拟煤与瓦斯突出现象，研究突出过程中地应力和瓦斯压
力的变化情况；徐涛等[10]基于RFPA2D(rockfailureprocessanalysis)流固耦合
理论研究煤与瓦斯突出机理，在球壳失稳假说[11]的基础上分析了瓦斯压力和地应
力对煤体的破坏变形。尽管煤与瓦斯突出机理的研究取得了多方面的进展，但突出
灾变过程中煤体破坏是由内部开裂的微裂纹等细观缺陷逐渐发展的结果[12]，宏-
细观参数选取的准确性很大程度上决定了煤与瓦斯突出数值模拟的精度。
Potyondy等[13]基于颗粒离散元理论研究了花岗岩的细观参数；Yang等[14]以
岩石样本作为研究对象，研究了平行粘结模型中颗粒细观参数和样本弹性参数(弹
性模量、泊松比)以及单轴抗压强度的关系；Hsieh等[15]研究了砂岩变形机制，
分析了细观组分含量与宏观力学特性之间的关系；Backstrom等[16]利用颗粒流
方法对花岗岩细观模拟结果与宏观应力-应变曲线的关系进行了研究。针对煤颗粒
细观参数与其宏观力学参数相关性的研究则鲜为人知，因此如何有效选取煤体宏-
细观参数是本研究的重点之一。
本文首先设计多组单轴压缩及巴西劈裂数值模拟试验，通过编制颗粒离散元
PFC3D程序并对其进行模拟计算，进而分析煤体宏观、细观力学参数之间的关系。
然后以寺河矿原始应力条件下取得的煤样为例，进行单轴压缩及巴西劈裂物理试验，
以期获得煤体宏观力学参数，根据所建立的煤体宏观、细观力学参数关系推算相对
应的细观参数，在此基础上利用PFC3D程序对煤与瓦斯突出进行数值模拟。
1.1颗粒流数值模拟基本原理
在PFC3D计算循环过程中，颗粒间的作用力主要遵循牛顿第二定律与力-位移定
律(如图1所示)，而颗粒的接触本构模型包含刚度模型、滑动模型和粘结模型三部
分[17-19]。
利用PFC3D进行单轴压缩和巴西劈裂数值模拟试验，在计算颗粒与墙体接触的平
均刚度基础上，通过设定顶部或底部墙体的速度来实现模型的加载过程，通过伺服
机制控制墙体的运动对试件进行压缩试验，将突出过程中瓦斯作用于煤颗粒的力等
效为有效应力并考虑压力梯度的影响[20]。
1.2单轴压缩及巴西圆盘劈裂数值模拟试验
根据物理力学试验提供的资料，结合相关的岩石力学性质进行取值，相关参数如表
1所示，通过PFC3D软件生成的试样如图2所示。
设计多组单轴压缩和巴西圆盘劈裂试验，研究杨氏模量、法向与切向刚度之比对单
轴抗压强度、泊松比、弹性模量以及单轴抗拉强度的影响。试样结果如图3所示，
浅色部分代表模拟过程生成的裂隙，深色部分在颗粒及裂隙剖面图中代表煤体颗粒，
在接触力及裂隙图中代表接触力。
1.2.1弹性模量的回归分析
由图4可以发现弹性模量与杨氏模量存在线性关系，与刚度比存在负对数相关关
系。基于上述结论，对各方案结果进行回归分析，得出弹性模量与细观参数关系式
(式(1))，其相关系数为R=0.9895，两者拟合效果较好。回归分析结果与模拟数
据对比如图5所示。
1.2.2泊松比的回归分析
由图6和图7可以发现泊松比与杨氏模量相关性较小，刚度比是影响泊松比的主
要参数。对各方案结果进行回归分析如图8所示，其函数关系如式(2)所示，相关
系数为0.9897。
1.2.3单轴抗压强度的回归分析
为了研究切向连接强度与法向连接强度不同时对单轴抗压强度的影响，设置杨氏模
量与刚度比均为定值，切向连接强度与法向连接强度比值取值0.1～5.5，杨氏模
量取值为10MPa，模拟结果与拟合曲线如图9所示。当≤2时，单轴抗压强度受
平行连接切向连接强度和法向连接强度共同影响；当>2时，平行连接切向连接强
度对σc的影响较弱，此时平行连接法向连接强度与σc成线性关系。对该方案结
果数据进行回归分析，得到单轴抗压强度与细观参数函数关系如式(3)所示，其相
关系数为0.9884，拟合效果较好。
1.2.4抗拉强度的回归分析
由图10可以得到与单轴抗压强度相类似的抗拉宏细观参数相关性结论，对该方案
结果数据进行回归分析，得到抗拉强度与细观参数函数关系如式(4)所示，其相关
系数R=0.9597，拟合效果较好。
1.3煤体细观参数的确定
在宏细观参数相关性的基础上，利用MTS815.03岩石伺服试验系统对实际煤样试
件进行单轴压缩和抗拉物理试验得出宏观参数的前提下，确定煤体的细观参数。试
验所用煤样取自山西省晋城无烟煤集团矿业有限公司寺河矿，试件详细信息如表2
所示，其中AY-1与BL-1、AY-2与BL-2分别取自矿井相同地点的同一钻孔。
单轴压缩试验及抗拉试验前后煤体试件见图11。对寺河矿煤样进行单轴压缩和抗
拉试验后，得出实际煤体的弹性模量、泊松比、单轴抗压强度和抗拉强度等宏观物
理参数，利用式(1)～式(4)在宏观物理参数的基础上推算得到相应的煤体细观参数，
如表3所示。
2.1数值模型
基于煤与瓦斯突出机理利用PFC3D程序进行颗粒流模拟，为了便于计算与分析将
煤颗粒理想化为圆形颗粒。综合考虑计算量和计算的精度，模型尺寸为8cm×15
cm×8cm，为了更好的模拟煤体颗粒分布的不均匀性，颗粒半径分布从
Rmin=0.35mm到Rmax=0.50mm采用了均匀分布。针对煤颗粒由于粘结强度
不够而四处逃逸的问题，在模型周围加入6面墙体，将模型的上下、前后边界设
置为不透气边界。
在地应力为20MPa，瓦斯压力梯度(左(3MPa)到右(0MPa))作用下的灾变过程，
在设置模型孔隙率为0.077的前提下，根据表3颗粒的宏-细观参数建立突出模拟
试样模型及网格划分如图12所示，模型包含颗粒数量8000个，网格划分
8×15×8个单元，黑色颗粒代表了煤颗粒，黑色线条表示所划分的网格。
2.2模拟结果
2.2.1突出前后过程及瓦斯压力分析
由图13(a)知，受到地应力和瓦斯压力的共同作用，使得部分煤体破碎并被带出孔
洞，堆积在孔洞外附近且造成抛出通道的面积逐渐减小，使得孔洞内的煤体破碎过
程缓慢，瓦斯流动的阻力增大且流速降低，此时瓦斯压力会急剧升高而突破阻止瓦
斯流动的障碍，继续将煤体破碎抛出并进一步堵塞通道，该过程反复进行，直至煤
颗粒将通道堵塞使得瓦斯压力不足以突破，达到基本平衡突出终止。
由图13(b)知，煤体内部最大主应力一般接近垂直方向且越靠近抛出通道时，这些
最大主应力越偏离垂直方向，其连线看起来更像一个凸向煤壁深部的圆弧。由于越
靠近抛出通道，在煤层上下部位的水平剪应力越大，使得最大主应力不能再用垂直
应力代替，所得结果与文献[21]中所提到的观点相一致。
图14为煤与瓦斯突出过程中，4×10×4单元(4#测点)、4×8×4单元(3#测点)、
4×4×4单元(2#测点)以及4×2×4单元(1#测点)处瓦斯压力变化情况(瓦斯压力设
置为3MPa)。图15显示了煤与瓦斯突出过程中，各测点处瓦斯压力变化情况(瓦
斯压力设置为5MPa)。
从图14～15可以看出，4#测点瓦斯压力测点随着充气过程的不断充入，瓦斯压
力最终稳定在2.4MPa左右。4#测点变化是最大的，突出完成后瓦斯压力稳定在
0.6MPa。原因是由于该点处在巷道突出口的内部，相当于内部瓦斯向外涌出，也
与现场实际情况相吻合。3#瓦斯压力测点随着充气过程的不断充入，瓦斯压力最
终稳定在1.4MPa左右。3#测点变化趋势和4#相似，突出完成后瓦斯压力稳定
在0.4MPa。但最终稳定的数值要明显低于1#测点，主要是由于该测点距离突出
口相对较近，突出发展的强度较低。1#和2#测点由于距离开采工作面较近，瓦斯
压力测点随着充气过程的不断充入，瓦斯压力最终稳定值略大于巷道大气压。突出
完成后瓦斯压力略微增大，主要是由于深部煤层及巷帮侧瓦斯向外涌出，巷道中瓦
斯压力及含量增大，使得1#和2#出现增长趋势。当瓦斯压力增大至5MPa时，
各压力曲线变化趋势与3MPa时相同。
因此：越靠近煤壁深部瓦斯压力越大，煤壁表面的瓦斯压力接近于巷道的大气压力；
突出准备阶段瓦斯压力基本保持不变。在发动以及发展阶段时，瓦斯压力在很短的
时间内急剧下降，导致煤体破坏及抛出。突出后煤壁深部的瓦斯压力略高于煤壁表
面的巷道压力。
2.2.2颗粒运动轨迹图
数值计算过程中3个监测颗粒的位置如图16所示。对颗粒轨迹进行分析，位于突
出口上部的5和下部的7号颗粒最先受到影响，由于实验初期瓦斯压力较大，并
且位于突出口附近，所以颗粒间的阻碍作用较小，两个颗粒沿试验槽运动，分别相
对于突出口(x=0.83m)产生了0.2和0.16m的水平位移；随着突出的发展，6号
颗粒被剥蚀出来，沿试验槽分别产生了0.2m的水平位移，由运动轨迹可以看出
试验槽上颗粒的淤积对其水平运动产生限制，运动轨迹也表现出了突出后中后部分
颗粒的堆积情况；所有监测的颗粒运动轨迹均向突出口移动且水平方向以及垂直方
向的移动位移较大。
2.2.3喷出口速度图
喷出口速度对研究突出过程有着重要意义。煤与瓦斯突出过程中颗粒抛出速度越大
时，颗粒所受到的力就越大，然而在室内试验中很难对某个颗粒速度进行监测追踪。
数值模拟过程中通过设置fish函数对喷出口处289、387、772、902号颗粒速度
进行监测，监测结果如图17所示。总体上看，颗粒最大抛出速度为30.2、66.2、
37.5和29.5m/s,平均速度为5.3、9.8、8.2和6.7m/s。颗粒刚被抛出时，速度
为整个过程中的最大值，且随着向前的移动，能量逐渐被耗散，速度较低。
在物理实验得到宏观参数及模拟实验得到宏-细观参数关系的基础上，利用PFC模
拟软件进行了煤与瓦斯突出的数值模拟，再现了煤体的堆积过程，得到如下结论：
1)突出的准备阶段，煤体内部最大主应力接近垂直方向且越靠近煤壁时，这些最
大主应力越偏离垂直方向，煤壁深部瓦斯压力较高且基本相等，突出孔洞附近瓦斯
压力梯度较高；在发动以及发展阶段，高压瓦斯对煤体进行破坏、抛出，并在很短
的时间内急剧下降；突出后孔洞周围存在一定范围的塑性区，并再现了突出煤体的
堆积过程。
2)突出发动和发展阶段，突出孔洞附近煤体的极限平衡状态被打破，煤岩体破坏
失稳并被抛出，形成了口小腔大的突出孔洞，随着突出过程的发展，掘进面前方煤
体的破坏范围和突出孔洞增大，破坏范围由巷道高度的中间位置向顶底板处扩展。
3)突出过程中，工作面前方煤体靠近煤壁的暴露面处的最大主应力向量逐渐偏离
垂直方向，且应力向量连线形成类似一个凸向煤壁深部的圆弧。突出口处颗粒最大
抛出速度为30.2、66.2、37.5和29.5m/s，平均速度为5.3、9.8、8.2和6.7
m/s。