徐维铮,吴卫国
(1. 武汉理工大学 高性能舰船技术教育部重点实验室,湖北 武汉 430063;2. 武汉理工大学 交通学院船舶、海洋与结构工程系,湖北 武汉 430063)
当爆炸发生在舰船舱室内部时,由于冲击波传播过程中受到舱室壁面限制,将产生冲击的多次反射和汇聚,以及持续时间较长的准静态压力,对内部结构和人员产生更加严重的毁伤效应。近年来,国内外学者针对内爆炸特性和计算方法方面做了大量的研究工作。
WU等[1]开展了封闭舱室内爆炸实验,研究了炸药形状和起爆位置对爆炸冲击波载荷的影响规律。侯海量等[2]实验研究了舱室内爆炸冲击波载荷特性。孔祥韶等[3]实验研究了药量、角隅结构形式对舱内爆炸载荷影响规律。侯海量等[4]采用MSC.DYTRAN软件数值研究了舱室内爆炸冲击波载荷特性。孔祥韶等[5]采用MSC.DYTRAN软件数值模拟了舰船舱室内部战斗部爆炸及爆炸毁伤效应。樊壮卿等[6]基于LS-DYNA软件,数值模拟了大体积复杂舱室内爆炸载荷传播特性。陈攀等[7]采用AUTODYN软件数值研究了舱室内爆工况下不同区域爆炸冲击波载荷特性及爆点位置对舱室内爆炸冲击波特性的影响规律。
炸药爆炸属于高压比、高密度比问题,其对激波捕捉格式提出较为严格的要求。传统的商用软件(MSC.DYTRAN,LS-DYNA,AUTODYN)主要采用低阶精度的黎曼求解器和FCT欧拉算法,针对爆炸载荷的计算比较粗糙,且当前发展比较成熟的激波捕捉算法还没有嵌入到商用程序中,因此自主开发高精度爆炸载荷数值计算程序显得更有意义。
目前,激波捕捉格式众多,如TVD格式[8]、ENO格式[9]、NND格式[10]等,本文主要关注被广泛采用的WENO格式(weighted essentially non-oscillation scheme),LIU等[11]于1994年提出了该格式,SHU等[12–14]发展了该格式。WENO在格式的有效性、通量的光滑性和收敛解的稳定性方面均优于ENO格式。
本文基于FORTRAN平台,采用3阶、5阶、7阶WENO有限差分格式,开发了舱室内爆炸载荷高精度三维数值计算程序。采用所开发的程序开展了封闭舱室、泄压舱室内爆炸载荷数值计算,研究了WENO格式精度对舱室内爆炸载荷影响规律,为结构抗爆设计提供一定的参考和指导。
假定炸药瞬时爆轰,将凝聚相炸药等效为高压、高密度气体,流场的控制方程采用三维可压缩欧拉方程进行描述,
其中:
程序中采用3阶、5阶、7阶WENO有限差分格式对欧拉方程进行数值离散和求解,下面给出不同精度WENO格式的离散过程。
3阶WENO格式(WENO-JS3)的数值离散和推导过程如下:
式中:
上式即为本文程序中采用的3阶WENO格式的最终形式,对于含激波间断流场,式中的根据下式求得:
五阶WENO格式(WENO-JS5)的数值离散和推导过程如下:
式中:
式(16)即为本文程序中采用的5阶WENO格式的最终形式,对于含激波间断流场,式中的根据下式求得:
7阶WENO格式(WENO-JS7)的数值离散和推导过程如下,单元面中心处数值通量的4种重构方式分别为[16]:
式中:
式(25)即为本文程序中采用的7阶WENO格式的最终形式,对于含激波间断流场,式中的根据下式求得:
时间项采用3阶TVD-RK法求解,格式如下[17]:
为了初步考察上述格式(WENO-JS3,WENO-JS5,WENO-JS7)的计算性能,选取3个经典一维黎曼算例进行计算。
该算例初始条件如式(30)所示[18],网格数为400,计算结束时间为0.18。图1给出计算结束时刻密度曲线图及局部放大图。
图 1 计算结束后密度曲线及局部放大图Fig. 1 Density curve and partial enlarged details
该算例初始条件如式(31)所示[19],网格数为400,计算结束时间为0.038。图2给出计算结束时刻密度曲线图及局部放大图。
该算例初始条件如式(32)所示,网格数为400,计算结束时间为1.8。图3给出计算结束时刻密度曲线图及局部放大图。
根据上述计算结果,评估各格式的计算性能发现:WENO-JS7格式计算性能最优,WENO-JS5格式次之,WENO-JS3格式计算性能最低。即高精度WENO格式对含激波间断流场具有更低的耗散特性和更高的分辨率。
图 2 计算结束后密度曲线及局部放大图Fig. 2 Density curve and partial enlarged details
图 3 计算结束后密度曲线及局部放大图Fig. 3 Density curve and partial enlarged details
为了考察不同精度WENO格式对舱室内爆载荷的影响规律,采用所开发的三维程序开展封闭舱室和泄压舱室内爆炸载荷数值计算。
4.1.1 初始条件设置
密闭舱室尺寸为800×800×800 mm,壁面上设置2个测点,分别编号为NO1,NO2对爆炸超压时间历程进行输出(见图4)。
图 4 舱室及测点分布Fig. 4 Cabin and gauging points
球形炸药(200 g)位于舱室中间,等效后的高压、高密度气体参数:半径为61.65 mm,密度为203.75 kg/m3,压力为3.791×108Pa。计算初始条件见图5(a)。考虑计算时间及精度的要求,经多次数值试验,网格数为40×40×40,见图5(b)。壁面边界条件设置为反射边界。
4.1.2 封闭舱室内爆炸载荷特性
图6给出封闭舱室内爆炸壁面测点NO1和NO2超压时间历程曲线。从图6可知,封闭舱室内爆炸载荷主要包含多峰值、强度逐渐衰弱的冲击波和持续时间较长、保持稳定的准静态超压。根据测点NO1和NO2的对比可知,测点的位置对爆炸前期爆炸波超压峰值有较大的影响,而对最终形成的准静态超压峰值几乎没有影响,说明封闭舱室内部爆炸形成的准静态超压峰值在空间上近似均匀。
4.1.3 不同精度爆炸载荷比较
为了探讨不同精度WENO格式对封闭舱室内爆炸载荷影响规律,选取壁面典型测点NO1对爆炸载荷进行输出。
图7(a)~图7(c)给出不同精度WENO格式,测点NO1爆炸超压时间历程曲线,通过对比发现,高精度WENO格式能分辨出更精细的爆炸载荷特征,即高精度WENO格式具有更低的耗散特性和更高的分辨率。
图 5 爆炸初场及网格分布Fig. 5 Initial condition and mesh distribution
图 6 封闭舱室内爆炸超压时间历程曲线Fig. 6 Overpressure time histories of all the gauging points in closed cabin
图7(d)和图7(e)分别给出3种精度WENO格式爆炸载荷对比图及局部放大图。从图7(d)和图7(e)可知,随着WENO格式精度的提高,爆炸冲击波载荷峰值得到较大的提高,即高精度WENO格式对于激波间断具有更低的耗散特性。然而,不同精度WENO格式对最终形成的准静态超压峰值影响较小,见图7(d)中的粗实线。
4.2.1 初始条件设置
计算初始条件同4.1.1小节,唯一的不同是,泄压舱室壁面上设置一个边长为160 mm的方形泄压口(见图8)。壁面边界条件设置为反射边界,泄压口处边界条件设置为透射边界。
图 7 封闭舱室内爆炸不同精度WENO格式测点NO1爆炸超压时间历程曲线Fig. 7 Overpressure time histories of gauging point NO1 for different schemes in closed cabin
图 8 舱室及测点分布Fig. 8 Cabin and gauging points
图 9 泄压舱室内爆炸超压时间历程曲线Fig. 9 Overpressure time histories of all the gauging points in venting cabin
4.2.2 泄压舱室内爆炸载荷特性
图9给出泄压舱室内爆炸壁面上测点NO1和NO2超压时间历程曲线。从图9可知,泄压舱室内爆炸载荷主要包含多峰值、强度逐渐衰弱的冲击波和持续时间较长、呈现指数衰减规律[20]的准静态超压。根据测点NO1和NO2的对比可知,测点的位置对爆炸前期爆炸波超压峰值有一定的影响,而对准静态超压时间历程几乎没有影响,说明泄压舱室内爆炸形成的准静态超压时间历程在空间上近似均匀。
4.2.3 不同精度爆炸载荷比较
为了探讨不同精度WENO格式对泄压舱室内爆炸载荷影响规律,选取壁面典型测点NO1对爆炸载荷进行输出。
图10(a)~图10(c)给出不同精度WENO格式,测点NO1爆炸超压时间历程曲线,通过对比发现,高精度WENO格式能分辨出更精细的爆炸载荷特征,即高精度WENO格式具有更低的耗散特性和更高的分辨率。
图 10 泄压舱室内爆炸不同精度WENO格式测点NO1爆炸超压时间历程曲线Fig. 10 Overpressure time histories of gauging point NO1 for different schemes in venting cabin
图10(d)和图10(e)分别给出3种精度WENO格式爆炸载荷对比图及局部放大图。从图10(d)和图10(e)可知,随着WENO格式精度的提高,爆炸冲击波载荷峰值得到较大的提高,即高精度WENO格式对于激波间断具有更低的耗散特性。然而,不同精度WENO格式对逐渐形成的、呈现指数衰减规律的准静态超压影响较小,见图10(d)中的粗实线。
基于自主开发的高精度舱室内爆炸三维数值计算程序研究了WENO格式精度对舱室内爆炸载荷影响规律,通过本文的研究主要得到以下4点结论:
1)高精度WENO格式对含激波间断流场具有更低的耗散特性和更高的分辨率。
2)WENO格式精度对舱室内爆炸冲击波载荷影响较大,高精度WENO格式给出更陡峭的激波峰值。
3)WENO格式精度对舱室内爆炸准静态超压(封闭舱室内爆炸形成的保持不变的准静态超压、泄压舱室内爆炸形成的呈指数衰减的准静态超压)载荷影响较小。
4)本文所开发的三维程序可用于舰船舱室内爆炸载荷数值计算,为结构抗爆设计提供载荷依据。
[1]WU C, LUKASZEWICZ M, SCHEBELLA K, et al.Experimental and numerical investigation of confined explosion in a blast chamber [J]. Journal of Loss Prevention in the Process Industries, 2013, 26(4): 737–750.
[2]侯海量, 朱锡, 李伟, 等. 舱内爆炸冲击载荷特性实验研究[J]. 船舶力学, 2010, 14(8): 901–907.
HOU Hai-liang, ZHU Xi, LI Wei, et al. Experimental studies on characteristics of blast loading when exploded inside ship cabin [J]. Journal of Ship Mechanics, 2010, 14(8): 901–907.
[3]孔祥韶, 吴卫国, 李俊, 等. 角隅结构对舱内爆炸载荷影响的实验研究 [J]. 中国造船, 2012(3): 40–50.
KONG Xiang-shao, WU Wei-guo, LI Jun, et al. Experimental research of influence of corner structure on blast loading under inner explosion[J]. Chinese Journal of Ship Research, 2012(3):40–50.
[4]侯海量, 朱锡, 梅志远. 舱内爆炸载荷及舱室板架结构的失效模式分析 [J]. 爆炸与冲击, 2007, 27(2): 151–158.
HOU Hai-liang, ZHU Xi, MEI Zhi-yuan. Study on the blast load and failure mode of ship structure subject to internal explosion [J]. Explosion and Shock Waves, 2007, 27(2):151–158.
[5]孔祥韶, 吴卫国, 李晓彬, 等. 舰船舱室内部爆炸的数值模拟研究 [J]. 中国舰船研究, 2009, 4(4): 7–11.
KONG Xiang-shao, WU Wei-guo, LI Xiao-bin, et al.Numerical simualtion of cabin structure under inner explosion[J]. Chinese Journal of Ship Research, 2009, 4(4):7–11.
[6]樊壮卿, 王伟力, 黄雪峰, 等. 典型舱室内爆炸仿真分析 [J].工程爆破, 2015, 21(3): 13–17.
FAN Zhuang-qing, WANG Wei-li, HUANG Xue-feng, et al.Simulation analysis on typical cabin internal explosion[J].Engineering Blasting, 2015, 21(3): 13–17.
[7]陈攀, 刘志忠. 舱室内爆冲击波载荷特性及影响因素分析[J]. 舰船科学技术, 2016, 38(3): 43–48.
CHEN Pan, LIU Zhi-zhong. Research on loading of explosive and influencing factors inside closed cabin[J]. Ship Science and Technology, 2016, 38(3): 43–48.
[8]HARTEN A. High resolution schemes for hyperbolic conservation laws[J]. Journal of Computational Physics, 1983,135(3): 357–393.
[9]HARTEN A, ENGQUIST B, OSHER S, et al. Uniformly high order accurate essentially non-oscillatory schemes, III [J].Journal of Computational Physics, 1987, 71(2): 231–303.
[10]ZHANG H, ZHUANG F. NND schemes and their applications to numerical simulation of two-and three-dimensional flows*[M]//JOHN W H, THEODORE Y W. Advances in Applied Mechanics. Elsevier. 1991: 193–256.
[11]LIU X D, OSHER S, CHAN T. Weighted essentially nonoscillatory schemes [J]. Journal of Computational Physics,1994, 115(1): 200–212.
[12]JIANG G S, SHU C W. Efficient implementation of weighted ENO schemes [J]. Journal of Computational Physics, 1996,126(1): 202–202.
[13]SHU C W. Essentially non-oscillatory and weighted essentially non-oscillatory schemes for hyperbolic conservation laws [M].Springer Berlin Heidelberg. 1998: 325–432.
[14]SHU C W. High order weighted essentially nonoscillatory schemes for convection dominated problems [J]. Siam Review,2009, 51(1): 82–126.
[15]谢昌坦. 气—水可压缩流物质界面的R-M不稳定性研究[D]. 哈尔滨: 哈尔滨工程大学, 2011.
XIE Chang-tan. The study on R-M instability of the material interface in gas-water compressible flow[D]. Harbin: Harbin Engineering University Master of Engineering, 2011.
[16]BALSARA D S, SHU C W. Monotonicity preserving weighted essentially non-oscillatory schemes with increasingly high order of accuracy [J]. Journal of Computational Physics, 2000,160(2): 405–452.
[17]SHU C W, OSHER S. Efficient implementation of essentially non-oscillatory shock-capturing schemes [J]. Journal of Computational Physics, 1988, 77(2): 439–471.
[18]SOD G A. A survey of several finite difference methods for systems of nonlinear hyperbolic conservation laws [J]. Journal of Computational Physics, 1978, 27(1): 1–31.
[19]WOODWARD P, COLELLA P. Colella, P. The numerical simulation of two-dimensional fluid flow with strong shocks[J]. Journal of Computational Physics, 1984, 54(1):115–173.
[20]ANDERSON C E, BAKER W E, WAUTERS D K, et al.Quasi-static pressure, duration, and impulse for explosions (e.g. HE) in structures [J]. International Journal of Mechanical Sciences, 1983, 25(6): 455–464.