三维椭球形缺陷启裂载荷随倾角和压缩度变化的规律1)

高经纬 肖 艳 范立峰

(北京工业大学城市建设学部, 北京 100124)

工程岩体内赋存着许多缺陷。岩体内的缺陷在外部载荷作用下会发生启裂，随后进一步扩展贯通形成宏观裂纹，这会使岩体的强度降低，进而可能危害工程安全[1]。缺陷的尺寸和倾角都会影响其启裂行为，研究缺陷尺寸和倾角对其启裂行为的影响对认识岩体内缺陷的启裂破坏机理以及防治工程灾害具有重要意义。

早期学者通常采用平面二维缺陷模型来简化研究三维缺陷启裂行为。李冰峰等[2]研究了二维缺陷的启裂及扩展行为，表明缺陷的存在会削弱矿物颗粒形态对岩石破坏行为的影响。Barquins等[3]将岩体内缺陷简化成二维椭圆形缺陷，研究了缺陷倾角和扁平程度（扁椭圆形缺陷短轴与长轴之比）对其启裂位置的影响。Liu等[4]研究了预制缺陷的倾角和长度对其启裂位置以及压缩强度的影响，研究表明缺陷倾角会影响其启裂位置，并且倾角对压缩强度影响程度与缺陷长度密切相关，随着缺陷长度增加，倾角对强度影响越大。郭少华等[5]研究了压应力作用下缺陷倾角对其启裂载荷的影响，研究表明椭圆形缺陷启裂载荷随着其倾角的增加先减小后增加。椭圆形缺陷存在一个临界倾角，当缺陷处于该倾角时，其启裂载荷最小。缺陷倾角对其启裂载荷的影响规律受缺陷形态影响。赵永等[6]研究发现条形缺陷启裂载荷随其倾角增加单调减小。Ma等[7-8]将岩体内缺陷简化成直裂纹、Z型和S型非直裂纹，研究了直裂纹和非直裂纹在轴向压缩载荷下的启裂位置，揭示了缺陷形态对其启裂行为影响的规律。辛宏举等[9]研究了圆形缺陷尺寸对压缩开裂及破坏模式的影响。Zhang等[10]将岩体内缺陷简化成圆形，矩形和马蹄形，研究了缺陷形态对试样裂纹演化和破坏模式的影响。平面二维缺陷的启裂行为受缺陷形态、尺寸以及倾角等因素的影响。前人对二维缺陷启裂行为的研究对于认识缺陷启裂破坏机理提供了很大帮助，开展三维缺陷启裂行为的研究将有助于进一步揭示岩体内复杂形态缺陷的启裂扩展破坏机理。

目前，前人对复杂的三维缺陷启裂扩展行为开展了相关研究。Germanovich等[11]将岩体内缺陷简化成三维币形缺陷，研究了三维币形缺陷在轴向载荷作用下的扩展行为。Fu等[12]将岩体内缺陷简化成三维椭圆形缺陷，研究了三维椭圆形缺陷在轴向载荷作用下的启裂扩展规律。他们侧重于研究三维缺陷在轴向载荷作用下的扩展模式。Yang等[13]将岩体内缺陷简化成三维椭圆形缺陷，研究了缺陷倾角对试样压缩破坏强度的影响，结果表明试样压缩破坏强度随着缺陷倾角增加先减小后增加。岩体内缺陷的启裂、扩展和破坏行为与缺陷形态、角度以及尺寸密切相关。前人关于三维缺陷启裂行为的研究通常将三维缺陷简化成三维币形和三维椭圆形缺陷，侧重讨论缺陷的倾角和尺寸对其扩展过程和试样力学强度的影响。实际工程岩体在宏观尺度和细观尺度上都含有类似椭球形的缺陷，如图1所示。椭球形缺陷三个方向上的尺寸对其启裂行为均会产生影响。另外，不仅需要关注岩体中缺陷的扩展过程以及岩体的强度，同样也需要关注缺陷的启裂行为，比如启裂载荷，这对于工程灾害预报具有重要意义。

图1 岩体内的椭球形缺陷Fig. 1 Ellipsoidal flaws in rocks

本文将岩体内缺陷简化成三维椭球形缺陷，基于扩展有限元法模拟三维椭球形缺陷的启裂行为，讨论三维椭球形缺陷倾角和尺寸对其临界启裂载荷的影响。

1 三维椭球形缺陷启裂行为数值模拟

1.1 模型的几何参数

图2 (a)呈现了三维椭球形缺陷模型的几何尺寸。模型为长、宽、高分别是20 mm，20 mm，40 mm的长方体，三维椭球形缺陷位于模型的正中心。o′x′y′z′坐标系是模型的整体坐标系，并且模型是关于x′o′z′平面对称的。oxyz坐标系是模型的局部坐标系，以三维椭球形缺陷的三向轴为坐标轴。缺陷倾角β为加载方向与x轴正方向的夹角。

