基于SHPB 的UHPC 冲击试验径向惯性效应分析*

2019-11-16贾永峰

爆炸与冲击 2019年10期

任亮，何瑜，王凯，李周，贾永峰

（1. 华东交通大学土建学院，江西南昌 330013；2. 中建钢构有限公司，广东深圳 518000；3. 江西省萍乡市公路管理局，江西萍乡 337000）

超高性能混凝土(ultra high performance concrete，UHPC)是20 世纪90 年代由法国Bouygues 实验室提出的一种新型的水泥基材料，具有损伤容限高、断裂韧性大、耐久性能优异等特点，可广泛应用于海上结构、地下空间、核废料容器和核反应堆防护罩等特殊工程和国防军事工程中。为适应其对抗冲击和抗爆性能的需求，许多学者通过分离式霍普金森压杆(split Hopkinson pressure bar，SHPB)对其冲击压缩性能开展试验研究，结果表明UHPC 动态力学性能指标随着应变率的增加而提高，UHPC 材料具有明显的动态增强效应，但在试验加载过程中应力波的弥散效应、试件与杆件端面摩擦效应和径向惯性效应将对试件应力平衡和一维应力波传播产生重要影响[1]，其中试件与杆件的端面摩擦效应，可通过端面涂抹凡士林等润滑剂减少其对试验结果的影响[2]。对于应力波的弥散效应，已有的研究表明可通过在撞击杆和入射杆之间设置波形整形器可来克服。其中：Hassan 等[3]通过冲击压缩试验研究指出，铝相对于铜更适合作为UHPC 材料的整形器；Xu 等[4]基于数值模拟，进一步指出整形器直径不宜大于杆件直径的0.4 倍，长径比不宜大于0.2。对于径向惯性效应，学者们通过理论分析和试验研究，提出了不同的普通混凝土试件惯性效应修正公式[5-7]，指出试件截面尺寸、长径比和恒应变率加载是影响惯性效应的重要因素。但对于UHPC 材料，目前少有关于减少径向惯性效应的报道，如参考传统混凝土材料的处理方法，其适用性和可靠性仍有待深入研究。

针对上述问题，本文中应用大型有限元分析程序LS-DYNA，对SHPB 试验中减少径向惯性效应的有效措施开展研究。通过对程序中KCC 损伤模型控制参数取值优化，拟合UHPC 材料动态性能，建立UHPC 试件应用SHPB 技术的数值仿真模型，并与实验验证。在此基础上，开展不同UHPC 试件直径、长径比以及有无整形器下的参数分析，探讨其对SHPB 试验中试件径向惯性效应的影响。

1 SHPB 装置

式中：c0表示应力波在压杆中的弹性波速，L表示试件的长度，A和As分别表示压杆和试件的截面面积，E表示压杆的弹性模量，εR(t)和εT(t)分别表示由应变片测得的反射波和透射波信号。

图1 分离式霍普金森压杆试验装置Fig. 1 Test device of split Hopkinson pressure bar

2 径向惯性效应评价指标

为使式(1)～(3)成立，需要保证一维应力传播和试件应力平衡假设成立，但由于试验过程中端面摩擦效应、惯性效应和波形弥散等因素的影响，易导致应力波偏离一维传播和试件应力不均匀。考虑到本文中主要研究径向惯性效应对SHPB 试验加载过程的影响，在数值模拟中拟通过合理设置摩擦系数消除端面摩擦效应的影响，并在撞击杆和入射杆之间增设整形器消除波形弥散效应。为评价径向惯性效应的影响，引入应力平衡因子和应力状态判定系数来衡量，应力平衡因子δ 的表达式为：

为进一步判定试件加载过程中的应力状态，采用试件主方向应力与侧向应力的比值Rσ和主方向应变与侧向应变Rε作为双重应力状态判定系数，相应的表达式为：

式中：σx、σy、σz、εx、εy、εz分别表示沿x轴、y轴和z轴方向的主应力和主应变，N表示选取的单元数量。考虑到一维应力状态时试件应力-应变需满足以下条件：

式中：E表示弹性模量， ν 表示泊松比。当UHPC 泊松比按0.2 取值时[4]，根据式(6)、(9)～(11)，可知满足一维应力状态的应力状态判定系数Rε取值约为5；而对于应力状态判定系数Rσ，根据文献[8]中的建议，认为主轴应力与非主轴应力比大于10 时能近似认为满足一维应力状态，为此本文中Rσ按大于10 取值。考虑到惯性效应对试样径向约束作用，同时将结合对不同试样动态增强因子(dynamic increase factor，DIF)的影响进行判定。

3 数值模拟

