黏弹介质槽波超前探测陷落柱三维数值模拟.pdf

Mine Engineering 矿山工程矿山工程,2023,11(3),416-424 Published Online July 2023 in Hans.https:/www.hanspub.org/journal/me https:/doi.org/10.12677/me.2023.113052 文章引用文章引用:何帝秀,于坤,张亚伟.黏弹介质槽波超前探测陷落柱三维数值模拟J.矿山工程,2023,11(3):416-424.DOI:10.12677/me.2023.113052 黏弹介质槽波超前探测陷落柱三维数值模拟黏弹介质槽波超前探测陷落柱三维数值模拟 何帝秀何帝秀，于于 坤，张亚伟坤，张亚伟 安徽理工大学地球与环境学院，安徽 淮南 收稿日期：2023年6月17日；录用日期：2023年7月18日；发布日期：2023年7月25日 摘摘 要要 工作面前方隐伏陷落柱的存在常引发突水等安全事故，给煤矿安全生产与经济收益造成影响。实际地下工作面前方隐伏陷落柱的存在常引发突水等安全事故，给煤矿安全生产与经济收益造成影响。实际地下煤层具有黏弹特性，会使地震波能量产生衰减。为研究实际煤矿中掘进巷道前方存在陷落柱条件下反射煤层具有黏弹特性，会使地震波能量产生衰减。为研究实际煤矿中掘进巷道前方存在陷落柱条件下反射槽波的传播规律，建立了三维黏弹性含陷落柱煤系地质模型，对三维黏弹一阶速度槽波的传播规律，建立了三维黏弹性含陷落柱煤系地质模型，对三维黏弹一阶速度应力方程进行了高应力方程进行了高阶有限差分数值模拟，分析了三分量的波场快照与合成地震记录。结果表明，阶有限差分数值模拟，分析了三分量的波场快照与合成地震记录。结果表明，Y分量的反射分量的反射Love槽波能槽波能量大、波形特征明显；量大、波形特征明显；X和和Z分量的反射分量的反射Rayleigh槽波能量较小、同相轴较模糊。因此，建议在实际煤矿槽波能量较小、同相轴较模糊。因此，建议在实际煤矿开采时利用开采时利用Love槽波超前探测陷落柱。槽波超前探测陷落柱。关键词关键词 煤层煤层，槽波槽波，黏弹介质黏弹介质，超前探测超前探测，波场特征波场特征 Three-Dimensional Numerical Simulation of Collapse Column in Viscoelastic Medium with In-Seam Wave Advance Detection Dixiu He,Kun Yu,Yawei Zhang School of Earth and Environment,Anhui University of Science and Technology,Huainan Anhui Received:Jun.17th,2023;accepted:Jul.18th,2023;published:Jul.25th,2023 Abstract The existence of hidden collapse columns in front of the working face often causes safety accidents such as water inrush,which affects the safety production and economic benefits of coal mines.The actual underground coal seam has viscoelastic properties,which will cause attenuation of seismic wave energy.In order to study the propagation law of reflected in-seam wave under the condition 何帝秀 等 DOI:10.12677/me.2023.113052 417 矿山工程 of collapse column in front of tunneling roadway in actual coal mine,a three-dimensional viscoe-lastic coal measure geological model with collapse column was established.The high-order finite difference numerical simulation of three-dimensional viscoelastic first-order velocity-stress equ-ation was carried out,and the three-component wave field snapshot and synthetic seismic record were analyzed.The results show that the reflected Love in-seam wave of Y component has large energy and obvious waveform characteristics.The reflected Rayleigh in-seam wave energy of the X and Z components is small and the in-phase axis is blurred.Therefore,it is suggested to use Love channel wave to detect collapse column in advance in actual coal mining.Keywords Coal Seam,In-Seam Wave,Viscoelastic Medium,Advance Detection,Wave Field Characteristics Copyright 2023 by author(s)and Hans Publishers Inc.This work is licensed under the Creative Commons Attribution International License(CC BY 4.0).http:/creativecommons.org/licenses/by/4.0/1.引言引言 煤矿开采工作中常遇到一些如断层、陷落柱、采空区等地质体的干扰。在陷落柱比较发育的地区，其是不容忽视的灾害性地质异常体，煤系地层中陷落柱的存在会破坏煤层的连续性，影响煤巷围岩的稳定性，影响巷道的掘进和煤矿的开采，还可能形成一种特殊的导水通道，引发突水等安全事故，严重时可能造成人员伤亡1 2。精准探明掌子面前方陷落柱的位置、范围与形态等参数尤其重要。而实际煤层的物性条件复杂，导致实际勘探的波场具有复杂性与不确定性。因此通过数值模拟研究含陷落柱煤系地层中的反射槽波传播规律具有重要意义。2011 年杨思通提出反射 Rayleigh 型槽波的 SV 分量可以作为超前探测陷落柱的有效波3。2012 年王勃等通过二维正演模拟分析了掌子面前方存在陷落柱时的地震波场特征，总结出基于绕射偏移的陷落柱边界探测方法4。目前，诸多学者的研究都是基于弹性各向同性的理想模型条件进行数值模拟，而实际地层具有黏弹性质，导致数值模拟结果与实际煤矿探测数据存在误差5。苑春方等6、刘瑞珣7发现Kelvin-Voigt 黏弹性介质模型相比其他黏弹模型考虑了介质对地震波能量的损耗，使用 Kelvin-Voigt 模型更加符合地下介质是粘弹性质的实际情况。目前基于黏弹性质的槽波超前探测理论研究很少。因此，笔者在前人的研究基础上，基于黏弹三维一阶速度应力方程，采用交错网格高阶有限差分方法进行槽波超前探测陷落柱的三维数值模拟，分析掌子面前方存在陷落柱时的波场特征与反射槽波的传播规律。2.原理原理 2.1.槽波超前探测原理槽波超前探测原理 超前探测指的是探测掘进巷道前方一定距离内的构造情况。超前探测的测线布置与侧帮探测相似，是布置在巷道侧帮煤壁上，检波器安装在浅孔内或锚杆头上，或者侧帮打深孔，在深孔内安装多个检波器。不同之处在于，为提高效率、减少探测施工对掘进工作的影响，超前探测的炮数较少，有时仅在测线的两端布置炮点8。槽波是一种仅沿煤层传播的制导波，具有传播距离远、能量强、波形特征易于识别、携带地质信息Open AccessOpen Access何帝秀 等 DOI:10.12677/me.2023.113052 418 矿山工程 丰富的特点，且具有明显的频散特征。槽波大部分能量不向围岩辐射，携带大量地质构造信息，加之煤层的波导性，槽波被应用于超前探测。图 1 为槽波超前探测原理图，槽波在沿超前探测方向传播过程中遇到陷落柱(存在波阻抗差异)时会产生反射槽波，对其进行处理、分析与解释，进而实现对掌子面前方地质构造异常区域的槽波超前探测9 10。Figure 1.Advance detection principle 图图 1.超前探测原理 2.2.三维黏弹性波动方程三维黏弹性波动方程 使用 vx、vy、vz变量表示位移的一阶导数，减少位移二阶导数的使用，简化整个计算过程。在没有受到外力影响或外力消失之后，根据 Kelvin-Voigt 模型，可获得三维黏弹一阶速度应力方程11：22211131311131322213111313111322131311131311=+=+=+yyxxxxzzyyyyxxzzyyxxzzzvvvvvvtxyzx ty tz tvvvvvvtxyzx ty tz tvvvvvtxyzx ty tCCCCCCCCC2224444224444224444=+=+=+zyzyyzzxzxxzzxyyyxxvz tvvvvtzyz ty tvvvvtzxz tx tvvvvtyxy tx tCCC (1)式中，xx、yy、zz为正应力；xy、xz、yz为剪应力；C为弹性矩阵，i、j 分别为矩阵的行号和列号；为黏滞矩阵；t 为传播时间，单位 s；vx、vy、vz分别为质点沿 x、y、z 方向的振动速度，x、y 和 z分别代表三维空间中的平面水平方向、平面竖直方向和纵向。2.3.有限差分方法有限差分方法 有限差分就是将求解区域进行网格划分，使用有限的网格节点取代连续的求解域，利用微商与差商的近似关系，将描述地震波传播的微分方程转化为差分方程。差分离散主要分两种，其一是将单变量的二阶波动方程直接转化为时间空间的二阶中心差分进行离散求解；其二是将用位移表示的二阶波动方程转化为用应力及质点速度表示的一阶方程组，即一阶速度应力方程，然后使用应力和速度的交错网格何帝秀 等 DOI:10.12677/me.2023.113052 419 矿山工程 进行离散求解12。低阶近似、规则网格的有限差分方法具有强频散和低精度的缺点13，为保证计算精度，本文使用交错网格高阶有限差分方法对公式 1 进行离散化。2L 阶空间差分精度、二阶时间差分精度三维各向同性介质交错网格高阶有限差分离散方程11，如下：()()()()()()()()()()()()()()()()(),2,2,nnnnxxxxxxyynzznnnnyyyyxxyyi j ki j ki j ki j ktLVij ki j ktLVi jki j ktLVi j ki j ki j ki j ktLVij ki j ki j ktLVi jki j kt+=+=+(),nzzLVi j k+()()()()()()()()()()()()()(),2,(,),+=+=+=nnnnzzzzxxyynzznnnnxyxyyxxynnxzxzi j ki j ki j ktLVij ki j ktLVi jki j ki j ktLVi j kijkijki j kt LVij kLVi jkij kij()()()(),+nnzxxzki j kt LVij kLVi j k (2)()()()()()()()()()()()()()()(),+=+=+=+nnnnyzyzzyyznnnnnxxxxxyxyzxznnnyyxxyi jki jki j kt LVi jkLVi j ktVij kVij kLi j kLijkLij ki j ktVi jkVi jkLijki j k()()()()()()()(),+=+nnyyyzyznnnnnzzxxzyyzzzzLi j kLi jktVi j kVi j kLij kLi jkLi j ki j k 式中，t 为时间差分间隔，上标 n 表示 n t 时刻，+n表示 n+1/2，n表示 n 1/2；i、j、k 表示空间差分格点位置；+xL等为差分算子。2.4.吸收边界与稳定性条件吸收边界与稳定性条件 地震波正演数值模拟不仅要考虑数值算法的构造，还需要考虑边界条件与稳定性条件。完美匹配层(Perfect Matched Layer，简称 PML)原理为：将求解对象的变量分为两个部分，平行于边界的分量和垂直于边界的分量，人为在计算区域四周加上完全匹配层，可以充分的吸收边界反射，且其可以在短时间内吸收衰减沿人工界面法向的各类平面波，从而降低人工界面所产生的各类反射波，除此之外，和人工界面彼此平行的平面波不会发生衰减。在理论和模型算例上证明该方法可以很好吸收来自各个方向、各种频率的波14。为避免产生较大的数值频散，在数值模拟时，时间步长、空间步长以及物理参数的选取均需满足差分计算的稳定性要求15。根据傅里叶分析，得到交错网格有限差分法各类差分精度的稳定性条件，如下11：()22211111=+pLLiitvxyzC (3)式中：vp指介质纵波速度，x、y 和 z 分别为 x、y 和 z 方向的网格步长，()LiC指差分系数。何帝秀 等 DOI:10.12677/me.2023.113052 420 矿山工程 3.模型与参数模型与参数 模型在 X、Y、Z 方向的大小分别为 300 m 100 m 25 m，图 2、图 3 分别是三维含陷落柱煤系模型的 XOY 切面和 XOZ 切面，其中 X 方向为平行巷道方向，Y 方向为垂直巷道方向，Z 方向为垂直方向。模型网格 dx、dy、dz 分别为 1 m 1 m 0.25 m，煤厚 5 m。煤巷长度 150 m，宽度 4 m，高度 4 m，起点坐标为(10,48,10.5)。掌子面位于 X=160 m 处，陷落柱中心点坐标为(240,50,12.5)，半径 D=20 m。使用主频为 150 Hz 的雷克子波作为激发震源，延迟时间为 20 ms，时间的采样间隔为 0.05 ms，震源坐标为(110,54,12.5)。检波器与震源同一水平线布置，布测于 X 方向 80160 m，道间距 1 m，共 81 个检波点。根据相关文献中煤层和围岩物性参数常见值的范围划分16，设计模型各物性参数如表 1 所示，煤层的顶底板岩性参数一致，巷道赋予空值。Figure 2.3D coal series model XOY section 图图 2.三维煤系模型 XOY 切面 Figure 3.3D coal series model XOZ section 图图 3.三维煤系模型 XOZ 切面 Table 1.Physical parameter table 表表 1.物性参数表 层位 vp/(ms1)vs/(ms1)/(kgm3)Qp Qs 顶板 3500 2000 2400 200 100 煤层 1900 1100 1300 100 50 陷落柱 2400 1400 1800 150 100 底板 3500 2000 2400 200 100 何帝秀 等 DOI:10.12677/me.2023.113052 421 矿山工程 4.正演模拟与分析正演模拟与分析 4.1.波场快照波场快照 图 4 为正演模拟的 X 分量 40 ms、80 ms、120 ms 与 160 ms 在 Z=12.5 m 的波场快照。图 4(a)为传播时间 40 ms，震源激发后，产生纵波、横波，体波传播速度快，向超前方向传播率先经过掌子面，一部分体波在掌子面发生反射。图 4(b)为传播时间 80 ms，直达纵波已经过陷落柱，且在陷落柱左界面(靠近掌子面一侧)产生反射波和透射波，反射波被直达横波干扰，波形不明显，陷落柱体内的透射波能量有衰减，直达纵波的波形变化较小，部分速度较快的波已传播至模型右侧壁(X=300 m)；直达横波紧跟直达纵波；在体波传播过程中，纵波与横波的 SV 波相互叠加、干涉形成 Rayleigh 槽波，与体波相比，直达Rayleigh 槽波能量较大、速度较慢；各直达波在掌子面处形成的系列掌子面反射波已沿巷道向后传播扩散。图 4(c)为传播时间 120 ms，直达横波经过陷落柱，透射进陷落柱的波能量有衰减，波形有变化。图4(d)为传播时间 160 ms，反射横波到达掌子面；直达 Rayleigh 槽波传播至陷落柱，在陷落柱界面产生反射波和透射波，透射入陷落柱的波能量衰减大，波形有明显变化，波形不连续；掌子面反射波沿巷道向后已传播至巷道中部；地震波传播至巷道头时，亦形成了巷道头反射波。(a)40 ms (b)80 ms (c)120 ms (d)160 ms Figure 4.X-component wave field snapshot(In the figure:wave group 1direct longitudinal wave,wave group 2direct transverse wave,wave group 3direct in-seam wave,wave group 4reflected wave on the palm surface,wave group 5reflected longitudinal wave,wave group 6reflected transverse wave,wave group 7reflected in-seam wave,and wave group 8reflected wave at the tunnel head)图图 4.X 分量波场快照(图中：波组 1直达纵波，波组 2直达横波，波组 3直达槽波，波组4掌子面反射波，波组 5反射纵波，波组 6反射横波，波组 7反射槽波，波组 8巷道头反射波)图 5 为正演模拟的 Y 分量 40 ms、80 ms、120 ms 与 160 ms 在 Z=12.5 m 的波场快照。观察发现，Y分量地震波传播的波场类似于 X 分量。图 5(a)为传播时间 40 ms，震源激发产生纵波、横波。图 5(b)为传播时间 80 ms，直达纵波已经过陷落柱，透射入陷落柱的波能量有衰减，遇陷落柱波形有变化；在传播过程中，SH 波相互干涉形成 Love 槽波。图 5(c)为传播时间 120 ms，直达横波已传播至陷落柱，透射进陷落柱的波能量衰减小，波形有变化。图 5(d)为传播时间 160 ms，直达 Love 槽波传播至陷落柱，产生反射Love 槽波，透射入陷落柱的波能量衰减小，波形变化不明显；反射横波到达掌子面；掌子面反射波沿巷道向后传播至巷道中部。214213421364378何帝秀 等 DOI:10.12677/me.2023.113052 422 矿山工程 (a)40 ms (b)80 ms (c)120 ms (d)160 ms Figure 5.Y-component wave field snapshot 图图 5.Y 分量波场快照 图 6 分别为正演模拟的 Z 分量 40 ms、80 ms、120 ms 与 160 ms 在 Z=12.5 m 的波场快照。观察发现，Z 分量中地震波的传播与 X 分量相似。图 6(a)为传播时间 40 ms，震源激发产生纵波、横波。图 6(b)为传播时间 80 ms，直达纵波已经过陷落柱，透射入陷落柱的波能量增大，波形变化明显；纵波与横波的SV 波在传播过程中相互干涉形成 Rayleigh 槽波。图 6(c)为传播时间 120 ms，直达横波传播至陷落柱，透射入陷落柱的波能量衰减较小，波形有变化。图 6(d)为传播时间 160 ms，反射纵波已快传播至巷道中部；反射横波到达掌子面；直达 Rayleigh 槽波在陷落柱界面波形不连续，产生的反射 Rayleigh 槽波紧跟反射横波。(a)40 ms (b)80 ms (c)120 ms (d)160 ms Figure 6.Z-component wave field snapshot 图图 6.Z 分量波场快照 不同分量上地震波受陷落柱的影响不同。X 分量的直达 Rayleigh 槽波受陷落柱影响较大，波形不连续，陷落柱的边界范围明显可见；遇陷落柱产生的反射 Rayleigh 槽波波形模糊。Y 分量的直达 Love 槽波受陷落柱影响较小；遇陷落柱产生的反射 Love 槽波波形较清晰。Z 分量的直达 Rayleigh 槽波受陷落柱影2134214312743268213214432647358何帝秀 等 DOI:10.12677/me.2023.113052 423 矿山工程 响较大；反射 Rayleigh 槽波波形较清晰。4.2.地震记录地震记录 图 7 为合成地震记录的三分量图，图中横轴为检波器编号，纵轴为传播时间。结合波场快照和波的时距曲线特征，图 7 各分量上按被接收的时间顺序可识别出：直达纵波、直达横波、直达槽波、掌子面反射波、反射纵波、反射横波、反射槽波和巷道头反射波，与波场快照中所识别的波组相对应。图 7(a)为 X 分量地震记录，反射 Rayleigh 槽波在能量较大掌子面反射波的干扰下难以识别，同相轴较模糊。图7(b)为 Y 分量地震记录，反射 Love 槽波能量较大，与其他波列在时间轴上间隔大，波列易识别。图 7(c)为 Z 分量地震记录，陷落柱的各反射波波组连续性好，反射 Rayleigh 槽波能量较小，与掌子面反射波振幅大小相似，波列较易识别。震源激发后，沿超前方向传播的体波和槽波遇到陷落柱表面产生了振幅较弱的反射纵波和反射横波，在 X 与 Z 分量产生了能量较小的反射 Rayleigh 槽波，在 Y 分量产生了振幅较强、同相轴明显的反射 Love槽波。因此可以利用反射 Love 槽波超前探测陷落柱。(a)X 分量 (b)Y 分量 (c)Z 分量 Figure 7.Three component map of seismic records 图图 7.地震记录三分量图 5.结论结论 通过数值模拟含陷落柱的黏弹性煤系地质模型，分析波场快照与地震记录，对槽波超前探测陷落柱得出如下结论：1)通过三个分量的波场快照图发现，地震波在沿超前方向传播时，遇陷落柱皆产生了反射槽波。X和 Z 分量的 Rayleigh 槽波受陷落柱影响大，Y 分量上 Love 槽波受陷落柱影响较小。2)三分量地震记录上，X 与 Z 分量接收的反射 Rayleigh 槽波能量较小，Y 分量接收到的反射 Love槽波能量较大，同相轴明显。建议利用 Y 分量的反射 Love 槽波超前探测陷落柱。限于篇幅，本文数值模拟模型仍不够多，后续还需建立更多模型更全面的研究陷落柱的反射槽波响应特征。6821347672134875213468何帝秀 等 DOI:10.12677/me.2023.113052 424 矿山工程 参考文献参考文献 1 李刚.透射槽波在煤矿陷落柱探测中的应用J.煤矿开采,2016,21(3):39-42.2 马欣,杨思通,李新凤,等.基于透射槽波的采煤工作面陷落柱探测模拟研究J.煤矿安全,2019,50(4):32-36.3 杨思通.矿井巷道地震超前探测三维全波场数值模拟与探测方法研究D:博士学位论文.青岛:山东科技大学,2011.4 王勃,刘盛东,胡泽安.陷落柱地震波超前探测数值模拟与应用J.中国煤炭地质,2012,24(3):56-59+65.5 李晓波,董良国.粘弹介质中可变网格地震波传播数值模拟J.石油物探,2012,51(1):1-10+111.6 苑春方,彭苏萍,张中杰,刘振宽.Kelvin-Voigt 均匀黏弹性介质中传播的地震波J.中国科学:D 辑,2005,35(10):957-962.7 刘瑞珣,张秉良,张臣.描述岩石粘弹性固体性质的开尔文模型J.地学前缘,2008,15(3):221-225.8 王季,叶红星,张广忠,等.煤矿反射槽波探测技术研究评述J.煤田地质与勘探,2023,51(2):292-300.9 程建远,江浩,姬广忠,吴海.基于节点式地震仪的煤矿井下槽波地震勘探技术J.煤炭科学技术,2015,43(2):25-28.10 高玉超,邓重青,李继路.槽波超前勘探技术在巷道前方探测断层的应用J.山东煤炭科技,2020(4):166-168.11 姬广忠.煤巷侧帮反射槽波成像方法及应用研究D:博士学位论文.北京:煤炭科学研究总院,2017.12 杨莹.二维地震波场有限差分法数值模拟研究D:硕士学位论文.北京:中国地质大学(北京),2009.13 陈可洋.高阶弹性波波动方程正演模拟及逆时偏移成像研究D:硕士学位论文.大庆:大庆石油学院,2009.14 姬广忠,程建远,朱培民.煤层 Love 型槽波数值模拟及其频散特征分析J.煤炭科学技术,2011,39(6):106-109.15 Liu,Q.-H.,Schoen,E.,Daube,F.,et al.(1996)A Three-Dimensional Finite Difference Simulation of Sonic Logging.The Journal of the Acoustical Society of America,100,72-79.https:/doi.org/10.1121/1.415869 16 李东会.煤储层各向异性波场模拟与特征分析D:博士学位论文.徐州:中国矿业大学,2012.