通过设置缺陷体积不变，改变椭球形缺陷三个方向上的尺寸，得到不同尺寸的三维椭球形缺陷。缺陷最大扩展长度大约是初始缺陷半径的1.50倍[14]。为了保证模型边界不影响缺陷扩展，设置缺陷最大半径不超过模型边界长度的五分之一。本研究设置缺陷z方向上的最小半径为0.2 mm，x或y方向上的最大半径为3.5 mm。表1呈现了本研究设置的三维椭球形缺陷尺寸。每一种尺寸的三维椭球形缺陷包括6种倾角，分别是0°，15°，30°，45°，60°和75°。

表1 三维椭球形缺陷尺寸Table 1 Geometric dimensions of the threedimensional ellipsoid flaw model

1.2 模型的力学参数

参考文献[15-17]中关于采用Vero Clear材料开展脆性岩石材料内缺陷启裂和扩展行为的数值模拟和试验研究，本研究在数值模拟材料属性赋值时，同样采用了Vero Clear材料的基本力学参数。Vero Clear材料的基本力学参数，包括单轴抗压强度、弹性模量和泊松比[18-19]。冻结温度下，Vero Clear材料表现出明显的弹脆性特征。本研究是基于线弹性断裂力学理论进行的，在缺陷启裂阶段材料仍然处于弹性阶段。冻结温度下Vero Clear材料的断裂韧度KIC= 60.4 kPa·m1/2，弹性模量E= 59.60 MPa，泊松比µ= 0.30，临界抗拉强度σt= 1.00 MPa。

1.3 有限元模型

模型关于x′o′z′面对称。为了提高计算效率，可以简化建立对称模型。在模型的对称面x′o′z′面上设置对称约束，约束y′方向的位移。模型的下端面设置固定约束，在模型上端面施加均布面载荷，如图2 (b)所示。采用单轴压缩方式进行加载。网格采用十结点二次四面体单元(C3D10)，模型网格划分结果如图2 (c)所示。

图2 三维椭球形缺陷模型Fig. 2 3D ellipsoid flaw model

本研究采用了最大主应力准则作为判别启裂依据。前人研究表明，采用最大主应力准则作为判断缺陷启裂的依据时，模拟的缺陷启裂情况与实际缺陷启裂情况相符，模型计算过程的稳定性较高且容易收敛[20]。最大主应力准则表达形式为

式 中，σm0ax为 最 大 容 许 主 应 力，σmax为 最 大 主应力。当断裂准则f在容差允许范围内达到1.0时，缺陷启裂。如果f> 1.0 +ftol，则满足缺陷启裂标准。为了使模型计算效率高且容易收敛，使计算更加精确，将容差修改为f= 0.01。

1.4 模型的启裂载荷

采用三维椭球形缺陷有限元模型，基于扩展有限元法进行计算，根据计算结果的最大主应力云图，观察缺陷启裂情况，记录模型在启裂时的外部载荷，即模型的启裂载荷。以尺寸x= 3.00 mm，y= 2.00 mm，z= 0.20 mm，倾角β= 30°的椭球形缺陷模型为例介绍启裂载荷的确定方法。

图3呈现了尺寸x= 3.00 mm，y= 2.00 mm，z= 0.20 mm，倾角β= 30°的椭球形缺陷模型在不同载荷下的最大主应力云图。从图中可以看到缺陷顶点附近在载荷作用下出现了应力集中现象。随着载荷增加至584.8 kPa时，缺陷顶点附近发生启裂[4]。随着载荷进一步增加，缺陷发生扩展。该模型的启裂载荷为584.8 kPa。

图3 x = 3.00 mm， y = 2.00 mm， z = 0.20 mm， β = 30°模型在不同载荷下的最大主应力云图Fig. 3 Maximum principal stress nephogram of model with dimensions x = 3.00 mm, y = 2.00 mm, z = 0.20 mm, and flaw inclination angle β = 30° under different loads

为了更直观地比较不同尺寸及倾角的三维椭球形缺陷模型的启裂载荷，将所有模型的启裂载荷归一化。尺寸为x= 3.50 mm，y= 1.75 mm，z= 0.20 mm，倾角β= 0°模型的启裂载荷P0cr最大，为 3.44 MPa。各模型归一化后的启裂载荷计算公式为

式中，Pk为归一化启裂载荷，Pcr为模型的启裂载荷。

2 结果与讨论

2.1 倾角对不同尺寸椭球形缺陷启裂载荷的影响