受试验方法、加载条件和试件数量等的限制，试验研究并不能完全反应UHPC 材料的冲击压缩性能，为此，开展数值模拟成为有效的补充。目前ABAQUS、AUTODYNA、MSC/DYTRAN、LS-DYNA 等均可应用于混凝土材料的非线性动力有限元分析，其中LS-DYNA 软件提供了Holmquist-Johnson-Cook(HJC)、KCC 等损伤模型，可模拟材料从加载至破坏全过程力学行为，而被广泛应用于混凝土材料SHPB 试验数值模拟中[1,9-10]。与HJC 模型相比，KCC 模型能够准确地捕捉材料在冲击荷载作用下的非线性行为，更适合中高等应变率下混凝土材料的动态力学性能模拟[11]。为此，本文中选用KCC 模型，但考虑到现阶段UHPC 动态损伤模型少见相关的文献报道，因此如何基于KCC 损伤模型中的输入参数，优化拟合UHPC 的动态损伤行为是开展数值模拟的关键。

3.1 KCC 模型描述

在KCC 模型中，混凝土材料复杂的力学行为通过定义的3 个剪切破坏面(初始屈服面、最大屈服面和残余应力面)来控制，图2 为对应的3 个剪切面破坏面和材料单轴本构关系。其中，Δσ 和p分别表示主应力差和静水压力；fc′表示静态抗压强度；ft表示拉伸强度；P1 表示屈服点，P2 表示峰值强度点，P3 表示残余强度点。

上述3 个剪切面对应的表达式为：

式中：Δσy表示初始屈服强度；Δσm表示最大屈服强度；Δσr表示等效残余强度；a0y、a1y、a2y、a0、a1、a2、a1f、a2f为常量，一般通过单轴和三轴试验确定。

图2 KCC 模型的剪切破坏面及单轴本构关系Fig. 2 Shear failure surfaces and uniaxial stress-strain relation in KCC model

在荷载作用下，材料应力先上升到屈服点P1，到达初始屈服面；随后材料进入强化阶段，上升到峰值强度点P2，到达最大屈服面；最后材料发生破坏，下降到残余强度点P3，到达残余应力面，此时材料完全破坏。当材料应力介于3 个剪切破坏面之间时，当前应力为：

式中：λ 和η(λ)分别表示等效塑性应变和损伤函数；λm表示等效塑性应变峰值。λ 从0 增大到λm时，表示材料进入屈服强化阶段，此时η(λ)由0 增加到1，但随着λ 值的继续增大，η(λ)值逐渐下降为0，表示材料进入破坏阶段。在KCC 模型中η(λ)和λ 可以数组形式输入。

与此同时，KCC 模型在使用时需结合状态方程来描述混凝土的抗压性能，KCC 模型使用Tabulated-Compaction 模型作为状态方程。相应表达式为：

式中：εV表示体积应变，γ 是温度常数，C(εV)为εV对应的体积压力值，T(εV)为εV对应的温度值，E0为初始体积内能，p为加载压力。

上述状态方程中，加载压力随体积应变的增大单向递增；当状态方程出现卸载时，相应的卸载刚度定义为卸载体积模量K，卸载以截止压力作为终点；再加载路径沿卸载曲线加载至卸载开始点，然后按状态方程骨架曲线继续加载。在KCC 模型中压力p、体积应变εV和卸载体积模量K需以数组形式输入。

3.2 KCC 模型参数校准拟合UHPC 动态性能

对于传统混凝土，KCC 模型在输入对应的静态抗压强度后，可自动生成上述参数来反映其动态力学性能。UHPC 作为一种新型水泥基材料，其力学性能与传统混凝土存在显著的差异，如采用KCC 模型自动生成的参数将产生较大的偏差，为准确模拟UHPC 的损伤行为，需通过修正控制参数，对KCC 模型中进行校准，相应的步骤如下。

3.2.1拟合不同围压下UHPC受压本构关系

应用KCC 模型拟合UHPC 动态力学性能的关键是如何确定UHPC 在不同围压下的受压本构关系，以确定控制剪切破坏面的8 个主要参数。目前UHPC 在不同围压下的受压本构关系少有相关文献报道，本文中参照文献[12]提出的不同围压下混凝土受压本构关系，代入UHPC 材料参数对其进行拟合，其中上升段为：

式中：f和f0分别表示应力和考虑围压下的峰值应力；A、B、C、D表达式详见文献[12]；X表示应变ε 与考虑围压下峰值应变ε0的比值。其中，ε0为：

式中：fc′分别表示单轴抗压强度；εc表示单轴抗压峰值应变；fr表示三轴试验中对应的围压值。

不同围压下抗压本构下降段为：

式中：frp表示受压残余应力；fic和εi分别表示本构曲线下降分支上拐点处的应力和应变；frp、fic和εi的表达式详见文献[12]。

