【化学工程与材料科学】
结构抗爆分析通常采用以数值分析为主、试验为辅的研究方法。通用动力学软件LS-DYNA被广泛地应用于抗爆数值分析中,LS-DYNA为用户提供了4种不同的爆炸荷载施加方法。为研究何种爆炸荷载施加方法最能反映真实的情况,首先选用具有代表性的钢筋混凝土板(以下简称RC板)抗爆试验文献作为研究对象,然后利用LS-DYNA软件建立与试验文献完全相同的数值分析模型并采用4种不同的方法施加爆炸荷载,最后通过与文献的试验结果进行比较选择出最合理的爆炸荷载施加方法。
图1为自由空气场爆炸产生的冲击波典型压力时程曲线,曲线由正压与负压作用两部分组成,并且负压峰值远小于正压峰值。文献[1-2]将图1的典型压力时程曲线简化为三角形荷载曲线,荷载曲线仅考虑正压作用,曲线由正压峰值Pmax与正压作用时间两参数确定,如图2所示。
图1 爆炸波冲击压力时程曲线
图2 简化后的三角形压力时程曲线
爆炸力学采用爆炸点至结构的距离R与药量W1/3的比值(定义为比例距离Z)作为衡量爆炸超压作用的单位。相同距离R的条件下,药量越大,则比例距离越小、产生的超压荷载越大。采用方法1施加爆炸荷载时需要根据比例距离Z计算出超压峰值Pmax与正压作用时间亨利奇(Henrgeh)在实验基础上得到了超压峰值Pmax与正压作用时间的计算公式为[3]:
Pmax=1.407 17/Z+0.553 97/Z2-0.035 72/Z3+
0.000 625/Z4 (0.05≤Z≤0.30)
(1)
Pmax=0.619 38/Z-0.032 62/Z2+0.213 24/Z3
(0.30≤Z≤1)
(2)
Pmax=0.0662/Z+0.405/Z2+0.3288/Z3
(1≤Z≤10)
(3)
(0.107+0.444Z+0.264Z2-0.129Z3+
0.033 5Z4)×10-3×W1/3
(4)
根据式(1)~(4)计算出已知比例距离Z的爆炸超压峰值Pmax与正压作用时间后,借助于LS-DYNA的关键字*DEFINE_CURVE将三角形超压时程曲线以均布面荷载的形式施加在爆炸波的接触面上。方法1不需要单独建立空气与炸药的数值分析模型,操作简单并且计算时间很短。
LS-DYNA中存有由试验得到的TNT空中爆炸超压数据,可以利用关键字按以下方法直接调用。首先采用*LOAD_BLAST_ENHANCED关键字定义TNT当量、爆炸点及爆炸类型,考虑地面反射作用时,还应借助关键字*DEFINE_VETOR定义地面与爆炸点之间的法向向量;然后采用关键字*LOAD_SEGMENT与*LOAD_BLAST_SEGMENT_SET定义爆炸波的接触面。文献[4-6]采用该方法施加爆炸荷载。方法2同样不需要建立空气与炸药的数值分析模型,施加方法简单,计算时间比方法1略长。
方法3在抗爆数值分析中应用最广。采用方法3的建模方法时,需要对炸药和空气分别单独划分网格并赋予其各自的状态方程,且炸药网格和空气网格应在交界面上共节点。文献[7-10]采用该方法施加爆炸荷载。炸药需要借助关键字*MAT_HIGH_EXPLOSIVE_BURN模型配合Jones-Wilkens-Lee(JWL)状态方程定义[10],以炸药爆轰产物的压力表示的JWL方程为:
(5)
式(5)中: P1为爆炸压力; V为相对体积; E0为初始内能密度;A、B、R1、R2、ω为状态方程参数。JWL方程参数与TNT材料参数取值如表1所示,表中ρ表示炸药密度,DT表示爆轰速度,Pcj表示爆轰压力。
空气采用*MAT_NULL配合*EOS_LINEAR_POLYNOMIAL状态方程定义[10],其线性多项式状态方程为:
P2=C0+C1μ+C2μ2+C3μ3+
(C4+C5μ+C6μ2)E1
(6)
式(6)中: μ=ρ/ρ0-1, ρ表示空气的当前密度, ρ0是初始时刻的空气密度; P2为空气压力;C0~C6为多项式方程系数; E1为内能密度; V0为初始相对体积。各相关参数取值见表2。
表1 TNT炸药参数[10]
ρ/(kg·m-3)A/GPaB/GPaR1R2ωE0/ (J·m-3)VDT/(m·s-1)Pcj/GPa1 631371.23.234.150.950.37×1091.06 93021.4
表2 空气参数[10]
ρ0/(kg·m-3)C0/MPaC1C2C3C4C5C6E1/(J·m-3)V01.293-0.10000.40.402.53×1051.0
方法3操作比以上方法1和方法2复杂,并且需要消耗较长的计算时间。
方法4是在方法3基础上改进而得到的一种爆炸荷载施加方法,其施加方式详见文献[11-13]。该方法首先需要定义炸药、空气的材料参数与状态方程;然后借助关键字*ALE_MUTI-MATERIAL_GROUP定义一个包含炸药与空气的材料组;最后借助关键字*INITIAL_VOLUME_FRACTION_GEOMTRY进行填充。以图3为例说明方法3与方法4的区别,图3所建模型为1 m×1 m×1 m的正方体空气与半径0.2 m的球形炸药,当采用方法3建模时需分别单独建立空气与炸药的网格。而采用方法4建模时仅需要建立空气网格作为背景网格,然后利用关键字*INITIAL_VOLUME_FRACTION_GEOMTRY在指定区域将炸药材料填充至建好的空气网格中,并且不要求空气网格和炸药网格共节点。操作完成后,在炸药几何模型所占据的空间中,原网格中的空气被炸药材料所替代。关键字*INITIAL_VOLUME_FRACTION_GEOMTRY的定义如表3所示,共3行。第1行的4个关键字定义初始的空气材料,其中FMSID表示空气的PART号;FMIDITYP表示PART类型;BAMMG表示空气材料在材料组中的编号;NTRACE表示填充采样点数量,可采用默认值。第2行的3个关键字定义炸药材料,其中CNTTYP表示炸药的形状,6为球形;FILLOPT表示填充方式,0表示内部填充;FAMMG表示炸药材料在材料组中的编号。第3行定义炸药的位置及体积,其中X0、Y0、Z0表示炸药球心坐标,R0表示球半径。
图3 方法3与方法4示意图
相对于方法3,方法4的优越性体现在以下2个方面:1) 不需要单独建立炸药单元,仅需在K文件中将炸药进行替换,可以非常方便地定义不同药量的爆炸荷载。2) 可以定义各种形状如圆形、方形、柱形等炸药数值分析模型,避免了不规则炸药体网格划分的困难。
表3 关键字*INITIAL_VOLUME_FRACTION_
GEOMTRY的应用
第1行FMSIDFMIDTYPBAMMGNTRACE数值111—第2行CNTTYPFILLOPTFAMMG—数值602第3行X0Y0Z0R0数值0000.2
孙文彬[14]对RC板进行了抗爆试验研究,笔者选取该文献作为参考对象。试验采用的2块RC板尺寸均为1 300 mm×1 000 mm×120 mm。板设计为单向板,长边简支,短边自由。采用双层双向配筋,受力的长边配筋为10@100,非受力的短边配筋为10@200。板保护层厚度取值25mm。试验钢筋采用HRB500,其屈服强度为560 MPa,极限强度为605 MPa,弹性模量为230 GPa。采用C40混凝土,抗压强度为48 MPa,弹性模量为48 GPa。试验在爱丁堡国防科技组织的爆炸仓内进行,爆炸点位于RC板中心正上部0.6 m处,爆炸Ⅰ、爆炸Ⅱ的TNT当量为0.079 kg、2.09 kg。爆炸Ⅰ比例距离Z=1.4 m/kg1/3,为了观测RC板在弹性阶段的动力响应。爆炸Ⅱ比例距离Z=0.47 m/kg1/3,为了观测RC板在塑性阶段的爆炸响应和破坏特征。试验布置及采用方法3建立的数值分析模型如图4所示。
图4 试验布置及数值分析模型示意图
结合试验提供的材料数据,采用LS-DYNA的72#材料模型模拟混凝土[15]。当采用国际单位时,72#模型定义C40强度混凝土材料的主要参数如表4所示。其中,MID为材料标识号; ρ为材料密度;A0为初始抗压强度,负值;RSIZE为单位转换数值:m/inches;UCF为单位转换数值:Pa/PSI;LCRATE为所要调用的曲线号。
表4 C40混凝土参数(72#模型)
ρ/(kg·m-3)A0/MPaRSIZE/无量纲UCF/无量纲LCRATE/无量纲2 320-4839.37145×10-6723
采用LS-DYNA的3#模型(*MAT_PLASTIC_KINEMATIC)模拟HRB500钢筋[15],钢筋的主要参数采用表5定义。其中, ρ表示密度; fy为屈服应力;HP为硬化参数;ν为泊松比;E和Etan分别为弹性模型和切线模量;C、P为Cowper-Symonds应变率参数;FS为失效应变。
添加关键字*MAT_ADD_EROSION定义混凝土单元的失效,定义最大抗压强度fc max≤40 MPa与最大抗拉强度 ft max≤4 MPa为失效准则,当材料强度大于以上任意数值时,混凝土单元将被删除。钢筋采用3#模型自带的失效准则定义,当钢筋的失效应变大于0.2时,钢筋单元将被删除。炸药与空气定义为多物质欧拉几何实体(EULER),结构体定义为拉格朗日几何实体(LAGRANGE),借助于关键字*CONSTRAINED_LAGRANGE_IN_SOILD实现爆炸波与结构体的接触,接触类型可采用无侵蚀的罚函数耦合。
表5 HRB500钢筋参数
ρ/(kg·m-3)E/GPaν/无量纲fy/MPaHP/无量纲C/无量纲P/无量纲Etan/GPaFS/无量纲7 8502300.356004052.30.2
爆炸Ⅰ的TNT当量为0.079 kg,爆炸点位于RC板中心正上部0.6 m处,比例距离Z=1.4 m/kg1/3。爆炸Ⅰ(Z=1.4 m/kg1/3)的RC板处于弹性阶段,相关文献给出了爆炸Ⅰ作用下RC板的跨中最大加速度amax、跨中最大挠度fmax及出现最大挠度的时间tmax。采用方法1施加爆炸荷载I时需要根据爆炸Ⅰ的比例距离Z=1.4 m/kg1/3并借助于式(3)、式(4)分别计算出Pmax=343 ms,然后将荷载曲线以面荷载的形式施加于爆炸波接触面上。表6为爆炸Ⅰ各种荷载施加方法的数据汇总,图5为4种爆炸荷载施加方法的RC板跨中位移时程曲线。
表6 爆炸I数据汇总
爆炸工况数据来源amax/(m·s-2)fmax/mmtmax/ms爆炸Ⅰ试验3101.1222爆炸Ⅰ方法13551.5330爆炸Ⅰ方法22620.3822爆炸Ⅰ方法34142.0533爆炸Ⅰ方法44402.6532
图5 爆炸I,RC板跨中位移时程曲线
通过对表6与图5分析可知,爆炸Ⅰ时:1) 方法1得到的数值分析结果与试验吻合最好,其最大加速度amax、最大位移fmax、最大位移出现时间tmax略大于试验结果。2) 方法2得到的最大加速度amax、最大位移fmax较试验结果偏小,而最大位移出现时间tmax与试验结合相吻合。3) 方法3、方法4得到的各数据均大于试验结果。由于试验结果的RC板在爆炸Ⅰ时最大位移fmax仅为1.12 mm,方法3、方法4比试验结果分别增大0.93 mm、1.53 mm。考虑到试验产生的误差,其计算结果仍然是合理的。(4)除方法4产生少量塑性变形以外,其余3种方法的RC板均处于弹性阶段,与试验结果相一致。以上说明,爆炸荷载较小时4种爆炸荷载施加方法均能较为准确地反映RC板在弹性阶段的动力响应。
爆炸Ⅱ的TNT当量为2.09 kg,爆炸点位于RC板中心正上部0.6 m处,比例距离Z=0.47 m/kg1/3。爆炸Ⅱ的RC板发生明显的破坏,图6为RC板顶与底的破坏图。由相关文献可知,板跨中表现为明显的弯曲破坏,支座未出现剪切破坏特征。板顶跨中出现120~150 mm的混凝土压碎区,板底出现2条宽度6~8 mm的主裂缝。板底形成近似正方形的坑,坑边长约为550 mm,周长约为2 200 mm,最大深度约为65 mm。板跨中产生不可恢复的永久变形,最大挠度为50 mm。
图6 爆炸Ⅱ,试验的RC板破坏图
采用方法1施加爆炸荷载Ⅱ时需要根据爆炸Ⅱ的比例距离Z=0.47 m/kg1/3并借助于式(2)、(4)分别计算出然后将荷载曲线以面荷载的形式施加于爆炸波接触面上。图7为爆炸Ⅱ,4种荷载爆炸施加方法的RC板破坏图,表7为爆炸形成的板底坑数据对比表,图8为4种爆炸荷载施加方法的RC板跨中位移时程曲线。通过对上述图表分析并结合试验数据可知:1) 方法1的RC板顶跨中出现宽度约200 mm的破坏区域,与试验结果吻合较好。板底产生大面积的塑性变形但未呈现明显的破坏特征,与试验结果不一致。板跨中最大挠度为118 mm,远大于试验结果,板破坏形式为跨中弯曲破坏与支座剪切破坏。2) 方法2的RC板顶破坏区宽度大于试验结果,板支座处产生通长裂缝,表现出明显的剪切破坏特征。与试验结果不同,板底坑的形状为圆形,坑最大深度及周长均小于试验结果。板跨中最大挠度略小于试验结果,板破坏形式为跨中弯曲破坏与支座剪切破坏。3) 方法3的RC板顶破坏区域宽度为150 mm,与试验结果一致。板底的纵横方向产生2条宽度为6~8 mm的主裂缝,与试验结果吻合。板底破坏形成接近正方形的坑,坑的周长、深度、最大宽度数值均非常接近试验结果。板跨中最大挠度为52 mm,与试验结果的误差仅为2%。板以跨中弯曲破坏为主,支座未发生明显的剪切破坏。4) 方法4的RC板顶破坏区域宽度为350 mm,大于试验结果。板支座处出现较为严重的塑性破坏,表现出剪切破坏特征。板底坑的形状为矩形,坑最大宽度、深度与周长均大于试验结果。板跨中最大挠度为58 mm,略大于试验结果,板破坏形式为跨中弯曲破坏与支座剪切破坏。
图7 爆炸Ⅱ,4种爆炸荷载施加方法的RC板破坏图
表7 爆炸Ⅱ深坑数据
爆炸工况数据来源形状最大长度/mm最大深度/mm周长/mm爆炸Ⅱ试验正方形550652 200爆炸Ⅱ方法1----爆炸Ⅱ方法2圆形-451 800爆炸Ⅱ方法3正方形540522 100爆炸Ⅱ方法4矩形950753 000
图8 爆炸Ⅱ,RC板跨中位移时程曲线
1) 当爆炸荷载较小时,RC板处于弹性阶段,4种方法均能较为准确地模拟结构构件的动力响应。
2)当爆炸荷载较大时,RC板处于塑性与破坏阶段。方法3可以较为精确地模拟RC板的破坏特征,方法2与方法4可以近似地模拟板的破坏特征,但与试验结果存在一定的误差,不应采用方法1施加爆炸荷载。
[1] 阎石,刘蕾,齐宝欣,等.爆炸荷载作用下方钢管混凝土柱的动力响应及破坏机理[J].防灾减灾工程学报,2011,31(05):477-482.
[2] 焦燏烽,赵果,侯延利.框架柱在爆炸冲击荷载作用下动力系数及破坏模式研究[J].建筑科学,2015,31(09):32-37.
[3] 亨利奇.爆炸动力学及其应用[M].北京:科学出版社,1987.
[4] 田力,付协伟.爆炸荷载作用下高层钢筋混凝土框架结构的连续倒塌机理研究[J].科学技术与工程,2016,16(03):232-238.
[5] 李忠献,师燕超,史祥生.爆炸荷载作用下钢筋混凝土板破坏评定方法[J].建筑结构学报,2009,30(06):60-66.
[6] 师燕超,李忠献,郝洪.爆炸荷载作用下钢筋混凝土框架的连续倒塌分析[J].解放军理工大学学报,2007,8 (06):652-658.
[7] SHA Y,HAO H.Nonlinear finite element analysis of barge collision with a single bridge pier[J].Engineering Strctures,2012(41):63-67.
[8] CUNG C,LEE J,JUNG R.Numerical simulations of missile impacts on reinforced concrete plates:IRIS-2010/2012 benchmark project[J].Nuclear Engineering & Design,2015,295:747-758.
[9] 索强,徐鹏,尤文斌.网格划分对冲击波波形的影响分析[J].兵器装备工程学报,2020,41(02):198-203.
[10] 陈铭,周云波,孙晓旺.TNT径高比对爆炸冲击波比冲量的影响研究[J].兵器装备工程学报,2018,39(06):62-66.
[11] 辛春亮,薛再清,涂建,等.True Grid和LS-DYNA动力学数值计算详解[M].北京:机械工业出版社,2019.
[12] 徐浩铭,顾文彬,刘建青,等.串联聚能炸药隔爆结构设计数值模拟和实验研究[J].兵工学报,2014,35(02):170-175.
[13] 毕程程,王志亮,石高扬,等.初始体积分数法在爆炸模拟中的应用[J].工程爆破,2017,23(04):26-33.
[14] 孙文彬.钢筋混凝土板的爆炸荷载试验研究[J].辽宁工程技术大学学报(自然科学版),2009,28(02):217-220.
[15] 李忠献,师燕超.建筑结构抗爆分析理论[M].北京:科学出版社.2015.