基于HJC模型高温后大理岩SHPB实验数值模拟研究

李 睿，刘 磊，2，张志华，王 远

(1.昆明理工大学国土资源工程学院，昆明650093；2.云南省中-德蓝色矿山与特殊地下空间开发利用重点实验室，昆明650093；3.中铁十九局集团轨道交通工程有限公司，广东 深圳518000)

在地下空间的开发随着浅部资源的日益减少和枯竭而逐渐转入岩体的深部过程中，随着采掘深度的不断增加，地温也在不断升高[1]。除了高温外，在深部地下工程施工的过程中，爆破施工、大型机械扰动等动载荷对围岩稳定性产生不同程度的影响，为了确保深部地下空间开发施工的安全及地下工程的稳定性，需要对岩石受高温和动载荷影响下物理力学性质进行探究。

近来，很多学者对不同岩石受高温影响的动态力学特性进行了不同方式的实验研究[2]，平琦等[3]通过对不同温度条件下砂岩试件进行动态力学实验，探究温度对砂岩试件的动态抗压强度、动态峰值应变、弹性模量等动态力学指标的影响；李明[4]对不同温度等级及不同冲击速率下煤系砂岩动态力学特性，并运用射线衍射分析及等实验手段对破坏后砂岩的细观力学机制进行研究，分析温度以及冲击载荷对煤系砂岩力学性质和行为的影响。

随着有限元分析技术发展,运用ABAQUS、ANSYS、MSC等有限元分析软件对爆破及岩石冲击动力学问题进行数值模拟分析并对实验研究结果与数值模拟结果进行对比分析已成为一种成熟的研究思路[5-7]。由Holmquist、Johnson和Cook对JC模型进行岩石材料适用性修改的HOLMQUIST-JOHNSON-CONCRETE(HJC)损伤本构模型[8]最初主要应用于混凝土的高应变率、大变形问题的研究，随着岩石动力学数值模拟研究技术的发展和对HJC材料模型研究的深入，HJC材料模型被不断的完善，巫绪涛等[9]通过数值模拟结果与实验结果相对比探究HJC材料模型的参数的确定方法；方秦等[10]将HJC材料模型推广到岩石材料中，探究相应参数确定的方法，并通过大理岩的弹体侵彻实验研究及数值模拟研究进行了验证。

综上所述，尽管目前，对于涉及温度影响的岩石材料动力学性能的研究已经取得了相当丰富的成果，但是对相关数值模拟的研究仍然鲜见报道，本文基于ANSYS/LS-DYNA数值模拟有限元分析软件和HJC材料模型，确定经不同温度等级处理并自然冷却后大理岩HJC本构模型参数，模拟不同温度等级处理后大理岩试件单轴动态压缩实验过程，并与高温后大理岩SHPB实验结果进行对比分析。

1 高温后大理岩SHPB实验

1.1 试件制备与实验装置

根据相关规范及文献，将大理岩加工φ50 mm×40 mm圆柱体并利用陶瓷马弗炉对加工好的标准大理岩试件进行不同温度等级(100～800 ℃)高温处理并自然冷却至室温，控制升温速率为5 ℃/min。

利用SHPB测试系统进行高温后大理岩动态冲击力学实验，该系统主要由加载驱动装置、荷载传递装置和测试系统3个部分组成。其中，入射杆和透射杆的相关参数如表1所示。

表1 杆件基本参数

1.2 实验方案与结果

分别对常温(25 ℃)及不同温度等级(100、200、400、600、800 ℃)处理大理岩试件进行恒定冲击气压单轴动态冲击压缩实验，主要实验结果如表2所示。

表2 不同温度等级大理岩动态压缩实验结果

常温及经不同温度加热并冷却大理试件动态应力应变曲线如图2所示，由图可以看出，高温后(200～800 ℃)大理岩单轴动态应力应变曲线与常温状态下单轴动态应力应变曲线特征基本相似，并且随着温度等级的升高，大理岩脆性增大，其大致可分为压密阶段、弹性阶段、破坏阶段3个阶段。

图1 大理岩试件

图2 应力应变

