含三维内裂纹试件单轴拉伸裂纹扩展数值模拟研究

彭 博

(西南交通大学 土木工程学院， 四川 成都 611756)

1 研究背景

裂隙与缺陷广泛存在于岩土和水利工程当中，与边坡失稳、坝体垮塌以及结构的使用寿命息息相关[1-4]，尤其是在水利岩土工程中的深埋地下结构，裂隙与软弱夹层在高应力作用下将成为结构稳定的关键因素[5-6]，因此，对于岩土及水利工程中的裂纹扩展演化规律一直是学术界的研究热点，同时也是难点问题。

对于工程中的岩体来说，抗拉强度要远远小于抗压强度[7]，含裂纹的岩体在受拉情况下更容易导致灾难性的后果，因此研究裂纹在拉应力下的裂纹扩展演化规律显得很有必要，针对于此，国内外许多学者进行了大量的室内试验与数值模拟研究。试验方面：喻勇等[8]对弹性模量以及泊松比对巴西圆盘抗拉强度的影响进行了室内试验的分析；于庆磊等[9]对含中心穿透型裂纹的巴西圆盘进行了室内试验，研究了不同角度预制裂纹在轴压下的裂纹扩展演化规律；张绪涛等[10]利用改进的拉伸试验装置对砂岩、灰岩和花岗岩的抗拉强度进行了测定。数值模拟方面：鲁国富等[11]对蒙皮材料试件进行了双轴拉伸数值模拟研究；黎坤运等[12]利用扩展有限元对含不同角度穿透型裂纹进行了数值模拟，得到了裂纹扩展规律。然而上述研究仅仅针对二维裂纹或者三维穿透型裂纹，在实际情况中，裂隙往往以三维形态存在于岩体中，因此，对三维内裂纹在轴拉情况下的裂纹扩展规律的研究意义重大。

本文利用Franc3D软件，对于含不同角度预制内裂纹的标准立方体试件进行了单轴拉伸数值模拟研究，并将数值模拟结果与室内试验对比，验证了模拟的合理性，同时给出了裂纹前缘的应力强度因子的变化规律，为相应的工程实际提供了一定的参考。

2 计算理论

2.1 M积分计算应力强度因子

J积分可以表示为以下形式[13]：

(1)

式中：δij为克罗狄克符号;q为赋值函数;Γ为计算积分域。

应力强度因子KⅠ，KⅡ，KⅢ，与J积分的关系可以表示为以下形式：

(2)

式中：KⅠ、KⅡ、KⅢ分别为Ⅰ型、Ⅱ型、Ⅲ型应力强度因子；E为弹性模量；ν为泊松比。

应力、应变及位移关系可以表述为：

(3)

(4)

(5)

Ki=Ki(1)+Ki(2)

(6)

其中应力、应变以及位移上的(1)与(2)为两种独立的线弹性状态，将公式(6)代入公式(2)，可得：

J=J(1)+J(2)+M(1,2)

(7)

其中：

(8)

(9)

(10)

其中：

(11)

(12)

(13)

(14)

因此利用M积分便可计算出3个应力强度因子，本文所使用的Franc3D软件正是利用这样的原理，计算出裂纹尖端的应力强度因子的分布，同时可以依据应力强度因子判断出裂纹扩展的方向。

2.2 最大拉应力准则(MTS)

裂纹扩展准则采用最大拉应力准则(MTS)，裂纹将沿着最大周向应力σθmax方向进行扩展，而环向应力可以表示为：

(15)

(16)

3 计算模型与计算参数

建立标准立方体模型，试样设计图与模型网格图如图1所示。其中试件尺寸为标准立方体试件尺寸，长×宽×高为50 mm×50 mm×100 mm，插入三维内裂纹取文献[11]室内物理试验中的实际尺寸，即2a×2b=20 mm×10 mm的椭圆形三维内裂纹，材料的弹性模量为E=17.92 GPa，泊松比μ=0.192，密度ρ=2.3 g/cm3。

计算边界条件为：模型底部施加固定3个方向的位移约束，模型上部施加固定x与z方向的位移约束，同时y方向施加应力边界。

计算工况一共设定为5组工况，分别为预制裂纹与水平方向夹角α呈0°、15°、30°、45°、60°情况，由Franc3D软件依据最大拉应力准则计算出裂纹的扩展。

图1 计算模型图

4 计算结果分析

4.1 不同工况下的裂纹扩展

预制裂纹角度分别为0°、15°、30°、45°和60°工况下的裂纹扩展过程数值模拟结果见图2。

