第12期 祁建东等：Bandis经验公式修正及其应用 ·1577· K=1.45a88. 表2岩石物理力学参数 Table 2 Physico-mechanical characteristics of rock [1-0.3050.+0.029j(0.360.+0.03》 岩石密度， 体积模量， 剪切模量， 摩擦角， 黏聚力， (7 p/(kg'm-3) K/GPa G/GPa /() c/MPa 表1剪切刚度参数拟合 2220 8.45 7.71 10.19 Table 1 Shear stiffness parameters fitting 利用3DEC内置函数GEN划分网格，共生成 Go/ Tp/ K./ Tulk/ 40930个有限差分域. MPa MPa (GPa-m-1)MPa 边界条件：在z=60及x=-30处施加不同应 0.2 0.10 0.55 0.104 力边界，其余边界面在各自法向量方向施加速度边 0.4 0.14 0.73 0.169 1.45 0.63 界v=0. 0.60.21 1.08 0.215 3.1基于Bandis经验公式数值实验 0.8 0.28 1.30 0.315 合理选取结构面刚度是数值分析结构面变形特 2.2.2法向刚度 性的保障.首先采用Bandis经验公式选取剪切刚度 利用实测数据通过共轭梯度法拟合得到lgσ。- 及法向刚度参数如表3所示，并据此进行了三组不 △V关系式： 同应力组合的数值实验.分析中利用History函数 lg0m=7.34△V-1.77. (8) 监测记录结构面中心点处的剪切位移. 结合式(6)及式(8)得到结构面法向刚度K。: 表3基于Bandis经验公式法选取的刚度参数 7.340m Table 3 Stiffness parameters for numerical experiments according to the K.=0.3434 (9) Bandis formula 3 Bandis经验公式修正 K/ K/ o。/MPa T/MPa (MPa“mm-l) (MPa'mm-) 通过基于室内试验结果的数值实验，利用反演 0.05 0.317 4.275 法对Bandis经验公式进行修正.本文中数值实验采 0.20 0.10 0.005 4.275 用三维离散元程序3DEC.该程序是基于粗糙结构 0.07 0.489 8.336 面的假设而开发的，虽然不能探究结构面不同性状 0.40 0.14 0.008 8.336 特征对其变形特性的不同影响，但对结构面的变形 0.11 0.595 12.611 特性进行定性半定量的分析是可行的. 0.60 0.21 0.009 12.611 综合考量边界效应影响及简化计算后建立模型 0.07 1.123 16.886 为60mm×60mm×120mm立方体，如图1所示. 0.14 0.755 16.886 0.80 0.21 0.304 16.886 0.28 0.011 16.886 当法向应力σ。=0.8MPa保持不变时，等差取 结构面 剪应力r直至r=T。=0.28MPa共4种组合，数值 实验结果如图2所示. 当法向应力σ。依次取为0.2、0.4、0.6和 0.8MPa时，选取各自对应的剪应力r=0.5r,共四种 组合，数值实验结果如图3所示 当法向应力o。依次取为0.2、0.4、0.6和 0.8MPa时，选取各自对应的剪应力T=T,共四种组 图1数值实验模型 Fig.1 Numerical experiment model 合，数值实验结果如图4所示. 从图2~图4可以看出：基于Bandis经验公式 选用刚性块体模型模拟砂岩试件，结构面选用 的刚度参数在剪应力相对较小（即，≤0.5r,)时，能 Coulomb滑移模型.砂岩物理力学参数见表2. 够得到与室内试验结果相近的剪切变形数值模拟第 12 期 祁建东等: Bandis 经验公式修正及其应用 Kst = 1. 45σ0. 63 n [ · 1 － τ 2 ( 0. 305σn + 0. 029) ( 0. 336σn + 0. 033 ] ) 2 ． ( 7) 表 1 剪切刚度参数拟合 Table 1 Shear stiffness parameters fitting σn / MPa τp / MPa Ksi / ( GPa·m － 1 ) τult / MPa Kj nj 0. 2 0. 10 0. 55 0. 104 0. 4 0. 14 0. 73 0. 169 1. 45 0. 63 0. 6 0. 21 1. 08 0. 215 0. 8 0. 28 1. 30 0. 315 2. 2. 2 法向刚度 利用实测数据通过共轭梯度法拟合得到 lgσn － ΔVj关系式: lgσn = 7. 34ΔVj － 1. 77． ( 8) 结合式( 6) 及式( 8) 得到结构面法向刚度 Kn : Kn = 7. 34σn 0. 3434． ( 9) 3 Bandis 经验公式修正 通过基于室内试验结果的数值实验，利用反演 法对 Bandis 经验公式进行修正． 本文中数值实验采 用三维离散元程序 3DEC． 该程序是基于粗糙结构 面的假设而开发的，虽然不能探究结构面不同性状 特征对其变形特性的不同影响，但对结构面的变形 特性进行定性半定量的分析是可行的． 综合考量边界效应影响及简化计算后建立模型 为 60 mm × 60 mm × 120 mm 立方体，如图 1 所示． 图 1 数值实验模型 Fig． 1 Numerical experiment model 选用刚性块体模型模拟砂岩试件，结构面选用 Coulomb 滑移模型． 砂岩物理力学参数见表 2． 表 2 岩石物理力学参数 Table 2 Physico-mechanical characteristics of rock 岩石密度， ρ /( kg·m － 3 ) 体积模量， K /GPa 剪切模量， G /GPa 摩擦角， φ/( °) 黏聚力， c/MPa 2220 8. 45 7. 71 43 10. 19 利用 3DEC 内置函数 GEN 划分网格，共生成 40930 个有限差分域． 边界条件: 在 z = 60 及 x = － 30 处施加不同应 力边界，其余边界面在各自法向量方向施加速度边 界 v = 0． 3. 1 基于 Bandis 经验公式数值实验 合理选取结构面刚度是数值分析结构面变形特 性的保障． 首先采用 Bandis 经验公式选取剪切刚度 及法向刚度参数如表 3 所示，并据此进行了三组不 同应力组合的数值实验． 分析中利用 History 函数 监测记录结构面中心点处的剪切位移． 表 3 基于 Bandis 经验公式法选取的刚度参数 Table 3 Stiffness parameters for numerical experiments according to the Bandis formula σn /MPa τ /MPa Kst / ( MPa·mm － 1 ) Kn / ( MPa·mm － 1 ) 0. 20 0. 05 0. 317 4. 275 0. 10 0. 005 4. 275 0. 40 0. 07 0. 489 8. 336 0. 14 0. 008 8. 336 0. 60 0. 11 0. 595 12. 611 0. 21 0. 009 12. 611 0. 07 1. 123 16. 886 0. 80 0. 14 0. 755 16. 886 0. 21 0. 304 16. 886 0. 28 0. 011 16. 886 当法向应力 σn = 0. 8 MPa 保持不变时，等差取 剪应力 τ 直至 τ = τp = 0. 28 MPa 共 4 种组合，数值 实验结果如图 2 所示． 当法 向 应 力 σn 依 次 取 为 0. 2、0. 4、0. 6 和 0. 8 MPa时，选取各自对应的剪应力 τ = 0. 5τp共四种 组合，数值实验结果如图 3 所示． 当法 向 应 力 σn 依 次 取 为 0. 2、0. 4、0. 6 和 0. 8 MPa时，选取各自对应的剪应力 τ = τp共四种组 合，数值实验结果如图 4 所示． 从图 2 ～ 图 4 可以看出: 基于 Bandis 经验公式 的刚度参数在剪应力相对较小( 即 τ≤0. 5τp ) 时，能 够得到与室内试验结果相近的剪切变形数值模拟 · 7751 ·