3D Mesoscopic Simulation of Shock Compression Behaviors of Reactive Materials
-
摘要: 为研究反应材料的细观冲击响应行为,采用细观数值模拟方法对Al/PTFE(铝/聚四氟乙烯)与Al/Ni反应材料的冲击压缩行为进行了分析研究。采用球磨混合和静压方法,制备了两种反应材料试样,并使用纳米CT设备获得了材料的细观图像,借助图像处理和网格映射方法,建立了基于材料真实细观结构的三维有限元模型,计算获得的冲击Hugoniot关系与理论结果吻合较好。细观数值模拟结果表明:反应材料的冲击波阵面在细观尺度内呈不平整状态,颗粒相在冲击波作用下主要表现为沿冲击压缩方向的体积压缩和运动;在较高撞击速度下,Al/PTFE反应材料的PTFE基体将发生熔化,而Al/Ni反应材料在本研究范围内保持固态。Abstract: In order to investigate the mesoscopic behavior of Al/PTFE reactive materials under shock loading, the shock response of Al/PTFE and Al/Ni reactive materials were investigated via numerical approaches. Two material samples were fabricated by ball mill mixing and cold isostatic pressing, and a Nano-CT system was employed to obtain their mesoscopic images. The 3D mesoscopic finite-element model based on the real configuration was established with the help of image processing and mesh mapping methods, and the numerical Hugoniot results agreed well with the theoretical results. The numerical modeling results indicated that the shock wave is uneven in the mesoscale, and the metal granules were compressed and moved along the impact direction. The PTFE matrix of Al/PTFE melted at high impact velocity, while the Al/Ni remained solid within the impact velocity range of this study.
-
Key words:
- reactive materials /
- mesoscopic model /
- shock temperature rise /
- image processing
-
不同流体间的扰动界面在冲击波作用下会发生不稳定性增长,称为Richtmyer-Meshkov(RM)不稳定性[1-2]。当RM不稳定性发展进入非线性阶段后,扰动界面会发展为“蘑菇”状的尖钉/气泡结构;之后由于界面两侧切向速度不匹配而产生剪切效应,导致界面上又发生Kelvin-Helmholtz(KH)不稳定性[3],并诱发更多小尺度涡。KH不稳定性的出现会破坏“蘑菇”状界面结构的对称性,使之破碎,并加速湍流混合。流体界面不稳定性和湍流混合涉及许多重要领域,如武器物理、惯性约束聚变、超音速燃烧、天体物理等,因此具有重大的研究意义。
较早之前,对于RM界面不稳定性和湍流混合,研究人员关注的重点是湍流混合区宽度的增长规律及影响因素[4-16]。随着实验诊断技术、高精度数值算法及并行计算技术的发展,研究重点逐渐深入到湍流混合的内在机制和统计规律[17-23]以及以此为基础的湍流混合建模[24]。然而,由于湍流混合的多尺度特性,实验中所能获取的流场信息依然十分有限,无法给出足够多的湍流混合区瞬态数据和精细复杂的混合区流场结构;而数值模拟基本上能够给出所需要的有关湍流混合区流场的所有信息,特别是借助超大规模的高性能计算,可以在不同尺度上研究湍流混合的物理机制及其演化的统计特性。Cabot等[25]就曾在IBM的BlueGene/L超级计算机上利用6万多个核实现了Rayleigh-Taylor(RT)不稳定性及湍流混合的近300亿网格(3 0723)的超大规模直接数值模拟,研究了湍流混合区增长的标度率问题和雷诺数效应。
AWE开展的“反尖端”RM不稳定性激波管实验中初始的SF6/空气右界面呈“∨”形,整体上类似于反的士兵臂章,称为“反尖端”界面。本研究采用自研的可压缩多介质黏性流动和湍流大涡模拟程序MVFT[6, 21],对“反尖端”界面不稳定性及湍流混合进行大规模并行数值模拟,在高分辨率的网格下细致地研究“反尖端”界面不稳定性的发展过程和规律,包括波系作用、湍流混合区结构的演化、湍流混合区统计量的变化等。
1. 计算方法
可压缩多介质黏性流动和湍流大涡模拟程序MVFT所求解的流场控制方程是经过Favre滤波后的可压缩多介质黏性流动Navier-Stokes方程组
{∂ˉρ∂t+∂ˉρ˜uj∂xj=0∂ˉρ˜ui∂t+∂(ˉρ˜ui˜uj+ˉpδij)∂xj=∂σij∂xj−∂τij∂xj∂ˉρ˜E∂t+∂(ˉρ˜uj˜E+ˉp˜uj)∂xj=∂(σij−τij)˜ui∂xj−∂(qlj+QTj)∂xj∂˜Y(s)∂t+˜uj∂˜Y(s)∂xj=∂∂xj(˜D∂˜Y(s)∂xj)−∂QYj∂xj (1) 式中:下标i、j分别代表x、y、z 3个方向,遵循张量运算法则;
ˉρ 、˜uk(k=i,j) 、ˉp 、˜E 分别为可解尺度流体的密度、速度、压强和单位质量总能量;N为流体种类;˜Y(s) 表示第s种流体的体积分数(s=1,2,⋯, N–1),且满足∑Ns=1˜Y(s)=1 ;˜D 为扩散系数,˜D=ν/Sc ,ν 为流体运动黏性,Sc为Schmidt数;qlj=−λl∂˜T/∂xj 是单位时间单位空间的可解尺度能量流,λl=μlcp/Prl 是可解尺度热传导系数,cp为比定压热容,Prl是Prandtl数,˜T 为流体温度;τij,QTi,QYi 为滤波后产生的亚格子输运项,分别表示亚格子应力张量、热流量和体积分数输运通量,需要通过亚格子尺度应力模型进行模化处理以使控制方程组封闭,本研究采用Vreman亚格子模型[26]。牛顿流体黏性应力张量σij 表示为σij=μl[∂˜ui∂xj+∂˜uj∂xi−23δij(∂˜uk∂xk)] (2) 式中:μl为流体黏性。气体介质的状态方程采用理想气体状态方程。
具体计算过程为:首先采用算子分裂技术将(1)式描述的物理过程分解为3个子过程进行计算,即整个通量分解为无黏通量、黏性通量和热通量3部分。无黏通量的计算采用多介质高精度PPM(Piecewise Parabolic Method)方法,两步Lagrange-Remapping型的PPM方法分4个步骤进行:(1)物理量分段抛物插值,(2)近似Riemann问题求解,(3)Lagrange方程组推进求解,(4)将物理量映射到静止的欧拉网格上。然后在无黏通量的基础上,采用二阶空间中心差分方法和两步Runge-Kutta时间推进方法求解黏性通量、热通量及标量输运通量。
2. 结果分析
本研究所采用的模型是英国AWE开展的“反尖端”RM不稳定性激波管实验[27],计算模型如图1所示,由细金属丝制成的网状支架支撑微米量级的硝化纤维薄膜将实验气体隔开,形成预设的初始界面。激波管内截面尺寸为20 cm×10 cm,“反尖端”区气体为SF6,其宽度为15 cm,两侧气体是空气,“反尖端”界面平衡位置距离激波管尾端20 cm。另外,在初始界面上设置振幅小于0.01 cm、波长为0.5~5 cm的小尺度随机扰动,以模拟支架的扰动效应。初始入射冲击波马赫数为1.26。计算域:(x, y, z)∈[0, 80 cm]×[–10 cm, 10 cm]×[0, 10 cm]。实验得到的结果十分有限,只是不同时刻混合区的演化图像及壁面气泡和尖钉位置(距激波管末端的距离),不足以深入分析界面不稳定性及湍流混合发展规律。为此,本研究借助大规模数值模拟,利用1 024个CPU,采用粗(网格大小1.00 mm)、中(网格大小0.50 mm)、细(网格大小0.25 mm)3种不同分辨率的网格,详细研究“反尖端”界面不稳定性及湍流混合的演化特性。细网格下的网格总量为1.024×109。气体初始参数如表1所示,其中γ为气体比热比。
表 1 空气和SF6的初始参数Table 1. Initial properties of air and SF6Gas ρ/(kg·m–3) p/MPa γ μl/(Pa·s) Diffusion coefficient/(m2·s–1) SF6 5.97 0.1 1.09 1.474 6×10–5 0.97×10-5 Air 1.18 0.1 1.40 1.852 6×10–5 2.04×10–5 冲击波从低阻抗介质(如空气)向高阻抗介质(如SF6)方向加载物质界面时,会产生一个透射冲击波和一个反射冲击波;冲击波从高阻抗介质向低阻抗介质方向加载物质界面,会产生一个透射冲击波和一个反射稀疏波。图2给出了以密度显示的“反尖端”计算模型轴线上的一维近似波谱,可以看出复杂的波与界面以及波与波之间的相互作用,冲击波加载使流场密度增大,稀疏波卸载使流场密度减小。首先,入射冲击波从左侧加载空气/SF6左界面,产生一个左行反射冲击波和一个右行透射冲击波;该右行透射冲击波加载SF6/空气右界面,在SF6气体中产生一个左行反射稀疏波,向右边空气中透射一个右行冲击波;右边空气中的右行冲击波会从激波管尾端反射回来再次与SF6/空气右界面作用,向SF6中产生一个左行透射冲击波(会在左界面上发生分解,向SF6中反射一个右行稀疏波),并向右边空气中产生一个右行反射冲击波,该右行反射冲击波会反复从激波管尾端反射回来加载SF6/空气右界面,在界面上发生分解;前面第一次在SF6/空气右界面上分解向SF6中产生的左行反射稀疏波和空气/SF6左界面作用后向SF6中产生一个反射压缩波,该压缩波和SF6/空气右界面作用分解产生一个左行反射稀疏波,之后在SF6左右界面之间反复分解成压缩波和稀疏波。如图2所示,在SF6气体区域的胞格结构显示了冲击波、压缩波、稀疏波、界面之间的复杂相互作用,使得界面不稳定性发展更复杂,而波和界面的每一次作用都会促使混合区快速发展[22]。
图3给出了“反尖端”界面及流场演化的实验图像和数值模拟结果比较,左列是实验结果,右3列从左至右依次是网格为1.00、0.50、0.25 mm时数值模拟的沿展向平均的密度场图像,清晰地显示了波和界面的相互作用及波系在流场中的相互作用和演化过程。实验和数值模拟二者之间,包括界面的形状和位置及波阵面的形状和位置均吻合很好。当冲击波与SF6区域左右界面作用后,界面上的小尺度随机扰动开始增长;而且在冲击波与SF6区域的右界面作用时,由于冲击波是从高阻抗介质向低阻抗介质方向加载界面,所以右界面的大尺度“反尖端”界面发生反相,从中心轴处生长出大尺度尖钉,靠近激波管上下壁处形成大尺度气泡结构,并且尖钉和气泡逐渐长大。对比3种不同网格分辨率数值模拟的展向平均密度场可以看出,界面大尺度结构和波系的演化差别很小,只是在中后期界面小尺度结构之间有较小的差别,这是由于展向平均抹去了一些细节结构造成的。图4、图5和图6分别给出了2.0、3.0、4.0 ms时刻不同网格分辨率下以SF6体积分数(YSF6)显示的湍流混合区三维图像。图中显示了复杂的三维湍流混合区演化,包括小尺度气泡/尖钉的长大、变形、融合等;不同分辨率网格计算所得的湍流混合区图像差别很大,其中高分辨率网格捕捉到了更精细的湍流混合区结构。图7给出了大尺度壁面气泡和中心尖钉位置随时间变化的曲线,该位置是相对于激波管尾端的距离,不同颜色曲线表示不同网格分辨率下的数值计算结果。可以看出,在不同网格分辨率下计算所得的大尺度壁面气泡和中心尖钉的位置差别较小,这是因为大尺度壁面气泡和尖钉的位置由湍流混合区的包络线确定,而湍流混合区的包络线主要由大尺度结构决定,湍流混合区大尺度结构的平均特征在不同网格分辨率下的差别很小;虽然差别较小,但是细网格下的模拟结果与实验结果吻合得更好。
图8和图9显示了1.0、2.0、3.0、4.0 ms时中心轴上流场密度和SF6体积分数在不同网格分辨率下的分布。图中清晰地显示了冲击波的传播及其在界面上的分解,以及界面的发展演化,比如1.0 ms时刻的密度分布显示,冲击波加载使流场密度增大,SF6气体中的右行透射冲击波与SF6/空气界面作用后向SF6气体中反射的左行稀疏波的卸载效应使流场密度减小。在前期(1.0~2.0 ms),不同网格分辨率计算的流场密度和SF6体积分数的分布差别很小,说明混合区流场湍流发展不充分,小尺度脉动还未完全发展起来;后期(3.0~4.0 ms)不同网格分辨率的计算结果差别增大,表明在冲击波的多次复杂加载作用下混合区流场已经进入完全湍流阶段,更多小尺度脉动出现,细网格对这些小尺度含能涡结构有更高的分辨率,而粗网格将大量的小尺度含能涡结构耗散掉。
湍动能和拟涡能可以用来表征流场湍流脉动和湍涡的强度和发展。图10和图11分别为1.0、2.0、3.0、4.0 ms时刻无量纲化的流场湍动能K和拟涡能Ω在不同网格分辨率下沿冲击波运动方向的分布,无量纲化速度参数取初始入射冲击波后的流场速度133.6 m/s,无量纲化长度参数取计算模型宽度0.1 m。可以看出:随着混合区的发展,小尺度脉动逐渐发展起来;随着网格分辨率的增大,对这些小尺度脉动的捕捉能力也逐渐提高,捕捉到了更强的湍流脉动和湍涡。所以更小尺寸的网格对流场有更高的分辨率,可以捕捉到更精细的混合区湍流结构,这对于研究RM不稳定性及湍流混合演化的动力学行为非常重要。
3. 结 论
利用可压缩多介质黏性流动和湍流大涡模拟程序MVFT,对“反尖端”界面不稳定性及其诱发的湍流混合进行了大规模并行三维数值模拟分析,揭示了界面不稳定性及湍流混合的复杂发展过程和规律。冲击波和空气/SF6/空气界面作用并发生分解而产生冲击波、稀疏波、压缩波,这些次生波在SF6气体中运动并相互作用,而且多次加载物质界面以及波和界面的每一次作用都会加速湍流混合区的发展和物质混合。高/低阻抗构型的“反尖端”界面受冲击加载后发生反相而发展为大尺度的壁面气泡和中心轴尖钉结构,该大尺度气泡和尖钉结构基本确定了湍流混合区的平均几何特征和包络范围而不依赖计算网格。采用高分辨率的计算网格能够捕捉到更精细的混合区小尺度湍涡结构及更强的湍流脉动,显示出湍流混合区的复杂结构和特征。
-
表 1 研究所用原材料及规格
Table 1. Details of the materials used in this research
Materials Sizes/meshes Manufacturers Al 100−200 Sinopharm Chemical Reagent Co., Ltd., Shanghai, China Ni 150−200 Aladdin PTFE 500 Shandong Fluorine Chemical Co., Ltd., Shandong, China Materials A/MPa B/MPa C m n Tm /K Al 265 426 0.015 1.000 0.34 775 PTFE 11 44 0.120 1.000 1.00 350 Ni 163 648 0.330 0.006 1.44 1728 -
[1] 周杰, 何勇, 何源, 等. 含能毁伤元冲击引爆模拟战斗部试验研究 [J]. 含能材料, 2016, 24(11): 1048–1056.ZHOU J, HE Y, HE Y, et al. Experimental study on shock initiation of simulative warhead by energetic kill element [J]. Chinese Journal of Energetic Materials, 2016, 24(11): 1048–1056. [2] 王海福, 郑元枫, 余庆波, 等. 活性破片引燃航空煤油实验研究 [J]. 兵工学报, 2012, 33(9): 1148–1152.WANG H F, ZHENG Y F, YU Q B, et al. Experimental research on igniting the aviation kerosene by reactive fragment [J]. Acta Armamentarii, 2012, 33(9): 1148–1152. [3] 徐松林, 阳世清, 张炜, 等. PTFE/Al含能复合物的本构关系 [J]. 爆炸与冲击, 2010, 30(4): 439–444.XU S L, YANG S Q, ZHANG W, et al. A constitutive relation for a pressed PTFE/Al energetic composite material [J]. Explosion and Shock Waves, 2010, 30(4): 439–444. [4] 徐松林, 阳世清, 徐文涛, 等. PTFE/Al反应材料的力学性能研究 [J]. 高压物理学报, 2009, 23(5): 384–388.XU S L, YANG S Q, XU W T, et al. Research on the mechanical performance of PTFE/Al reactive materials [J]. Chinese Journal of High Pressure Physics, 2009, 23(5): 384–388. [5] 赵鹏铎, 卢芳云, 李俊玲, 等. 活性材料PTFE/Al动态压缩性能 [J]. 含能材料, 2009, 17(4): 459–462.ZHAO P D, LU F Y, LI J L, et al. The dynamic compressive properties of PTFE/Al reactive materials [J]. Chinese Journal of Energetic Materials, 2009, 17(4): 459–462. [6] 周杰, 何勇, 何源, 等. Al/PTFE/W 反应材料的准静态压缩性能与冲击释能特性 [J]. 含能材料, 2017, 25(11): 903–912.ZHOU J, HE Y, HE Y, et al. Quasi-static compression properties and impact energy release characteristics of Al/PTFE/W reactive materials [J]. Chinese Jouranal of Energetic Materials, 2017, 25(11): 903–912. [7] XIONG W, ZHANG X F, WU Y, et al. Influence of additives on microstructures, mechanical properties and shock-induced reaction characteristics of Al/Ni composites [J]. Journal of Alloys and Compounds, 2015, 648: 540–549. doi: 10.1016/j.jallcom.2015.07.004 [8] XIONG W, ZHANG X F, TAN M T, et al. The energy release characteristics of shock-induced chemical reaction of Al/Ni composites [J]. Journal of Physical Chemistry C, 2016, 120(43): 24551–24559. doi: 10.1021/acs.jpcc.6b06530 [9] EAKINS D E, THADHANI N N. Mechanistic aspects of shock-induced reactions in Ni+Al powder mixtures [C]//Shock Compression of Condensed Matter–2007. Waikoloa, HI, 2007. [10] AYDELOTTE B B, THADHANI N N. Mechanistic aspects of impact initiated reactions in explosively consolidated metal+aluminum powder mixtures [J]. Materials Science & Engineering A, 2013, 570: 164–171. [11] QIAO L, ZHANG X F, HE Y, et al. Mesoscale simulation on the shock compression behaviour of Al-W-Binder granular metal mixtures [J]. Materials & Design, 2013, 47: 341–349. [12] GE C, MAIMAITITUERSUN W, DONG Y X, et al. A study on the mechanical properties and impact-induced initiation characteristics of brittle PTFE/Al/W reactive materials [J]. Materials, 2017, 10(5): 452. doi: 10.3390/ma10050452 [13] SUN Y, LI Q M, LOWE T, et al. Investigation of strain-rate effect on the compressive behaviour of closed-cell aluminium foam by 3D image-based modelling [J]. Materials & Design, 2016, 89: 215–224. [14] 褚巧龙. 基于Otsu的图像阈值分割算法的研究[D]. 秦皇岛: 燕山大学, 2011.CHU Q L. Research on image threshold segmentation arithmetic based on Otsu [D]. Qinhuangdao: Yanshan University, 2011. [15] 袁小翠, 黄志开, 马永力, 等. Otsu阈值分割法特点及其应用分析 [J]. 南昌工程学院学报, 2019, 38(1): 89–94.YUAN X C, HUANG Z K, MA Y L, et al. Analysis of characteristics and application of Otsu threshold method [J]. Journal of Nanchang Institute of Technology, 2019, 38(1): 89–94. [16] JOHNSON G R, COOK W H. A constitutive model and data for metals subjected to large strains, high strain rates and high temperatures [C]//7th International Symposium on Ballistics. The Hague, Netherlands, 1983. [17] SYSTÈMES D. Abaqus analysis user’s guide [M]. Providence, RI: Dassault Systèmes Simulia Corp, 2014. [18] HERBOLD E B, NESTERENKO V F, BENSON D J, et al. Particle size effect on strength, failure, and shock behavior in polytetrafluoroethylene-Al-W granular composite materials [J]. Journal of Applied Physics, 2008, 104(10): 103903. doi: 10.1063/1.3000631 [19] 汤文辉, 张若棋. 物态方程理论及计算概论[M]. 2版. 北京: 高等教育出版, 2008.TANG W H, ZHANG R Q. Introduction to theory and calculation of equation of state [M]. 2nd ed. Beijing: Higher Education Press, 2008. [20] 周杰. 典型氟聚物基活性材料冲击反应特性研究[D]. 南京: 南京理工大学, 2018.ZHOU J. Study on the impact-induced reaction characteristics of typical fluoropolymer-matrix reactive materials [D]. Nanjing: Nanjing University of Science and Technology, 2018. [21] GUO J, ZHANG Q M, ZHANG L S, et al. Reaction behavior of polytetrafluoroethylene/Al granular composites subjected to planar shock wave [J]. Propellants Explosives Pyrotechnics, 2017, 42(3): 299–306. [22] EAKINS D E, THADHANI N N. Mesoscale simulation of the configuration-dependent shock-compression response of Ni plus Al powder mixtures [J]. Acta Materialia, 2008, 56(7): 1496–1510. doi: 10.1016/j.actamat.2007.12.009 [23] EAKINS D E, THADHANI N N. Shock compression of reactive powder mixtures [J]. International Materials Reviews, 2009, 54(4): 181–213. doi: 10.1179/174328009X461050 [24] ZHANG X F, SHI A S, ZHANG J, et al. Thermochemical modeling of temperature controlled shock-induced chemical reactions in multifunctional energetic structural materials under shock compression [J]. Journal of Applied Physics, 2012, 111(12): 123501. doi: 10.1063/1.4729048 [25] ZHOU J, HE Y, HE Y, et al. Investigation on impact initiation characteristics of fluoropolymer-matrix reactive materials [J]. Propellants, Explosives, Pyrotechnics, 2017, 42(6): 603–615. doi: 10.1002/prep.201700003 [26] ZHANG X F, QIAO L, SHI A S, et al. A cold energy mixture theory for the equation of state in solid and porous metal mixtures [J]. Journal of Applied Physics, 2011, 110(1): 013506. doi: 10.1063/1.3603018 [27] 何源, 何勇, 张先锋, 等. 疏松金属材料冲击温度理论分析 [J]. 爆炸与冲击, 2012, 32(2): 143–149.HE Y, HE Y, ZHANG X F, et al. Theoretical analysis on shock temperature of porous metal [J]. Explosion and Shock Waves, 2012, 32(2): 143–149. [28] SONG P, CAI L C, WANG Q S, et al. Sound velocity, temperature, melting along the Hugoniot and equation of state for two porosity aluminums [J]. Journal of Applied Physics, 2011, 110(10): 103522. doi: 10.1063/1.3662193 [29] QIAO L R, ZHANG X F, HE Y D, et al. Multiscale modelling on the shock-induced chemical reactions of multifunctional energetic structural materials [J]. Journal of Applied Physics, 2013, 113(17): 173513. doi: 10.1063/1.4803712 [30] SRAJ I, SPECHT P E, THADHANI N N, et al. Numerical simulation of shock initiation of Ni/Al multilayered composites [J]. Journal of Applied Physics, 2014, 115(2): 023515. doi: 10.1063/1.4861402 [31] 谈庆明, 黄风雷. 超高速碰撞动力学引论[M]. 北京: 科学出版社, 2000.TAN Q M, HUANG F L. Introduction to the dynamics of hypervelocity impact [M]. Beijing: Science Press, 2000. [32] 夏睿全, 张小平. 聚四氟乙烯废料的热解实验 [J]. 化工进展, 2008, 27(1): 98–103.XIA R Q, ZHANG X P. Experimental study on waste polytetrafluoroethylene pyrolysis [J]. Chemical Industry and Engineering Progress, 2008, 27(1): 98–103. [33] 李玲琴. 金属/氟聚物反应材料性能的研究[D]. 太原: 中北大学, 2015.LI L Q. Research on detonation properties of metal-fluoride reactive materials [D]. Taiyuan: North Central University, 2015. -