2. 中国船舶集团有限公司第七〇九研究所,湖北 武汉 430205;
3. 武汉理工大学,湖北 武汉 430063
2. The 709 Research Institute of CSSC, Wuhan 430205, China;
3. Wuhan University of Technology, Wuhan 430063, China
舰船结构的抗爆性能是舰船生命力的重要保障,目前国内外学者针对结构抗爆已开展了较多研究。战斗部产生的爆炸载荷主要包括冲击波和破片。其中破片载荷多以高速破片群的形式出现,较多文献研究了破片的质量、速度在空间上的分布。Mott[1]采用统计学方法,从微观角度分析金属壳体破碎机理,推导了破片平均质量、数量的累积分布表达式。在Mott公式的基础上,Cohen[2]进一步给出了Weibull分布式,为前者提供了一个合理的统计理论。爆炸冲击波的作用过程为准静态过程,其传播符合质量、能量和动量守恒定律。Taylor[3],Bach[4],Henrych[5],Brode[6]和张玉磊[7]等采用理论或试验的方法,先后提出了估算爆炸冲击波的多个经验公式,为爆炸冲击波的预估提供了参考。李营等[8]设计了多舱室缩比模型,开展了舱内充水和无水的对比试验研究,结果表明水对准静态压力有衰减作用,进而会对舱室毁伤产生影响。
舱内爆炸下舱室结构响应的研究多以加筋板结构为基础,大量文献通过理论推导[9]、模型试验[10]、数值仿真[11]等方法研究了加筋板结构的失效模式,探讨不同的变形模式,建立简化模型为舱内结构变形进行合理预报。刘土光等[12]、黄震球[13]通过假定板架变形的模态,从理论上分析了加筋强弱对板架变形的影响。梅志远等[14]使用仿真手段研究了爆炸载荷下船用加筋板架毁伤效应的非线性问题,同时结合试验进行验证,结果表明大尺寸加筋板架分散了冲击波的破坏作用,有利于结构安全性的保障。
综合来看,已有文献主要集中于舱内爆炸下船体板架结构响应的研究,给出了冲击波和破片载荷的经验公式与预估方法。但在内爆后船体结构毁伤效应方面,多数基本只考虑了冲击波单独作用下的情况,考虑冲击波和破片联合作用还相对欠缺。本文以船舶舱段为研究对象,运用数值方法,考虑在舱室不同位置爆炸时,舱内爆炸冲击波与破片联合的毁伤效应。
1 舱段有限元模型本文的研究对象为某船的三舱段模型,以中间舱段作为目标舱段来进行舱室内爆数值仿真计算。相邻2个舱段被横舱壁隔开,每个舱段设置了3层甲板和双层底。舱段典型横剖面示意图如图1所示。
利用ANSA建立舱段的有限元模型,如图2所示。考虑到AUTODYN不能识别BEAM单元,在ANSA建模时仅使用SHELL单元。为提高计算效率和保障单元质量,对舱段模型局部细节进行了简化处理。舱段有限元模型的基本网格尺寸为135 mm,最大网格尺寸为150 mm,最小网格尺寸为40 mm。
舱段结构采用船用945钢,失效应变取为0.25[15],在AUTODYN中选用双线性弹塑性本构模型,船用945钢的材料参数如表1所示。
材料的应变率效应使用Cowper-Symonds(C-S)模型进行计算。C-S模型可表示为:
$ \mathop {{\varepsilon _p}}\limits^ \bullet = D{(\frac{{{\sigma _d}}}{\sigma } - 1)^q} 。$ | (1) |
其中:
战斗部壳体的材料为Q235钢,失效应变取为0.25,选用Johnson-Cook本构模型来描述材料行为。孔祥韶[16]通过准静态拉伸试验确定了Q235钢的初始屈服应力即参数
战斗部内装有TNT炸药,状态方程选用JWL状态方程。其爆炸速度为6 930 m/s,爆压为2.1×107 kPa,密度ρ=1630 kg/m3,比内能e=4.17×106 J/kg,表达式为:
$ P = A \left(1 - \frac{{\omega \eta }}{{{R_1}}}\right){e^{\frac{{ - {R_1}}}{\eta }}} + B \left(1 - \frac{{\omega \eta }}{{{R_2}}}\right){e^{\frac{{ - {R_2}}}{\eta }}} + \omega \eta \rho e 。$ | (2) |
式中:
空气采用Gamma方程,其表达式为:
$ p = (\gamma - 1)\rho e。$ | (3) |
式中:
本文数值仿真时选取的炸药当量为150 kg。战斗部形状选用圆柱形,战斗部高为980 mm,壳体厚为20 mm,TNT的质量为150 kg。战斗部的有关参数如表3和图3所示。根据战斗部的尺寸,舷侧破口设置为圆形,破口直径选取为400 mm,假定战斗部在01甲板和1甲板之间发生爆炸,战斗部穿甲的位置位于01甲板和1甲板之间的纵舱壁,如图4所示。
TNT药柱和战斗部金属壳体都通过SPH粒子进行定义,由式(4)和式(5)计算可得SPH粒子大小为20 mm,TNT的SPH粒子个数为11520,战斗部外壳的SPH粒子个数为3964,共计15484个粒子;圆柱形战斗部的起爆位置位于圆柱体的中心处。战斗部计算模型如图5所示。
$ {\mu ^{0.5}} = \frac{{K{t_0}{{({t_0} + {d_i})}^{1.5}}{{\left(1 + \dfrac{M}{{2{m_t}}}\right)}^{0.5}}}}{{{d_i}}} ,$ | (4) |
$ {m_p} = 2\mu。$ | (5) |
式中:
为研究冲击波和破片联合作用下的舱段内爆载荷特性和舱段结构响应,需要计算战斗部爆炸时的TNT等效裸装药量。利用Fisher公式[17]计算等效裸装药量:
$ \frac{{{C_{EB}}}}{C} = 0.6 + \frac{{0.4}}{{1 + M/C}} 。$ | (6) |
式中:
将战斗部的有关参数代入式(6)可知,计算冲击波和破片联合作用下战斗部的等效裸装TNT药量为115 kg。
2.4 爆炸位置圆柱形战斗部爆炸时起爆点位于中间舱段的中横剖面内且战斗部的轴向与Y轴平行。本文只考虑中心爆点工况,即将起爆点设置在01甲板与1甲板之间舱室中纵剖面1/2舱室高度处,如图6所示。
欧拉域建成长方体,欧拉域沿X、Y、Z方向的长度分别为13 300 mm、24 000 mm和10 000 mm,完全覆盖中间舱段的01甲板、1甲板、2甲板和横舱壁。采用全耦合方式建立舱室的耦合面,分别建立中间舱段01甲板与1甲板之间,第一层舱室内的空气区域、1甲板与2甲板之间第二层舱室内的空气区域和舱外空气区域3个空气区域。欧拉域中的欧拉网格数目近似为817万个,不同的欧拉区域互相采用耦合面即拉格朗日单元进行隔离。舱段结构与欧拉域的耦合模型如图7所示。
为了研究舱室内部冲击波压力的变化情况,在欧拉域中布置了5个测点,用于测量舱室内爆后的压力变化,1 ~ 5号测点与起爆点处于同一横剖面内且位于同一水平线上,其与起爆点的距离分别为1 600 mm、3 600 mm、7 600 mm、8 660 mm和9 710 mm,1 ~ 5号测点的具体位置如图8所示。
为研究舱室结构的响应,在舱室结构上布置了8个测点,这些测点与起爆点在同一横剖面内。其中6 ~ 9号测点位于01甲板上,离中纵剖面的距离分别为0 mm、1000 mm、2000 mm和3000 mm;10 ~ 13号测点位于1甲板上,离中纵剖面的距离分别为0 mm、1000 mm、2000 mm和3000 mm。6 ~ 13号测点的具体位置如图9所示。
研究裸装药爆炸冲击波作用下,舱内冲击波传播和压力特性、舱室结构等效塑性应变及位移。TNT炸药质量为150 kg,形状为圆柱体,直径为360 mm,高度为900 mm。计算时间为75 ms。
3.2.1 冲击波传播和压力特性由AUTODYN数值仿真得到的不同时刻下,舱内冲击波压力云图如图10所示。
可知,在0.3 ms时刻爆炸产生的冲击波首先接触01甲板与1甲板,接着发生反射;在1.5 ms时刻入射冲击波与甲板的相互作用基本已完成,产生的反射冲击波互相叠加并向舷侧方向传递;在11 ms时刻冲击波到达舷侧,并在角隅位置发生汇聚,角隅处出现了很大的峰值压力;之后,舱室内部冲击波不断发生反射,冲击波能量逐渐减弱,在75 ms时刻舱室内部的压力已经较为稳定,维持一种准静态状态。
由AUTODYN数值仿真得到的测点1 ~ 测点5的冲击波压力时间历程曲线如图11所示。
可知,测点1的入射冲击波压力峰值大于测点2的入射冲击波压力峰值,其原因是压力峰值随着离爆点距离的增大而减小。从入射冲击波压力峰值的出现时间来看,距离爆点越远,入射冲击波峰值超压的出现时间越晚。在测点1处,在入射冲击波压力峰值出现后,又出现了2个非常明显的压力峰值,且后一个压力峰值大于前一个压力峰值,这种多个压力峰值的现象也体现了不同反射冲击波在传播过程中,发生了追赶和叠加。在测点2处,在入射冲击波压力峰值出现后,很快又出现了一个更高的压力峰值。测点3和测点4均离爆点较远,测点4的入射冲击波压力峰值大于测点3的入射冲击波压力峰值,且在入射冲击波压力峰值出现后,又连续出现了3个越来越高的压力峰值。这是由于反射冲击波的追赶和叠加,再加上测点靠近舷侧,在入射冲击波到达后,测点3和测点4不仅接受到来自01甲板与1甲板的反射冲击波,还受到来自舷侧的反射冲击波影响。但由于测点3离舷侧较测点4更远,所以其压力峰值小于测点4。测点5的冲击波压力峰值大于测点3的冲击波压力峰值,其原因是测点5位于破口位置,破口对于舱内爆炸具有泄爆作用。当数值仿真结束时,除舷侧预先设置的破口外,爆炸载荷没有对舱室结构造成明显的破口,由于预置破口位置离爆点较远且面积较小,破口的泻爆作用并不明显。爆炸当舱维持了一种准静态压力状态,准静态压力的大小约为0.35 MPa。
3.2.2 舱室结构等效塑性应变及位移响应由AUTODYN数值仿真得到的01甲板上测点6 ~测点 9以及1甲板上测点10 ~ 测点13的等效塑性应变时间历程曲线分别如图12和图13所示。
将各测点的等效塑性应变时间历程曲线进行对比,发现01甲板和1甲板上的4个测点中,距离爆点较近的前2个测点在爆炸发生后很快产生了塑性应变,而距离爆点最远的两个测点直到10 ms后才开始产生明显的塑性应变。因此,冲击波载荷对01甲板和1甲板的影响半径在2 000 ~ 4 000 mm之间,超过4000 mm后冲击波载荷就已经对甲板产生的等效塑性应变影响较小,该处的塑性应变主要由准静态压力造成。由此可得结论,甲板整体的大变形更多取决于准静态压力而非冲击波压力。
由AUTODYN数值仿真得到的75 ms时刻爆炸当舱各板架的等效塑性应变云图和位移云图,如图14所示。
从图14可以看出,冲击波单独作用下的爆炸当舱除了预先设置的破口以外,没有产生破口。内爆载荷加载可以分为2个阶段,即冲击波载荷加载和准静态载荷加载,这2种载荷先后作用于爆炸当舱的各板架上,造成各板架最终的变形。01甲板板架和1甲板板架呈现典型的板架大变形模式,有明显的颈缩现象。01甲板板架的中心处和1甲板板架与舱壁连接处的等效塑性应变较大,01甲板板架和1甲板板架的最大位移分别为0.772 m和0.769 m,1甲板的最大位移要小于01甲板是因为爆炸邻舱对1甲板的抗爆能力有加强作用。横舱壁和纵舱壁的最大位移相对于甲板来说较小,且横舱壁的位移峰值位于板架中心,而舷侧板架的位移峰值位于舷侧板架的中心处。甲板的位移较大,可能的原因是爆点位置距离甲板较近,甲板受到的冲击波载荷很大,且舱室宽度较大,能够变形的程度也就更大,产生的垂向位移也越大。
4 舱内爆炸冲击波和破片联合作用下的数值仿真 4.1 测点布置为研究舱室内部冲击波压力的变化,在欧拉域中布置了6个测点,用于测量舱室内爆后的压力变化,1 ~ 5号测点与起爆点处于同一横剖面内且位于同一水平线上,其与起爆点的距离分别为1 600 mm、3 600 mm、7 600 mm、8 660 mm和9 710 mm。6号测点与起爆点位于同一垂直线上且位于起爆点的下方,其与起爆点的距离为2 700 mm,布置6号测点的目的是考虑到1甲板可能会出现破口,冲击波可能会泄漏到爆炸邻舱。1 ~ 6号测点的具体位置如图15所示。
为了研究舱室结构的响应,在舱室结构上布置了11个测点,这些测点与起爆点在同一横剖面内。其中7 ~ 11号测点位于01甲板上,离中纵剖面的距离分别为0 mm、2000 mm、4000 mm、6000 mm和8000 mm;12 ~ 16号测点位于1甲板上,离中纵剖面的距离分别为0 mm、2000 mm、4000 mm、6000 mm和8000 mm;17号测点位于2甲板上,离中纵剖面的距离为0 mm。7 ~ 17号测点的具体位置如图16所示。
在爆点位置同时布置了战斗部和TNT裸装药,TNT等效裸装药药量为115 kg,形状为圆柱体,直径为320 mm,高度为880 mm。战斗部高度为980 mm,壳体厚度为20 mm,装填150 kg直径360 mm、高900 mm的TNT,计算时间为75 ms。
4.2.1 冲击波传播和压力特性、破片飞散及舱室的毁伤过程由AUTODYN数值仿真得到的不同时刻下舱内冲击波压力云图,如图17所示。
从图17可已看出,在0.38 ms时刻爆炸产生的冲击波接触01甲板与1甲板,接着发生反射;在2.35 ms时刻爆炸当舱的01甲板以及1甲板的中心位置在冲击波和破片的联合作用下出现变形和破口,冲击波通过破口向舱外和爆炸邻舱传播,并且01甲板以及1甲板产生的反射冲击波发生叠加并向舷侧方向传播;在7.7 ms时刻冲击波到达了舷侧,冲击波在角隅位置发生汇聚,角隅处出现了很大的峰值压力。之后,舱室内部冲击波通过01甲板以及1甲板中心位置的大破口继续泄漏,并在舱室内部不断发生反射,冲击波能量渐渐减弱,在75 ms时刻舱室内部的压力已较为稳定,维持一种准静态状态。
由AUTODYN数值仿真得到的不同时刻下,舱内冲击波和破片联合作用对舱室的毁伤过程如图18所示。
从图18可以看出,在0.11 ms时刻战斗部外壳已经在爆轰产物作用下破碎;在3.84 ms时刻01甲板、1甲板及横舱壁在冲击波和战斗部破片联合作用下遭到毁伤并产生破口,冲击波和破片开始向舱外和爆炸邻舱(1甲板与2甲板之间舱室)发生泄漏和飞散;在8.78 ms时刻爆炸当舱(01甲板与1甲板之间舱室)内破片先于冲击波抵达2甲板及舷侧,并对其结构造成毁伤;在75 ms时刻大量破片通过破口向舱外飞散或汇集在角隅处,舱内冲击波和破片对舱室结构的联合毁伤作用趋于稳定。
由AUTODYN数值仿真得到的测点1 ~ 测点6的冲击波压力时间历程曲线如图19所示。
由图19可知,冲击波和破片联合作用下,测点1 ~ 测点5的冲击波压力曲线前期趋势与冲击波单独作用下的大致相同。这是因为冲击波单独作用下甲板不会产生破口,冲击波只在爆炸当舱(01甲板与1甲板之间舱室)内传播,而冲击波和破片联合作用下,01甲板、1甲板以及横舱壁均出现大破口。冲击波通过破口向舱外和爆炸邻舱(1甲板与2甲板之间舱室)泄漏,所以爆炸当舱在30 ms后就慢慢趋于一个准静态压力状态,准静态压力的大小约为0.25 MPa。由于1甲板出现了大破口,冲击波向爆炸邻舱(1甲板与2甲板之间舱室)传播,由测点6的压力时程曲线可看到测点6的冲击波压力峰值没有超过0.5 MPa,并在冲击波的反射作用下,随后出现了一个较小的峰值。最后逐渐趋于一个准静态压力状态,准静态压力的大小约为0.2 MPa。
4.2.2 舱室结构等效塑性应变及位移响应由AUTODYN数值仿真得到的01甲板上测点7 ~ 测点11、1甲板上测点12 ~ 测点16以及2甲板上测点17的等效塑性应变时间历程曲线,分别如图20、图21及图22所示。
将01甲板、1甲板以及2甲板上各测点的等效塑性应变时间历程曲线进行对比,可以发现测点7和测点12的等效塑性应变曲线在达到0.25以后就消失了。这是由于在AUTODYN中定义的材料失效应变为0.25。此外,冲击波和破片联合作用下,01甲板和1甲板上各测点的等效塑性应变曲线变化趋势与冲击波单独作用下的大致相同,但各测点的等效塑性应变明显小于冲击波单独作用下相应测点的等效塑性应变。这是因为01甲板、1甲板以及横舱壁很早出现了大破口,冲击波很早就向外泄漏,冲击波载荷受到削弱,舱室内的准静态压力也较早达到平衡。01甲板的等效塑性应变的最大值明显大于1甲板,这是因为爆炸邻舱对1甲板的抗爆性起到了加强作用。位于2甲板上的测点17的等效塑性应变峰值很小,只有0.001左右。这是由于通过1甲板破口传播而来的冲击波载荷,以及1甲板与2甲板之间舱室内准静态压力较小。
由AUTODYN数值仿真得到的75 ms时刻爆炸当舱和爆炸邻舱各板架的等效塑性应变云图和位移云图,如图23所示。
从图23可以看出,冲击波和破片联合作用下的毁伤区域主要为爆炸当舱的01和1甲板以及横舱壁,其中01和1甲板的损伤主要表现为甲板出现大破口以及角隅处出现的较大变形,横舱壁的损伤则是其边缘处出现面积较大的断裂。而冲击波单独作用下爆炸当舱除了预先设置的破口以外,不会产生明显破口。爆炸当舱的01甲板、1甲板及横舱壁分别出现了1.6 m×1.7 m、1.5 m×1.6 m和1.3 m×2.1 m大小的破口,01甲板及1甲板的最大位移分别为1.4 m和1.3 m,比冲击波单独作用下01甲板及1甲板的最大位移分别大了82.9%和69.8%。爆炸当舱的纵舱壁虽因离起爆点较远而未出现破口,但其最大位移和等效塑性应变较冲击波单独作用下大得多。由于冲击波单独作用下1甲板并未破损,故爆炸邻舱的2甲板及横舱壁并未受到太大影响,而在冲击波和破片联合作用下,爆炸邻舱的各板架均受到了冲击波和破片的毁伤作用,2甲板及横舱壁虽然未出现破口,但均出现了一定塑性变形,2甲板的最大塑性变形为0.2 m。
5 结 语1)只考虑冲击波的单独作用时,爆炸当舱除了预先设置的破口以外,没有产生其他破口。各板架最终的变形主要是由冲击波载荷和准静态载荷造成的,舱内结构整体的塑性应变并不大。
2)考虑冲击波和破片的联合作用时,爆炸当舱的01甲板、1甲板和横舱壁均出现了大破口。冲击波压力曲线的前期趋势与冲击波单独作用下大致相同。但由于大破口的存在,冲击波发生泄漏,舱室内压力较早达到准静态压力状态。
3)冲击波和破片联合作用下的毁伤区域主要为爆炸当舱的01和1甲板以及横舱壁,其中01和1甲板的毁伤主要表现为大破口以及角隅大塑性变形,横舱壁的毁伤主要表现为其边缘大面积撕裂。爆炸当舱舷侧的最大位移和等效塑性应变,较冲击波单独作用下大得多,爆炸邻舱的2甲板和横舱壁在冲击波和破片的联合作用下也出现了一定塑性变形。
[1] |
MOTT N F. Fragmentation of shell cases[J]. Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences, 1947, 189(1018): 300–308.
|
[2] |
COHEN JR E A. New formulas for predicting the size distribution of warhead fragments[J]. Mathematical Modelling, 1981, 2(1):19–32.
|
[3] |
TAYLOR G. The formation of a blast wave by a very intense explosion. II. The atomic explosion of 1945[J]. Proceedings of the Royal Society of London Series A, 1950, 201(1065): 175-186. |
[4] |
BACH G G, LEE J H S. An analytical solution for blast waves[J]. AIAA Journal, 1970, 8(2): 271-275. DOI:10.2514/3.5655 |
[5] |
HENRYCH J. The dynamic of explosion and its use[M]. Amsterdam: Elsevier, 1979.
|
[6] |
BRODE H L. Numerical solutions of spherical blast waves[J]. Journal of Applied Physics, 1955, 26(6): 766-775. DOI:10.1063/1.1722085 |
[7] |
张玉磊, 苏健军, 李芝绒. TNT内爆炸准静态压力特性[J]. 爆炸与冲击, 2018, 38(6): 252-257. |
[8] |
李营, 任广为, 张玮, 等. 水介质对舱内爆炸抑制作用的实验研究[J]. 爆炸与冲击, 2017, 37(6): 1-7. |
[9] |
刘敬喜, 刘尧, 汤皓泉, 等. 爆炸载荷作用下单向加筋方板的大挠度塑性动力响应分析[J]. 振动与冲击, 2011, 30(4): 182-187. |
[10] |
CHUNG KIM YUEN S, NURICK G N. Experimental and numerical studies on the response of quadrangular stiffened plates. Part I: subjected to uniform blast load[J]. International Journal of Impact Engineering, 2005, 31(1): 55-83. DOI:10.1016/j.ijimpeng.2003.09.048 |
[11] |
KADID A. Stiffened plates subjected to uniform blast loading[J]. Journal of Civil Engineering and Management, 2008, 14(3): 155-161. DOI:10.3846/1392-3730.2008.14.11 |
[12] |
刘土光, 胡要武, 郑际嘉. 固支加筋方板在爆炸载荷作用下的刚塑性动力响应分析[J]. 爆炸与冲击, 1994, 14(1): 55-65. |
[13] |
黄震球. 固支加筋方板的大挠度塑性动力响应[J]. 固体力学学报, 1995, 16(1): 5-12. |
[14] |
梅志远, 朱锡, 刘润泉. 船用加筋板架爆炸载荷下动态响应数值分析[J]. 爆炸与冲击, 2004, 24(1): 80-84. |
[15] |
曹广博. 空中爆炸载荷作用下船体结构动态响应及剩余强度研究[D]. 镇江: 江苏科技大学, 2019.
|
[16] |
孔祥韶. 爆炸载荷及复合多层防护结构响应特性研究[D]. 武汉: 武汉理工大学, 2013: 137-154.
|
[17] |
DUNNETT J, FLYNN D, WHARTON J. Blast algorithm development: definition of modified blast algorithms for PBX based explosives[C]//Insensitive Munitions and Energetic Materials Technical Symposium IMEMTS, 2006.
|