大理岩单轴动态抗压强度与温度关系如图3所示，由图可以看出，随着温度的增加，25、200、400、600、800 ℃温度等级下大理岩试件单轴动态抗压强度呈现先不断增大后逐渐减小的现象，可见温度对大理岩试件单轴动态抗压强度有着显著的影响，并且这种影响与温度有着密切的联系。在25～400 ℃范围内，随着温度作用于大理岩试件，温度表现出对大理岩试件动态抗压强度的强化作用，并且这种强化作用随着温度的升高不断增强;当温度升高至400 ℃时，温度对大理岩试件动态抗压强度的增强作用达到最大;在400～600 ℃范围内，温度对大理岩单轴动态抗压强度仍有强化的作用(相对于常温25 ℃时)，但是这种强化作用随着温度的升高而不断减小;当温度大于600 ℃时，温度对大理岩动态抗压强度产生一定弱化作用(相对于常温25 ℃时)且这种弱化作用随着温度的升高不断增大。

图3 峰值应力变化

不同温度等级典型大理岩破坏形态如图4所示，从图可以看出，大理岩破碎后的细粒量随温度的升高而增多。这是由于在高温作用下，岩石内部的水分蒸发，矿物成分改变，导致黏结力逐渐变小，试样破碎后的碎块尺寸由粗粒向细粒变化，试样更加破碎，进而导致随着温度等级的升高大理岩试件破碎程度向更剧烈的趋势发展。

图4 不同温度等级典型大理岩破坏形态

2 模型建立与参数确定

2.1 有限元模型建立

选取ANSYS/LS-DYNA数值模拟软件Solid164模块建立大理岩SHPB实验有限元模型(见图5)，用直接在入射杆入射端加载半正弦波的方式代替设置子弹冲击。对于杆件的精度要求不高，设置入射杆与透射杆分别沿径向网格划分18份、沿轴向网格划分200份；试件沿径向网格划分50份，沿轴向网格划分35份。

图5 有限元模型

定义试件与杆件接触为面面自动接触,通过删除单元的形式体现大理岩试件的失效，采用HJC材料模型自带的失效准则定义大理岩试件在单轴动态压缩下的失效。

2.2 杆件基本物理力学性质

在SHPB实验模拟中，不仅要考虑应力波在杆件中的传播，同时需要考虑杆件因受到冲击而发生的的变形，选择ELASTIC材料模型作为杆件的材料模型。杆件的材质为40Cr合金钢，参照表1中SHPB测试系统杆件参数在软件中对杆件进行设置。

2.3 常温大理岩HJC模型参数

HJC 模型主要由强度模型、损伤模型和状态方程3部分组成，共涉及21个不同的参数，按照性质共分为基本物理力学参数、强度参数、损伤参数、压力参数4类。借鉴文献[10-13]所介绍的方法，采取在直接对大理岩试件进行动、静态物理力学实验的基础上，同时参考工程经验公式并不断试算的方法确定常温状态下大理岩HJC材料模型参数。

1)常温下大理岩基本物理力学参数。根据相关实验规范，利用YAW-2000 型微机控制全自动压力实验机对常温状态下大理岩试件进行基本的岩石物理力学实验(准静态单轴抗压、巴西圆盘劈裂等),并通过式(1)、式(2)进行计算，确定常温状态下大理岩试件基本物理力学参数(见表3)。

表3 常温状态下大理岩基本物理力学参数

(1)

(2)

式中：G为剪切模量；E为弹性模量；K为体积模量；μ为泊松比。

2)强度参数。HJC材料模型强度模型如式(3)所示，其中D为损伤变量，其取值范围为0～1，表示岩石试件损伤程度，在自然状态下默认大理岩试件为无损伤(D=0)且忽略应变率与温度对大理岩试件的影响后可改写如式(4)所示：

(3)

σ*=A+Bp*N

(4)

式中：A为黏性强度系数；B为压力硬化系数；N为标准化静水压力的指数；p*为无量纲的静水压力；σ*为无量纲的等效应力。

(5)

(6)

(7)

对大理岩试件进行不同围压等级下的静态三轴压缩实验，并结合M-C准则[14]对实验结果进行拟合，拟合曲线如图6所示，得到拟合公式(8)，将拟合结果代入式(9)可得黏聚力C为18.73 MPa，可得到黏性强度系数A=0.25。

图6 M-C准则拟合

参照文献[11]的方法对无量纲的等效应力σ*和标准化静水压力p*进行拟合，拟合曲线如图7所示，并得到拟合公式(10)，得到压力硬化系数B=1.8，标准化静水压力指数N=0.83，Smax=5，并确定应变率系数C=0.001 84。

