·1314 北京科技大学学报 第33卷 向刚度，k为平行链接的切向刚度，N.为平行链接 行数值模拟试验.岩石样品加工成高80mm、直径 数，A为链接的截面积，I为链接的转动惯量，k、为 40mm的圆柱形试样，其弹性模量为42.6GPa,泊松 接触法向刚度 比为0.21.PFC2D数值模型建立与常用的连续介 当岩石模型受到外荷载的作用后，颗粒处会产 质程序不同，连续介质程序建模使用材料参数是宏 生局部应力，当其值超过链接强度时，链接断裂，链 观参数，可以由实验得到，但目前还没有理论能够反 接的应变能释放.在岩石模型中，颗粒是以非常紧 映模型中的细观参数与模型最终所反映的宏观力学 密的方式来排列，颗粒能量释放的表现方式是从链 性能的联系（这是今后研究的重点），因而PFC2D无 接断裂处以地震波的形式向外传播，这和真实岩石 法直接给模型设定细观力学参数和初始应力条件， 中声传播的过程非常相似 这需要采用逆向模拟的方法反复进行数值模拟试 地震波在岩石中传播时，由于岩石对波的吸收 验，先做一系列物理试验，得到材料的各种宏观参 特性，会有部分机械能转变成热能，从而使波在岩石 数，然后不断调整构成模型介质的基本粒子级配组 中的传播不断衰减，通常由地震波品质因子Q来表 成、接触方式和相应的细观力学参数，直至使数值模 示，是在波传播一个波长入的距离后，原来储存的 拟结果与物理试验结果相吻合).岩样内部的不 能量E与消耗能量△E之比的2π倍☒，即 均匀性通过设置颗粒链接强度空间分布的不均匀性 E 来实现，设定颗粒平行链接的平均法向强度G。及平 Q=2AE (2) 均切向强度T。后，再分别设置两者的标准差来体现 地震波品质因子和岩石的微观性质，如岩石内 出岩体在细观层面的不均匀性.当强度的标准差设 部裂纹的密度、分布和构造等相关.Q值越大，波通 定的比较小，材料宏观上表现比较均匀，以线性性质 过岩石的衰减越小，在PFC2D中通过定义阻尼系数 为主，而当标准差较大时，材料宏观的非线性特征显 α来模拟岩石的这一特征.对于单自由度系统，α在 著，达到峰值应力后期曲线出现明显的弱化行为. 一个循环中的能量损失可由下式来表示： 根据Hazzard和Young的经验，花岗岩的原位地 asAb 震波品质因子值大约为200，但对于室内试验的岩 (3) 4E 样而言，其地震波品质因子会有很大的衰减，因此本 因此，由式(2)和式(3)，可得 文中取Q为100，相应的阻尼系数a取值为0.0157. 岩样由6984个颗粒组成，粒径为0.25~0.45mm, T a=20 (4) 服从均匀分布，计算模型的参数以Lac du Bonnet花 通过式(4)可以在P℉C2D中定义阻尼系数来确 岗岩的P℉C2D力学参数为基础，进行单轴压缩 定所需要的弹性波衰减水平 试验，以其所得宏观力学参数弹模、泊松比及单轴抗 压强度等与室内物理试验结果相比较，然后不断对 2计算模型和结果分析 细观力学参数进行调整，再进行数值单轴压缩，直至 2.1计算模型及参数 结果和室内物理试验基本相符为止，具体参数见 根据上述原理，对取自某煤矿巷道的花岗岩进 表1. 表1PFC2D模型细观力学参数表 Table 1 Meso-structure parameters of PFC2D model 颗粒接触 平行链接半平行链接刚颗粒摩擦 平均法向强度， 法向强度标 平均切向强度，切向强度标 阻尼 模量，E。/CPa 径系数，入度比，瓜 因数μ G/MPa 准差，△d。/MPa 7/MPa 准差，△F。MPa 系数，a 5.7×100 1.0 2.5 0.5 157 36 157 36 0.017 数值模型建立后，通过岩样上端的墙来施加荷 2.2试验结果及分析 载，荷载施加速度为6.7×10-6mm·step-1,以保证 从试验结果可以看出，岩样的破坏形式与围压 在整个试验过程中试样都能保持静力平衡状态.根 的大小有关，围压为5MPa时，岩样的峰值抗压强度 据工程中实测地应力数据和应力分布特征，确定了 较低，破坏形式属于突发失稳.随围压的增加，岩样 三种不同的应力水平（σ3=5、20和50MPa)进行双 的峰值抗压强度有明显提升，破坏前的应变增加，岩 轴压缩试验 石的塑性增大，破坏形式逐渐转变为渐进破坏，见北 京 科 技 大 学 学 报 第 33 卷 向刚度，ks 为平行链接的切向刚度，Npb为平行链接 数，A 为链接的截面积，I 为链接的转动惯量，kN 为 接触法向刚度． 当岩石模型受到外荷载的作用后，颗粒处会产 生局部应力，当其值超过链接强度时，链接断裂，链 接的应变能释放． 在岩石模型中，颗粒是以非常紧 密的方式来排列，颗粒能量释放的表现方式是从链 接断裂处以地震波的形式向外传播，这和真实岩石 中声传播的过程非常相似． 地震波在岩石中传播时，由于岩石对波的吸收 特性，会有部分机械能转变成热能，从而使波在岩石 中的传播不断衰减，通常由地震波品质因子 Q 来表 示，是在波传播一个波长 λ 的距离后，原来储存的 能量 E 与消耗能量 ΔE 之比的 2π 倍［12］，即 Q = 2π E ΔE ( 2) 地震波品质因子和岩石的微观性质，如岩石内 部裂纹的密度、分布和构造等相关． Q 值越大，波通 过岩石的衰减越小，在 PFC2D 中通过定义阻尼系数 α 来模拟岩石的这一特征． 对于单自由度系统，α 在 一个循环中的能量损失可由下式来表示: α = ΔE 4E ( 3) 因此，由式( 2) 和式( 3) ，可得 α = π 2Q ( 4) 通过式( 4) 可以在 PFC2D 中定义阻尼系数来确 定所需要的弹性波衰减水平． 2 计算模型和结果分析 2. 1 计算模型及参数 根据上述原理，对取自某煤矿巷道的花岗岩进 行数值模拟试验． 岩石样品加工成高 80 mm、直径 40 mm 的圆柱形试样，其弹性模量为 42. 6 GPa，泊松 比为 0. 21． PFC2D 数值模型建立与常用的连续介 质程序不同，连续介质程序建模使用材料参数是宏 观参数，可以由实验得到，但目前还没有理论能够反 映模型中的细观参数与模型最终所反映的宏观力学 性能的联系( 这是今后研究的重点) ，因而 PFC2D 无 法直接给模型设定细观力学参数和初始应力条件， 这需要采用逆向模拟的方法反复进行数值模拟试 验，先做一系列物理试验，得到材料的各种宏观参 数，然后不断调整构成模型介质的基本粒子级配组 成、接触方式和相应的细观力学参数，直至使数值模 拟结果与物理试验结果相吻合［13］． 岩样内部的不 均匀性通过设置颗粒链接强度空间分布的不均匀性 来实现，设定颗粒平行链接的平均法向强度 σc 及平 均切向强度 τc 后，再分别设置两者的标准差来体现 出岩体在细观层面的不均匀性． 当强度的标准差设 定的比较小，材料宏观上表现比较均匀，以线性性质 为主，而当标准差较大时，材料宏观的非线性特征显 著，达到峰值应力后期曲线出现明显的弱化行为． 根据 Hazzard 和 Young ［14］的经验，花岗岩的原位地 震波品质因子值大约为 200，但对于室内试验的岩 样而言，其地震波品质因子会有很大的衰减，因此本 文中取 Q 为 100，相应的阻尼系数 α 取值为0. 015 7． 岩样由 6 984 个颗粒组成，粒径为 0. 25 ～ 0. 45 mm， 服从均匀分布，计算模型的参数以 Lac du Bonnet 花 岗岩的 PFC2D 力学参数为基础［11］，进行单轴压缩 试验，以其所得宏观力学参数弹模、泊松比及单轴抗 压强度等与室内物理试验结果相比较，然后不断对 细观力学参数进行调整，再进行数值单轴压缩，直至 结果和室内物理试验基本相符为止，具体参数见 表 1． 表 1 PFC2D 模型细观力学参数表 Table 1 Meso-structure parameters of PFC2D model 颗粒接触 模量，Ec /GPa 平行链接半 径系数，λ 平行链接刚 度比，kn /ks 颗粒摩擦 因数，μ 平均法向强度， σc /MPa 法向强度标 准差，Δσc /MPa 平均切向强度， τc /MPa 切向强度标 准差，Δτc /MPa 阻尼 系数，α 5. 7 × 1010 1. 0 2. 5 0. 5 157 36 157 36 0. 017 数值模型建立后，通过岩样上端的墙来施加荷 载，荷载施加速度为 6. 7 × 10 － 6 mm·step － 1 ，以保证 在整个试验过程中试样都能保持静力平衡状态． 根 据工程中实测地应力数据和应力分布特征，确定了 三种不同的应力水平( σ3 = 5、20 和 50 MPa) 进行双 轴压缩试验． 2. 2 试验结果及分析 从试验结果可以看出，岩样的破坏形式与围压 的大小有关，围压为 5 MPa 时，岩样的峰值抗压强度 较低，破坏形式属于突发失稳． 随围压的增加，岩样 的峰值抗压强度有明显提升，破坏前的应变增加，岩 石的塑性增大，破坏形式逐渐转变为渐进破坏，见 ·1314·