混凝土材料SHPB主动围压实验的数值模拟

张 华， 郜余伟， 武守锋

（河海大学 土木与交通学院，江苏 南京 210098）

在自然界及实际工程中有相当一部分混凝土工程在承受动荷载作用之前，已经处于三向应力状态，而大量的研究表明，脆性材料在围压作用下会表现出特殊的力学行为。在高应变率下，随着围压的增大，围压和应变率的耦合作用将使得混凝土材料的应力峰值增加，破坏模式由脆性破坏转变为延性破坏［1－2］。因此，有必要选取合适的方法对脆性材料在三向应力状态下的动态力学性能进行研究，目前较为多用的实验方法是SHPB主动围压实验。随着SHPB研究的深入，主动围压技术也得以不断发展，中国科学技术大学［3］、空军工程大学（与洛阳立特公司合作）［4］、宁波大学［5］都推出了主动围压装置。但是霍普金森压杆及与其配套的主动围压装置存在以下不足：

（1）对于脆性材料，受应力不均匀性及波形弥散效应的影响，SHPB装置数据测量困难、实验精度低；SHPB装置很难达到每秒104～106范围内更高应变率的加载；由于实验过程瞬间结束，很难观测材料内部体貌的演变过程。

（2）主动围压装置实验技术和测试手段较为复杂；实验过程中试件的变形将导致围压值的波动，影响实验精度；目前主动围压装置，多为二维加载，只是近似模拟工程实际中的三向受力。

数值模拟技术能够提供详细而且丰富的信息量，因此可以直接获取试件各单元的应力、应变、位移等信息，能有效分析混凝土在动态荷载下的力学响应，且数值模拟可以避免电磁干扰、端面摩擦效应、试件平行度公差等外界因素的影响，进而提高了测试精度。通过在试件周围施加恒定压力可以对主动围压试验进行模拟，且压力值不受试件变形的影响，该方法有望进一步发展，对实验中难以达到的三维围压进行模拟。

综合以上分析，数值模拟手段可以弥补霍普金森压杆实验中的不足，因此，本文利用LS－DYNA有限元软件和 HJC（Holmquist－Johnson－Concrete，简称HJC）混凝土动态本构模型对Φ74mm分离式变截面霍普金森压杆主动围压实验进行模拟。

1 SHPB实验装置

1.1 Φ74mm直锥变截面霍普金森压杆

Φ74mm直锥变截面霍普金森压杆装置，是由中国科学技术大学胡时胜教授研制［6］。入射杆采用锥形过渡的变截面形式，直径由37mm增大到74mm，可以利用截面较小的子弹驱动装置对大尺寸试样进行SHPB实验。其基本实验原理与常规等截面SHPB装置相同。SHPB装置中，子弹、入射杆、混凝土试件、透射杆长度分别为：400、3 200、32、1 600mm。Φ74mm直锥变截面SHPB装置简图如图1所示。

图1 Φ74mm直锥变截面SHPB装置图

1.2 主动围压装置

目前较为常用的SHPB围压实验装置，是通过在密封的油室内泵入液压油，待油室内压力达到所需压力值时，停止供油，压力通过弹性很好的橡胶隔层作用在试件周围，围压值大小可以在实验前自由调节，该类围压技术属于二维加载。实验过程中试件纵向受到冲击荷载作用，径向受到恒定围压的作用，试件在SHPB主动围压实验中的受力示意图如图2所示。虽然是通过二维动静组合加载的方式来近似考察材料三向应力状态下的动态响应，但是这比只考虑动荷载作用或只考虑静荷载作用更能反映工程的实际情况［7］。三向围压的主动围压装置，由于技术上较难实现，目前报道较少。

图2 主动围压下试件受力示意图

2 SHPB实验模拟与数据处理方法

2.1 SHPB有限元模型的建立

模拟对象为Φ74mm直锥变截面分离式霍普金森压杆，为减少计算量，采用1／4压杆模型建模，并在对称面上施加对称边界约束。对几何模型采用映射网格划分方法，入射杆、透射杆沿径向划分成16份，轴向网格长度为15mm。混凝土试样直径为74mm，厚度为32mm，网格划分时对试件做加密处理，沿轴向划分成8份，径向划分成32份，单元总数为58 176。入射杆、透射杆应力波采集位置分别位于距试件打击端1.2m和0.4m处。网格划分图，如图3所示。

图3 SHPB杆件及试件网格划分图

压杆、试件均采用Solid164三维实体单元，入射杆和透射杆采用相同的线弹性材料模型：密度ρ0＝7 850kg／m3，弹性模量E＝210GPa，泊松比μ＝0.25。混凝土试件采用HJC动态本构模型，该模型是文献［8］给出的混凝土材料在高应变、高应变率、高压下的动态本构模型，模型中考虑了材料损伤、应变率效应以及静水压力对于屈服应力的影响，能较好地描述混凝土在高速撞击与侵彻下所产生的损伤、破碎及断裂（或层裂）等行为。文献［9］对HJC模型进行修正，并且编制了计算软件，研究了混凝土靶体的侵彻过程，得到了满意的结果。