图4 (a)呈现了z= 0.20 mm椭球形缺陷启裂载荷随倾角β的变化规律。尺寸x= 3.50 mm，y= 1.75 mm的椭球形缺陷启裂载荷随缺陷倾角β增加先快速减小后缓慢增加。当缺陷倾角β=30°时，模型启裂载荷最小。其他尺寸的椭球形缺陷模型启裂载荷随倾角的变化趋势与之类似。

图4 (b)呈现了z= 0.30 mm椭球形缺陷启裂载荷随倾角β的变化规律。尺寸x= 2.00 mm，y= 2.00 mm的椭球形缺陷启裂载荷随缺陷倾角β增加先减小后缓慢增加。同样，缺陷倾角β=30°时，模型启裂载荷最小。其他尺寸的椭球形缺陷模型启裂载荷随倾角的变化趋势与之类似。

图4 不同尺寸椭球形缺陷模型启裂载荷随倾角 β 变化图Fig. 4 Variation of the crack initiation load with inclination angle β for different models

图4(c)和图4(d)分别呈现了z= 0.40 mm和z= 0.50 mm椭球形缺陷启裂载荷随倾角β的变化规律。同样地，椭球形缺陷模型的启裂载荷随倾角β增加先减小后增加。

图5呈现了近球形缺陷模型启裂载荷随倾角β的变化规律。近球形缺陷模型启裂载荷随缺陷倾角β增加变化轻微。尺寸x=y=z= 1.06 mm的球形缺陷模型启裂载荷随倾角增加无变化。

图5 近球形缺陷模型启裂载荷随倾角 β 变化图Fig. 5 Variation of the crack initiation load of sub-spherical flaw models with inclination angle β

2.2 倾角对不同压缩度缺陷启裂载荷的影响

缺陷的形态从扁平的椭球形变化到近球形，其启裂行为会发生显著变化[21]。选取表1中编号为3，6，10，13，15和17的模型进一步分析椭球形缺陷的扁平程度对模型启裂载荷的影响。在这些模型中，x方向上的尺寸等于y方向上的尺寸。引入压缩度参数来描述缺陷的扁平程度。三维椭球形缺陷的扁平程度随着压缩度增大而增大。压缩度的计算公式为

式中B为缺陷压缩度，x为椭球形缺陷x方向上的半径，z为椭球形缺陷z方向上的半径。

图6呈现了不同压缩度的缺陷模型启裂载荷随缺陷倾角β的变化规律。当缺陷压缩度B=12.25时，模型启裂载荷随着缺陷倾角增加先快速减小后缓慢增加，在缺陷倾角β= 30°时，模型启裂载荷最小，缺陷最容易启裂。缺陷压缩度B=6.67，4.33和3.10的模型启裂载荷随缺陷倾角增加先缓慢减小后缓慢增加，同样缺陷倾角β=30°时，模型启裂载荷最小。当缺陷压缩度B=1.53和1.00时，模型启裂载荷随缺陷倾角增加变化较小。另外从图6中看到，曲线的凹陷程度随着缺陷压缩度增加而增加，这表明随着缺陷压缩度增加，缺陷倾角β对启裂载荷的影响逐渐增大。

图6 不同压缩度的缺陷模型启裂载荷随缺陷倾角 β 变化图Fig. 6 Variation of the crack initiation load with inclination angle β for the models with different compressibility

2.3 缺陷压缩度对不同倾角缺陷启裂载荷的影响

图7呈现了缺陷压缩度B对不同倾角缺陷启裂载荷影响的规律。当缺陷倾角β= 0°时，模型启裂载荷随缺陷压缩度B增加先减小后逐渐增加。当缺陷倾角β= 15°，30°，45°，60°和75°时，模型启裂载荷随压缩度B增大而减小，并且减小速率逐渐降低。

图7 不同倾角的缺陷模型启裂载荷随压缩度B变化图Fig. 7 Variation of the crack initiation load with compression degree B for the models with different inclination angles

3 结论

本文研究了体积相同的三维椭球形缺陷的启裂载荷随其倾角和压缩度的变化规律，得到以下结论。

（1）三维椭球形缺陷的倾角对其启裂载荷影响程度与缺陷压缩度有关。对于压缩度较小的近球形缺陷，其启裂载荷随倾角增加变化轻微。对于压缩度较大的三维椭球形缺陷，其启裂载荷随倾角增加先减小后增加，在倾角大约为30°时，启裂载荷最小。随着三维椭球形缺陷压缩度增加，其启裂载荷受倾角影响程度越大。

（2）三维椭球形缺陷的压缩度对其启裂载荷影响规律与缺陷倾角有关。对于平行加载方向的三维椭球形缺陷，即倾角为0°时，其启裂载荷随压缩度增加先减小后增加。然而当倾角为15°，30°，45°，60°和75°时，启裂载荷随压缩度增加单调递减，并且减小速率逐渐降低。