由图2可见：当预制裂纹角度为0°时，裂纹面承受拉应力，此时裂纹尖端Ⅰ型应力强度因子最大，而Ⅱ、Ⅲ型应力强度因子几乎为0(见4.2节应力强度因子分析)，预制裂纹呈现自相似扩展的规律，即裂纹沿着原来裂纹面进行扩展，最终裂纹贯穿裂纹面，试件破坏；当预制裂纹角度为15°时，预制裂纹尖端裂纹呈现“翼裂纹”状扩展，即裂纹前缘的扩展方向与原预制裂纹方向呈现一个较小的夹角。翼裂纹扩展与原预制裂纹方向呈18°夹角，最终裂纹扩展逐渐趋向于垂直最大主应力方向，随后裂纹贯穿试件表面，试件破坏；当预制裂纹角度为30°时，预制裂纹扩展形态与15°较为类似，即裂纹尖端呈现“翼裂纹”状扩展，而翼裂纹与原预制裂纹尖端呈现49°夹角，与预制裂纹角度为15°试件相比较大，最终裂纹扩展逐渐趋向于垂直最大主应力方向，随后裂纹贯穿试件表面，试件破坏；当预制裂纹角度为45°时，预制裂纹尖端仍然呈现“翼裂纹”状扩展，翼裂纹与原预制裂纹尖端呈现64°夹角，最终裂纹扩展逐渐趋向于垂直最大主应力方向，随后裂纹贯穿试件表面，试件破坏；当预制裂纹角度为60°时，翼裂纹与原预制裂纹方向角度更大，为80°，最终裂纹扩展逐渐趋向于垂直最大主应力方向，随后裂纹贯穿试件表面，试件破坏。

4.2 应力强度因子分析

为直观描述裂纹尖端的应力强度因子随预制裂纹角度的变化规律，定义相对应力强度因子不同计算方案下的裂纹尖端应力强度│Ki│/│Kmax│为相对应力强度因子，从而可以得到裂纹尖端的相对应力强度因子随裂纹尖端相对位置的变化规律，预制裂纹的前缘相对位置定义为以裂纹右端为起点，逆时针旋转一周回到起点一周长定义为1，如图3所示。

由图3(a)可见，随着预制裂纹与水平方向的角度增大，相对Ⅰ型应力强度因子逐渐减小，可见预制裂纹角度越大，裂纹越难扩展，同时对于预制裂纹从相对位置0到0.5的距离上，相对Ⅰ型应力强度因子呈现先增大后减小的规律，可见预制裂纹短轴顶点上的相对Ⅰ型应力强度因子较其他方向要大，因而试件最终破坏是以裂纹先贯穿试件前后表面，进而使试件破坏。由图3(b)可见，随着预制裂纹角度的增大，相对Ⅱ型应力强度因子整体上先增大后减小，预制裂纹30°与预制裂纹60°的相对Ⅱ型应力强度因子几乎一致。值得注意的是，预制裂纹角度为0°时相对Ⅱ型应力强度因子为0，表明裂纹面处于纯拉状态，与数值模拟中裂纹沿原裂纹面扩展的规律较为一致。对于同一种预制裂纹角度来说，预制裂纹的左端与右端(对应的位置为0、0.5和1.0处)相对Ⅱ型应力强度因子达到最大，而预制裂纹短轴处顶点(对应的位置为0.25和0.75处)相对Ⅱ型应力强度因子达到最小为0。

5 与室内试验结果的对比

为验证本文数值模拟的合理性，取李术才等[14]、张波等[15]和张敦福等[16]的试验结果进行对比，不同裂纹角度的裂纹扩展形态与本文的数值模拟结果对比见图4～6。

图2 不同角度预制裂纹试件裂纹扩展过程

图3 裂纹前缘相对应力强度因子

图4 预制裂纹角度为0°裂纹形态对比

图5 预制裂纹角度为30°裂纹形态对比

图6 预制裂纹角度为45°裂纹形态对比

由图4可见，对于预制裂纹角度为0°的情况，由于相对Ⅱ型应力强度因子为0，裂纹面处于纯拉状态，则裂纹沿着原裂纹面呈自相似扩展，数值模拟结果也验证了室内试验的结论。

由图5和6可见，由于Ⅱ型应力强度因子的存在，使得裂纹的扩展偏离原预制裂纹面，室内试验的断口呈现出明显的翼裂纹扩展，这与数值模拟的结果较为吻合。

6 结 论

本文利用Franc3D软件，对含不同角度(0°、15°、30°、45°和60°)三维内裂纹标准立方体试件单轴拉伸下的裂纹扩展过程、扩展形态以及应力强度因子进行了数值模拟分析，并将数值模拟结果与相关参考文献中的室内试验结果进行了对比分析，得出了以下结论:

(1)当预制裂纹角度为0°时，预制裂纹呈现自相似扩展，当预制裂纹为其他角度时，预制裂纹尖端出现翼裂纹扩展，预制裂纹角度为15°、30°、45°和60°时，翼裂纹分别与原预制裂纹角度呈18°、49°、64°与80°。

(2)随着预制裂纹与水平方向的角度增大，相对Ⅰ型应力强度因子逐渐减小，对于同一预制裂纹角度来说，预制裂纹短轴顶点上的相对Ⅰ型应力强度因子较其他方向要大。

(3)随着预制裂纹角度的增大，相对Ⅱ型应力强度因子整体上先增大后减小，对于同一种预制裂纹角度来说，预制裂纹的左端与右端(对应于位置为0、0.5和1处)相对Ⅱ型应力强度因子达到最大，而预制裂纹短轴处顶点(对应位置为0.25和0.75处)相对Ⅱ型应力强度因子达到最小为0。

(4)将数值模拟结果与室内试验结果进行对比，证明了本文数值模拟的合理性。