以静态抗压强度为200 MPa 的UHPC 为例[3]，取密度为2 291 kg/m3、单轴抗压强度为200 MPa、单轴抗压峰值应变为0.004、单轴抗拉强度为7.89 MPa，可得到UHPC 单轴受压和不同围压下的受压本构曲线如图3～4 所示，为便于对比，在图3 中加入了UHPC 单轴受压试验曲线和未修正的KCC 模型生成曲线。

图3 UHPC 本构关系拟合值与试验值对比Fig. 3 Comparison of constitutive relation of UHPC between fitted and experimental values

图4 不同围压下的应力应变曲线Fig. 4 Stress-strain curves under different confining pressures

从图3 可以看出，未修正的KCC 模型生成曲线与试验值偏差较大，且在下降段由于普通混凝土较小的压应变极值导致快速下降，而修正后KCC 模型生成曲线与试验值吻合较好，表明修正后的KCC 模型能较好拟合UHPC 单轴受压本构曲线。

从图4 可以看出，修正后KCC 模型生成的有围压本构曲线在下降段通过拐点后均出现了近似水平的平台，这是由于式(19)中当ε/ε0比值达到一定程度后，第2 项趋近于零，导致应力变化幅度减少，从而出现图中近似水平的平台。

图5 剪切破坏面拟合Fig. 5 Fitted shear failure surfaces

3.2.2 确定剪切破坏面

根据式(1)确定的不同围压下的UHPC 受压本构曲线，可分别获得等效屈服强度点、等效破坏强度点和等效残余强度点对应的轴向应力σ1和侧向围压σ2、σ3，进而得到剪切破坏面的主应力差Δσ 和静水压力p，通过数值拟合即可得到3 个剪切破坏面，图5 为静态抗压强度200 MPa 的UHPC 数值拟合得到的剪切破坏面。

3.2.3 确定剪切破坏面参数

根据拟合的剪切破坏面曲线，结合式(12)～(14)，可确定剪切破坏面8 个控制参数，见表1。

表1 剪切破坏面确定参数Table 1 Parameters for determining shear surface

3.2.4 确定损伤方程和状态方程

KCC 模型参数校准拟合UHPC 材料动态性能时，损伤方程中参数η、λ 及状态方程中p、εV和K参照文献[13-14]取值。

3.3 数值仿真模型

应用大型有限元程序LS-DANA 建立SHPB 数值仿真模型。模型中整形器、杆件和UHPC 试件均采用三维solid 164 实体单元，单元横向网格均按5 mm 进行网格划分，纵向网格试件和整形器按1 mm 划分，杆件按20 mm 划分。为防止负体积的产生，在杆件与试件、杆件与整形器接触部分单元，杆件纵向网格按1 mm 划分。为消除端面摩擦效应的影响，所有界面接触均采用自动面面接触，且不激活静/动态摩擦系数。考虑到整形器和试件在试验过程中都将进入塑性，甚至发生破坏，为此整形器选用各向同性、随动硬化或两者混合的Plastic-Kinematic 模型，UHPC 试件选用参数校准后的KCC 动态损伤模型。为节省计算时间和占用的空间，应用对称性仅选取截面的四分之一进行建模，图6 为建立的分析模型。

图6 SHPB 数值仿真模型Fig. 6 Numerical model of SHPB

3.4 试验验证

Hassan 等[3,15]分别对无整形器和有整形器的UHPC 试件开展了SHPB 试验。其中，整形器材料选用铝和铜，厚度均为2 mm，直径均为10 mm，撞击杆速度铝质整形器为10 m/s，铜质整形器为11.9 m/s。表2列出了铝和铜材对应的材料参数；SHPB 设备中撞击杆、入射杆和透射杆特征参数见表3；UHPC 试件静态抗压强度为200 MPa，长度和直径分别为9.2 mm 和23.8 mm。应用上述建立的SHPB 数值仿真模型，本文中对文献[3-15]中实测的无整形器、铝质整形器、铜质整形器3 种情况下的入射波、反射波和透射波分别进行数值模拟，图7 为实测波形和计算波形的对比。

表2 整形器材料参数Table 2 Material parameters of pulse shapers

表3 杆件的特征参数Table 3 Characteristic parameters of bars

图7 实测波形与计算波形对比Fig. 7 Comparison of waveforms between experiment and calculation

从图7 可以看出，SHPB 试验实测波形与计算波形大致吻合。与无整形器波形相比，SHPB 试验在使用铝质和铜质整形器后，入射波均出现了明显的上升平台，有利于UHPC 试件到达应力平衡，且透射波取得峰值前，反射波出现了近似的水平段（恒应变率加载平台）。计算波形中较好地反映了波形的这种变化，说明本文中提出的数值仿真模型有较好的精准性。

4 参数分析