图7 屈服面拟合

σ1=4.21σ3+76.81 (R2=0.882)

(8)

(9)

σ*=0.25+1.8p*0.83(R2=0.955)

(10)

(11)

D1(p*+T*)D2≥EFmin

(12)

4)压力参数。HJC材料模型压力参数主要与HJC材料模型的状态方程有关(见式13、式14)，参考文献[10]介绍的方法,按照式(15)、式(16)分别计算pcrush、μcrush并得到K1、K2、K3及μlock。整理常温状态下大理岩HJC材料模型参数如表4所示。

表4 常温状态下大理岩HJC材料模型参数

(13)

(14)

(15)

(16)

2.4 高温后大理岩HJC参数确定

相对于实时高温需要对岩石试件进行温度场和动态冲击加载的耦合计算，高温后的SHPB实验过程中并没有热量的变化，温度对大理岩试件造成的温度损伤直接改变大理岩部分物理力学性质，通过材料模型参数的变化体现高温对大理岩试件的损伤，结合实验结果并参照常温状态下大理岩HJC材料模型参数确定方法确定不同等级高温处理并自然冷却至常温大理岩HJC材料模型参数如表5、表6所示。

表5 不同等级高温后大理岩HJC材料模型部分参数

3 数值模拟结果

3.1 应力应变

对所提取的波形利用三波法[18]进行计算并重构应力-应变曲线如图8所示,峰值应力变化趋势如图9所示。由图8、图9可知，随着温度等级的提升，大理岩动态抗压强度随着温度等级的提升并不是简单的上升或是下降而是先上升并在接近阈值时开始下降，当温度等级超过600 ℃后，大理岩动态抗压强度开始低于常温时动态抗压强度，这与室内实验过程中大理岩动态应力应变曲线及所体现性质(见图2、图3)总体趋势相一致。

图8 数值模拟应力应变

图9 数值模拟峰值应力变化

3.2 破坏过程与破坏形态

HJC材料模型自带的失效准则对满足失效条件的单元进行了删除从而形成裂隙，裂隙发育最终导致试件破坏。温度等级为200 ℃大理岩试件单轴冲击压缩模拟破坏过程如图10所示。

图10 Temp=200 ℃ 破坏过程

从图10可以看出，在动态冲击压缩过程中，大理岩试件表面应力经历了均匀-集中-分散的过程，呈现裂隙由边缘向中心发展的趋势，表现为先是以劈裂拉伸破坏产生裂隙进而进一步完全破坏，这与室内实验过程中以压碎破坏为主的破坏形态有一定差异，有学者[19]认为造成这种现象的原因是因为岩石试件抗拉强度远小于抗压强度，当应力波传递到大理岩试件表面时，先发生了劈裂拉伸破坏，相对于数值模拟对破坏形态的观测可以任意缩短观测时间区间，室内实验持续时间短，不容易观察试件完整破坏过程，只能通过最后的破坏形态分析破坏的过程。

不同温度等级大理岩试件模拟破坏结果如图11所示，通过与室内实验结果相比较，数值模拟破坏形态与实验破坏形态体现的规律具有一致性，温度的损伤效应对于大理岩试件动态冲击破坏形态的影响体现在破碎后细粒量，随着处理温度的升高，大理岩破碎后细粒也不断增多。

图11 不同温度等级数值模拟破坏结果

4 结论

1)通过与室内实验结果对比，利用ANSYS/LS-DYNA数值模拟软件对经不同温度等级处理并自然冷却后大理岩进行SHPB实验模拟所得实验结果与室内实验结果具有一致性。

2)HJC材料模型参数可以有效地表征经高温处理并自然冷却大理岩试件各项物理力学性质，在确定HJC材料模型参数过程中，“参考-计算-试算-调整”的方法具有一定的可参考性，通过直接修改HJC材料模型参数的方法可以有效地对经高温处理并自然冷却的大理岩试件物理力学参数描述。

3)大理岩试件在SHPB试件过程中，呈现裂隙由边缘向中心发展的趋势，表现为先以劈裂拉伸破坏产生裂隙进而进一步完全破坏。温度的损伤效应对于大理岩试件动态冲击破坏形态的影响体现在破碎后细粒量，随着处理温度的升高，大理岩破碎后细粒也不断增多。