混凝土HJC动态本构模型在生成的关键字文件中修改、定义，关键字定义方式为“＊MATJOHNSON－HOLMGUIST－CONCRETE”，采用文献［8］中所给出的模型参数，主要参数见表1所列。

表1 HJC动态本构模型参数

考虑到变截面压杆波形弥散和应力波波头过冲严重，且采用波形整形器后需不断调整整形器材料参数及整形器尺寸才可以获得理想波形，过程较为繁琐，故采用弥散较小的三角形速度波加载。三角形速度波的加载原则：加载的速度波形长度应为子弹长度的2倍，直接加载于压杆端面节点上的速度波形的波峰值等于撞击速度的1／2［10］。模型中的接触面为：入射杆和混凝土试件接触、混凝土试件和透射杆接触。考虑到混凝土试件在加载过程中的单元失效，试件和压杆之间的接触面采用 ERODING－SURFACE－TOSURFACE，即侵蚀接触，侵蚀接触用于1个或2个表面的单元在接触时发生材料失效，接触依旧在剩余的单元中进行，它通常用于实体单元表面发生失效贯穿问题。通过引入“MAT－ADD－EROSION”侵蚀失效准则来控制单元的失效，且拉应变和切应变失效为单元的主要失效方式［11］。

2.2 SHPB模拟数据处理方法

对于大直径霍普金森压杆，反射波传播的距离较远，波形弥散较为严重，因此采用两波法即利用波形良好的入射波εi（t）、透射波εt（t）进行数据处理［12］，处理公式为：

其中，E、c0、A分别为压杆的弹性模量、弹性波速、横截面面积；As、L0分别为混凝土试样的初始横截面面积和长度。由于数值模拟可以避免实验中电噪声的干扰，因此入射波的起点位置可以取入射波与基线相交的位置［13］。透射波相对于入射波的延迟时间τ由（3）式确定，即

分别对子弹速度为10、20m／s时的混凝土SHPB实验进行了模拟，并利用（1）～（3）式所描述的两波法重构了混凝土的应力应变曲线，如图4所示。

图4 不同子弹冲击速度下的混凝土应力－应变曲线

由图4可知，在冲击速度较低时，即应变率较低条件下，混凝土峰值应变小于破坏应变，由于损伤演化不充分，应力－应变曲线的卸载段表现为黏弹性机制的回滞。当冲击速度较高时，峰值应变大于破坏应变，应力－应变曲线具有明显的弹性上升段、损伤演化段和结构破坏段。由模拟结果体现出的应变率强化效应、未出现宏观裂纹时的黏弹性回滞等，均与混凝土SHPB实验体现出的力学行为本质非常相近［14］。因此，采用上述模拟方法和数据处理手段可以对混凝土SHPB实验进行有效的模拟。

3 主动围压条件的数值模拟

在常规SHPB数值模拟的基础上，对SHPB主动围压实验做了模拟，即在试件周围加上固定围压值。该项模拟存在的难点主要体现在：圆柱形试件径向压力施加困难；实验过程中压力值难以保持稳定，因此对该模型做了以下处理：

（1）对试件外围单元按轴向划分为均等的2份，定义为2个不同的单元组，调试出适合各单元组的Surface Key。

（2）围压值曲线由零开始定义，在0.4ms时达到所需围压值，而后压力值保持恒定，压力作用时间超出压杆冲击作用的时间。

（3）调整合适的网格大小及刚度罚因子k，在原有的网格划分条件下设置刚度罚因子值为0.4。

利用LS－PREPOST后处理器，在32MPa围压条件下，随机选取试件表面单元，绘制出单元压力－时间历程曲线，如图5所示。

图5 试件外围单元压力－时间曲线

由图5可以看出，在围压条件施加初期，压力产生微弱振荡，随着围压的加载，压力值逐渐稳定。此处的压力为各向主应力的平均值，计算公式为：

在该模型中，试件单元的压力值p＝64／3MPa。随机选取多个试件表面及内部单元，将所选各单元压力值平均后，得到均化的p－t曲线，如图6所示。

图6 试件单元平均后的压力－时间曲线

由图6可以看出，混凝土试样表面及内部单元应力均匀性较好，且模拟得到的围压值与所需定义的围压吻合较好，即采用LS－DYNA有限元程序可以有效模拟SHPB主动围压实验。

在上述基础上，对混凝土试样进行了围压大小为：4、8、12、16、20、24MPa条件下的数值模拟，子弹冲击速度为24m／s，重构出各围压等级下的应力－应变曲线，如图7所示。

图7 不同围压下混凝土应力－应变曲线