基于上述数值分析模型，本文中选取试件的直径、长径比和整形器等关键参数，对UHPC 试件在SHPB 试验中应力均匀性和单轴应力传播进行分析，探讨上述参数对径向惯性效应的影响。分析时SHPB 杆件尺寸、UHPC 静态抗压强度与文献[3]一致，试件直径和长径比分析时波形整形器参照文献[4]选用直径8 mm、厚度2 mm 的铝。

4.1 试件直径

为探讨试件直径对试验结果的影响，分别选取试件直径d和杆件直径DB的比值为0.8、0.85、0.9、0.95、1.0 和1.05 情况下的透射波波形、应力平衡因子、DIF 和应力状态判定系数进行分析。为剔除应变率的影响，分析时通过改变撞击杆速度将应变率控制在约90 s-1，相应的撞击杆速度分别为9.0、9.2、10.0、10.1、10.5 和11.1 m/s；分析时试件的长径比参照文献[3]按0.4 取值。图8 为对应的分析结果。

图8 试件直径对径向惯性效应影响对比Fig. 8 Comparison of the influence of specimen diameters on inertia effect

从图8(a)～(b)可以看出，在加载应变率基本保持不变的情况下，随着试件与杆件直径比的增加，在透射波达到峰值之前，应力平衡因子δ 均能持续趋近于零，试件破坏前能较好地实现应力平衡。

从图8(c)～(e)中可以看出，当d/DB小于1.0 时，Rσ均大于10，且DIF 变化不大，Rε经过短时间振荡后趋近于5 并在试件破坏之后发生偏离，表明此时试件在破坏前达到了一维应力状态；当d/DB达到1.0 时，Rσ明显小于10，而Rε在试件破坏前明显大于5，且DIF 相对于d/DB为0.95 的试件显著提高，试件不再保持一维应力状态。在试件保持恒应变率加载和应力平衡的前提下，小直径的试件较大直径试件能更好地维持一维应力状态，受惯性效应影响更小，建议试件d/DB不大于0.95。考虑到小直径试件的波阻抗相差较大且制作相对困难，建议试件d/DB不低于0.9。

图9 试件长径比对惯性效应影响对比Fig. 9 Comparison of the influence of specimen aspect ratios on inertia effect

图10 整形器设置对惯性效应影响对比Fig. 10 Comparison of the influence of pulse shaper on inertia effect

4.2 试件长径比

为探讨试件长径比对试验结果的影响，分别选取试件长度L和直径d的比值为0.30、0.35、0.40、0.45、0.50 和0.55 情况下的透射波波形、应力平衡因子、DIF 和应力状态系数进行分析，如图9 所示。分析时试件直径按杆件直径的0.95 倍取值，撞击杆速度为10 m/s。

从图9(a)～(b)可以看出，不同长径比试件的透射波在达到峰值之前，应力平衡因子δ 均具有明显趋近于零的平台，表明试件直径与杆件直径在匹配的情况下，试件长径比对试件加载过程中的应力平衡影响较小，但是随着试件长径比的减少，透射波的上升沿持续时间逐渐减少，试件破坏的时间不断提前。考虑到试件长径比过小不仅难于制作，而且UHPC 材料中钢纤维均匀性难以控制，建议试件长径比不宜低于0.35。

从图9(c)～(e)可以看出，当试件长径比小于0.5 时，应力状态判定系数Rσ和Rε虽然在前期存在一定的震荡，但在破坏前均能满足一维应力传播判定准则，且DIF 相对平稳；当长径比达到0.5 后，Rσ明显小于10，Rε在破坏前一直保持震荡，DIF 相对于长径比为0.45 的试件显著提高，试件在破坏前应力传播难以满足一维应力状态。为确保试验的有效性，建议试件长径比不超过0.45。

4.3 整形器

已有的研究表明，恒应变率加载是减少SHPB 冲击实验中径向惯性效应的有效措施[6]，而Xu 等[4]指出，整形器是实现恒应变率加载的有效途径。为此，对不同撞击杆速度下无整形器、使用铝质和铜质整形器下的透射波波形、应力平衡因子和应力状态系数进行分析，图10 为对应的分析结果。其中，试件直径参考前述分析结果，按杆件直径的0.95 倍取值（2.41 mm），长径比按0.4 取值。

从图10(c)～(h)中可以看出：在未使用波形整形器的UHPC 试件SHPB 试验中，即使试件直径和长径比满足上述取值，在不同的撞击杆速度下，试件在破坏前均难以实现应力平衡和一维应力传播；而在使用铝质和铜质整形器后，在不同的撞击杆速度下，试件在达到破坏时刻前，应力平衡因子δ 均具有明显趋近于零的平台，且应力状态判定系数Rσ和Rε均能满足一维应力传播判定准则。这表明，使用波形整形器实现恒应变率加载是UHPC 材料应用SHPB 开展冲击试验减少径向惯性效应的重要前提。