由图7可以看出，在子弹冲击速度一定的条件下，混凝土抗压强度随主动围压值的增加而增加；峰值应变随围压的增大而增加，即材料延性增加；当材料充分损伤演化时，混凝土破坏时的应变值同样随围压的增大而增加，即材料韧性提高。由于围压的作用限制了材料的横向变形，试样要达到压溃体积应变，则需提高轴向应力，即表现出材料的峰值应力增加。当主动围压值较大时，受应变率的限制，材料没有发生屈服和塑性变形即卸载，此时的体积模量与线弹性段相差不大，应力应变曲线上表现为与线弹性段几乎平行的卸载。例如当围压值达到24MPa时，试件没有产生宏观裂纹，导致材料损伤演化不充分，应力应变曲线卸载段表现出一定程度的黏弹性回滞，即围压产生的强化效应使得试件抗冲击能力提高，需进一步提高应变率（此处表现在冲击速度的增加），才能使得试件在较高围压下损伤演化充分。由于损伤的不充分，材料单元发生不均匀失效，透射波的下降段出现震荡，表现在应力应变曲线上为卸载前产生一段较为明显的强化过程，该过程为材料卸载时应力不均匀产生的波形震荡。

围压值与对应的试件峰值应力呈现出较好的线性关系，如图8所示，但是围压值的增大将导致应变率的减弱，考察相同应变率条件下不同围压对试件力学性能的影响，需提高子弹的冲击速度，经LS－DYNA计算可知，24MPa围压比4MPa围压的子弹冲击速度高1.5m／s左右。相同应变率下围压效应对混凝土材料力学性能影响情况的数值模拟有待于进一步研究。

图8 各围压值与峰值应力的函数拟合图

4 结束语

（1）模拟结果显示，LS－DYNA 作为SHPB实验的辅助分析软件，同样可以较好地模拟SHPB主动围压实验，与SHPB形成有效的互补，下一步将会对于实验中难以达到的三维主动围压进行数值模拟。

（2）为了确保主动围压值在模拟过程中保持稳定，本文推荐的解决方法如下：调整网格大小，保证网格划分匹配；调整罚因子k，如果接触刚度过大将会导致计算不稳定；添加刚度阻尼，合理的刚度阻尼值将会有效减少压力振荡；调整接触类型，面面接触不仅可以保证计算稳定且可以减少计算机时。

（3）混凝土试样损伤演化不充分时，应力应变曲线在卸载段中后期会表现出一定程度的黏弹性回滞，此时的应变率表现为负值，影响平均应变率计算的有效性。且卸载段中后期，数值模拟和实验有较大不同，数值模拟采用了拉式算法及失效准则删除单元，单元删除后形成的空洞不会被填补，试样声阻抗大大降低，表现在应力－应变曲线上为曲线卸载段较陡，进一步影响了平均应变率计算的有效性。因此，平均应变率的计算有待进一步规范。

［1］何满潮，谢和平，彭苏萍，等.深部开采岩体力学研究［J］.岩石力学与工程学报，2005，24（16）：2803－2813.

［2］郭文兵，李小双.深部煤岩体高温高围压下力学性质的研究现状与展望［J］.河南理工大学学报：自然科学版，2007，26（1）：16－20.

［3］薛志刚，胡时胜.水泥砂浆在主动围压下的动态力学性能［J］.爆炸与冲击，2008，28（6）：561－564.

［4］刘军忠，许金余，吕晓聪，等.大直径分离式Hopkinson压杆主动围压加载实验技术［C］／／冯夏庭，李海波.第十一次全国岩石力学与工程学术大会论文集，2010：161－166.

［5］施绍裘，喻 炳，王礼立.PP／PA共混高聚物在多轴压缩应力下的冲击响应［C］／／谭 华.中国力学学会第八届全国冲击动力学学术讨论会论文集，2008：54－61.

［6］刘孝敏，胡时胜.应力脉冲在变截面SHPB锥杆中的传播特性［J］.爆炸与冲击，2000，20（2）：110－114.

［7］左宇军，李夕兵，唐春安，等.二维动静组合加载下岩石破坏的实验研究［J］.岩石力学与工程学报，2006，25（9）：1809－1820.

［8］Holmquist T J，Johnson G R，Cook W H.A computational constitutive model for concrete subjective to large strains，high strain rates，and high pressures［C］／／Jackson N，Dickert S.The 14th International Symposium on Ballistics.A－merican Defense Prepareness Association，1993：591－600.

［9］张凤国，李恩征.大应变、高应变率及高压强条件下混凝土的计算模型［J］.爆炸与冲击，2002，22（3）：198－202.

［10］董 钢.分离式霍普金森压杆实验技术数值模拟［D］.合肥：合肥工业大学，2005.

［11］李 耀，李和平，巫绪涛.混凝土HJC动态本构模型的研究［J］.合肥工业大学学报：自然科学版，2009，32（8）：1244－1248.

［12］巫绪涛，杨伯源，李和平，等.大直径SHPB装置的数值模拟及实验误差分析［J］.应用力学学报，2006，23（3）：431－434.

［13］李英雷，胡昌明，王 悟.SHPB实验数据处理的规范化问题讨论［J］.爆炸与冲击，2005，25（6）：553－557.

［14］巫绪涛，孙善飞，李和平.用 HIC本构模型模拟混凝土SHPB实验［J］.爆炸与冲击，2009，29（2）：137－142.