吕海成,黄孝龙,李 宁,翁春生
(南京理工大学瞬态物理国家重点实验室,江苏 南京 210094)
混合气体爆轰作为一种特殊的气体燃烧现象,在航空发动机和燃气泄漏爆炸等场景广泛出现。由于气体爆轰在极短的时间内释放出大量能量,当爆轰波冲击周围固体介质(如容器壁,活塞)时,会产生较大的冲击力。对于发动机等存在工作循环的动力机械,这样的冲击将造成机械零件的应力疲劳、热疲劳并最终导致失效[1];在地下矿井、燃气管道等场合,混合气体爆轰会对周围的结构和环境产生巨大的破坏[2-3];这种冲击也在诸如轻气炮,爆炸喷涂等特殊场景中具有应用价值。因此,研究混合气体爆轰冲击气固界面的动力学响应有很强的工程意义。
国内外针对相关问题展开了很多研究。在爆轰波冲击界面理论研究方面,Lieberman 等[4]利用实验观察了爆轰波从可燃气体透射到不可燃气体混合物的过程,分析了入射波曲率、隔膜和二次燃烧等对该过程波结构的影响。于明等[5]针对凝聚炸药爆轰波在高声速材料界面上的折射,利用改进的爆轰波极曲线理论展开理论分析,结合数值模拟研究了典型的爆轰波折射过程的相关现象。Tang 等[6]研究了爆轰波在不同密度界面上的传播,并利用数值模拟讨论了密度差别大小和反应动力学对爆轰波传播的影响。Peace 等[7]对氢氧爆轰波与界面的一维相互作用进行了分析和实验研究。发现透射激波的放大或衰减取决于界面反射类型和介质的声阻抗比,Damazo 等[8]实验观察了气体冲击波撞击固体的现象,并对其中波速变化和反应区的关系给出了理论解释。
在气体爆轰波冲击固体介质的工程应用方面,张洋溢等[9]将爆轰波斜冲击金属介质理论引入聚能药型罩的罩高参数设计中,结合数值计算得出了不同形状药型罩在点起爆条件下罩高参数的工程算法。周宁等[10]实验研究了闭口管道中甲烷/空气混合气爆炸产生的反射波造成的管壁动态响应。郑龙席等[11]建立了爆震燃烧室的有限元模型并通过加载实验获取真实爆震载荷,计算分析了多循环工作条件下不同因素对爆震燃烧室等效应力的影响,并以此为基础对等壁厚爆震燃烧室进行了壁厚优化。Li 等[12]实验分析了对氢气/空气混合物由管道向圆柱容器泄爆的现象,并探究了通风口爆破压力和管道内障碍物对爆轰波传播的影响。Portaro 等[13]利用气体爆轰作为动力,提出了一种新型的无针注射器的动力装置,并利用实验研究了装置工作过程中气体爆轰作用在注射装置过程的相关规律。Wang 等[14]实验研究了利用甲烷/氧气爆轰作为地震波激发源冲击岩土体,分析了由混合气爆轰激发的地震波的频率和激发距离。Du 等[15]对气体爆轰与管道动态裂纹扩展的过程进行了数值模拟,讨论了爆轰过程的应力变化等参数对裂纹扩展程度的影响。高健等[16]通过实验研究了乙炔/氧气混合物爆轰波对活塞材料的破坏,探究了初始压力及当量比变化对破坏力的影响。
目前,针对气体爆轰冲击气固界面的工程研究多是将爆轰载荷直接输入固体进行模拟,由于缺少理论模型以及气固界面双向耦合的数值模拟方法,因此没有考虑到爆轰波反射到达界面后产生的高压。本文中,针对此问题,使用冲击波相关理论建立气相爆轰波冲击气固界面的一维简化理论模型,探究该过程特征点物理量的规律,采用CE/SE 方法进行气体爆轰冲击气固界面模拟,得到该过程中气体部分压力变化规律及固体中应力波的相关特征,并通过实验验证简化模型和数值模拟的正确性,以期为气体爆轰有关的应力分析以及相关诊断和探测领域提供一定参考。
图1 所示为气体激波冲击气固界面的一维模型示意图,左侧的理想气体柱和右侧的无限长固体杆紧贴在一起,气体柱和固体杆垂直x方向面积为δA。设左侧气体和右侧固体杆初始时刻均处于自然状态,初始速度U0=0,固体中应力和气体气压相同,即σ0=p0。由左侧传来一道激波,到达气固界面时,会瞬间产生两道波:反射波和透射波,由连续性条件可知,气固界面两侧的物理量满足
图1 气体激波冲击气固界面一维模型Fig. 1 A one-dimensional model of gas shock wave impacting solid
式中:p1为气体侧压力,σ1为固体侧应力,U1、U2分别为气体和固体侧的物质速度。
金属、陶瓷等材料中的应力波特性和应力大小有关。区分材料弹性阶段和塑性阶段的应力阈值被称为Hugoniot 弹性极限(Hugoniot elastic limit, HEL)[17]。应用在工程中的可燃混合气爆轰波撞击壁面的压力一般在5~200 MPa 之间,而金属,陶瓷等固体材料的HEL 都远超过这个值,如铝的HEL 约为0.6 GPa,铁的HEL 为1~1.5 GPa。因此,固体内的应力在气体爆轰波冲击界面的过程中均处在弹性阶段,可认为气相爆轰波冲击气固交界面产生的冲击波为弹性波。
根据一维应力波的动量守恒条件,在固体内部和气固边界处满足如下p-U关系:
由激波管产生的激波和由混合气体在一端封闭管道产生的爆轰波有很多不同之处,主要表现在:(1)激波波后为一稳态区,一端封闭管道产生的爆轰波波后会紧跟一个稀疏波;(2)爆轰波在到达稳定状态(C-J 状态)之前存在化学反应过程,前导激波后存在一个稳定的高压区,被称为von Neumann 尖峰。
爆轰波的von Neumann 尖峰同燃烧后的C-J 状态相比,其温度、压力和成分存在很大差别,但爆轰波化学反应诱导时间很短,只有微秒量级[18],在工程应用中可以忽略不计。因此,可以将一端封闭管道产生的爆轰波看作一道波后为C-J 状态的激波紧跟一道稀疏波,爆轰波冲击气固界面的问题也得到了相应简化。面,反射膨胀波经过不断发展,当入射膨胀波完全到达交界面后,反射膨胀波对应图中曲线,其波尾状态为状态4。实际情况中,爆轰波到达界面后,化学反应区域消失,反射稀疏波会追赶上反射激波并对其进行卸载。但该卸载过程时间很长,在界面到达稳定状态之前进行程度较低,对以上过程影响较小。
图2 爆轰波冲击界面的x-t 图和p-U 图Fig. 2 The x-t diagram and p-U diagram of detonation wave impact the interface
固体介质初始时刻同未燃混合气接触且处于静止状态,其p-U关系满足
由两曲线交点可求出状态1。
爆轰波后的膨胀波,其波头状态为C-J 状态,其p-U关系满足
其波尾状态2 的物质速度为0,由此可求出状态2 的相关参数。
其波尾位于气固边界,满足固体介质p-U关系,由此可求得最终气固边界状态4。
由于以上激波及膨胀波全都位于已燃区,化学反应结束,且温度变化不大,因此以上p-U曲线中的绝热指数γ 基本相同,均可以C-J 状态下的γ 进行计算。
为进一步探究可燃混合气爆轰波冲击气固界面过程中各物理量的变化,建立一维流固耦合数值模型。气体部分模型采用欧拉网格,计算方法采用时空守恒元和求解元(CE/SE)方法,将时间与空间统一处理,具有计算格式简单,精度高等优点,能够捕捉激波等强间断。为更加确切描述可燃混合气的爆轰波的产生和发展过程,本文采用基元反应描述可燃混合气的燃烧过程[20]。固体部分采用拉格朗日网格,采用有限元方法进行计算,气体和固体网格之间采用浸没边界法(immersed boundary method,IBM)[21]进行相互耦合。
气体部分以氢氧可燃混合气爆轰为例,采用7 组分16 反应的基元反应模型,模型已被证实可用于模拟氢氧爆轰[22]。根据理论分析的结论,在气相爆轰波冲击气固交界面过程中固体材料基本表现为线弹性特性,不同材料性质不同主要表现为材料的声阻抗不同。工程中参与爆轰波冲击的材料一般为金属,其声阻抗基本都位于1×107~4×107kg/(m2·s)之间,相差不大,金属材料不同对反射激波峰值压力变化影响较小。因此本部分模拟固体部分选择弹塑性材料模型,材料选用铝,密度为2 700 kg/m3,弹性模量为68 GPa,屈服极限为0.325 GPa。
数值模型如图3 所示,左侧流体区域长0.05 m,为了浸入边界法计算的需要,在气固界面处,流体区域网格和固体区域网格有部分重合,图中灰色重合区长度0.01 m。右侧固体区域长0.05 m,计算域总长0.11 m。
图3 一维数值计算模型Fig. 3 The one-dimensional numerical simulation model
初始条件和边界条件如上图所示,为产生稳定的混合气爆轰波,气体段初始填充化学恰当比的氢氧混合气,温度300 K,压力0.1 MPa。为引爆混合气,在气体段左侧封闭端处取一小区域设为点火区域,温度3 000 K,压力1 MPa。左侧气体封闭端采用固壁边界条件,右侧固体边界采用无反射边界条件。
为确定网格大小对于模拟结果的影响,以氢氧可燃混合气爆轰波产生并稳定为例设计5 组算例。设定网格大小从0.1 mm 到0.004 mm依次递减,模拟结果如表1 所示。可以看到,网格尺寸到达0.01 mm 后,网格尺寸再减小基本不对爆轰波波后压力和波速产生影响,为减小计算量,选取气体部分网格尺寸为0.01 mm。一般运用浸入边界法进行流固耦合计算时,固体网格尺寸一般为气体网格尺寸的2 倍。因此,数值模拟选取的固体网格尺寸为0.02 mm。
表1 网格收敛性模拟结果Table 1 Simulation results of mesh convergence
初始时刻,气体区域的可燃混合气在左端高温高压段的作用下开始反应,形成爆轰波并向右传播。图4选取了爆轰波形成和传播过程中不同时刻气体部分压力分布。可以看到,燃烧产生的激波开始时压力迅速上升,到t=1.5 µs 后压力缓慢增大,在t=9 µs 时刻之后爆轰波保持了波后压力稳定,此时,波速也稳定为2 833 m/s,和NASA CEA 软件[23]计算的C-J 理论波速2 835 m/s 基本一致,因此可认为气体段已形成稳定爆轰。
图4 气体部分爆轰波的形成与传播Fig. 4 Formation and propagation of detonation wave in gas
在t=18.1 µs 时刻,爆轰波波面到达气固界面。图5 所示为t=18.1,20.5,25.5,31.5,42.5 µs时刻,气体部分的压力分布和固体部分的主应力分布。从图5 中所示的气体侧可以看到,爆轰波到达气固界面(x=0.05 m)后,在界面处产生一高压点,紧接着气体部分反射激波左传,在固体部分弹性波右传。反射激波和原本爆轰波后的泰勒稀疏波相交,波后压力迅速下降。稀疏波在界面上也发生反射,两反射稀疏波相互作用,在反射激波后产生部分压力相对均匀的区域。反射激波完全通过稀疏波后,以基本不变强度左传,其后出现相交后的反射稀疏波。同时,随着原稀疏波波尾右传,稀疏波相交区域缩小,在t=42.5 µs时刻,原稀疏波波尾已完全到达气固界面,气固界面变为稳定状态,由于固体的声阻抗远大于气体的,与固壁性质相近,气固界面稳定状态压力同封闭端压力基本相同,但会有一个向右运动的速度。在该过程中,气固界面两侧始终保持压力相等,由于气固边界处产生压力峰值较小,气体压力变化透射入气固界面后在固体中产生弹性波。气固界面在开始时应力突变至峰值,对应气体侧反射激波波后到达压力峰值状态,随后应力以指数形式衰减,最后应力到达稳定状态保持不变,对应气体侧反射膨胀波波尾达到稳定状态。固体内应力波在传播过程中波形基本保持不变,如图5 固体侧所示。
图5 不同时刻气体段的压力分布和固体段的应力分布Fig. 5 Pressure distribution in gas section and stress distribution in solid section at different times
图6 为爆轰波波面到达气固界面后,气体中反射激波和固体中应力波波速随时间的变化,其中,气体中反射激波波速根据反射激波波面移动距离得到。可以看到,形成的反射激波在t=18.6 µs时刻与原稀疏波相交,此时,反射激波波速约为1400 m/s。随着该过程进行,反射激波波前压力下降,但其实际波速上升。这是因为该区域为反射激波与原稀疏波交会区域,激波与当地气体速度反向。初始交汇时刻,稀疏波波头位置气体速度大,随着反射激波向左侧传播,波前气体速度逐渐降低,反射激波速度逐渐上升。在t=25.3 µs时刻,反射激波到达原稀疏波波尾,此时激波波速到达最大值,约为2100 m/s,接下来反射激波由于受反射稀疏波的追赶卸载,波速缓慢下降。固体中应力波波速恒定为5000 m/s,符合金属铝的理论一维弹性波速5018.5 m/s。
图6 气体反射激波和固体弹性波波速曲线Fig. 6 Velocity of the reflected shock wave in gas and stress wave in solid
为验证理论及数值模拟结果,设计气相爆轰波冲击活塞实验装置,通过测量靠近活塞位置的压力及活塞受到的冲击力,得到气相爆轰波冲击气固界面过程的实际情况。分析气体段和固体段的压力变化特征,并与理论计算和数值模拟结果进行进一步印证。实验系统由爆轰管、供气系统、点火系统和采集系统等组成。如图7 所示。
图7 中爆轰管直径30 mm,点火点位置后爆轰管长1.2 m,利用分压填充控制填充气体比例,实验采用化学恰当比的氢气/氧气混合物,总填充压力为0.12 MPa。气体填充完成后,静置2 h 进行充分混合,由于爆轰管气密性等问题,经2 h 预混后,爆轰管内压力与大气压基本一致,可认为试验为0.1 MPa 爆轰。之后,高能点火器点火,通过扰流片等结构实现爆燃转爆轰过程,形成稳定爆轰波传播至爆轰管末端。在爆轰管末端使用铝制活塞封闭,在靠近活塞的爆轰管壁布置2 个PCB113B24 型高频动态压力传感器监测爆轰管末端交界面压力,其满量程为6.895 MPa;在活塞底部布置PCB208C05 型力传感器监测活塞底部冲击,其满量程为22.24 kN。活塞通过螺栓与力传感器相连,力传感器固定于基座上。2 个压力传感器间距5 cm,靠近底部压力传感器距活塞顶部平面12 mm,活塞长48 mm。实验过程中,3 个传感器同步进行测试,采集频率为1 MHz。
图7 实验系统示意图Fig. 7 Schematic diagram of the experimental system
图8 为爆轰波到达交界面时刻前后的压力数据及活塞冲击数据。爆轰波在t=421.12 ms时刻到达传感器P1,波后压力为2.03 MPa;在t=421.14 ms 时刻到达传感器P2,波后压力为1.98 MPa。经过计算可得爆轰波波速为2 941.17 m/s。可以看到,爆轰波波后压力基本保持不变,且波速与该条件下C-J 理论爆轰波速2 835 m/s 相近,由此可以认为爆轰波到达爆轰管末端时以发展为稳定爆轰状态。
图8 管内气体压力和活塞冲击力曲线Fig. 8 Variation of gas pressure in the pipe and impact force on the piston
爆轰波到达界面后,在气体内反射一道激波并在固体中透射一道应力波。反射激波在t=421.15 ms 时刻到达传感器P2,波后压力为4.11 MPa;在t=421.19 ms 时刻到达传感器P1,波后压力为3.15 MPa。经过计算可得反射激波波速1 515 m/s,与数值模拟得到的反射激波初期的速度基本一致。另一方面,反射激波通过2 个压力传感器后,可以看到两压力传感器测得的波形基本重合,呈指数形态。这符合数值模拟结果中,反射激波后的原稀疏波和反射稀疏波相交区压力均匀分布的现象。随着原稀疏波波尾到达交界面,稀疏波相交区消失,气固界面到达稳定状态,表现在压力传感器的数据中则是压力逐渐下降到约0.7 MPa 并保持稳定。
从力传感器测得的数据中可以看出,应力波在t=421.15 ms 时刻到达,在t=421.19 ms 时刻形成峰值,峰值为3.345kN,由此可得到活塞中应力波压力峰值为4.73 MPa,远小于铝材料的HEL。根据压力传感器得到的数据可估算出,由爆轰波到达界面到弹性波到达力传感器时间间距为0.01 ms,由此得到活塞中应力波波速为5 333 m/s,与铝材料声速5 103 m/s 基本符合。可以判定,气相爆轰波冲击铝活塞,在活塞中产生了弹性波。力传感器的数据到达峰值后不是以指数形式下降而是呈现不规则波动,其可能原因主要有:(1)爆轰波冲击活塞表面后,由于活塞通过力传感器固定于基座上,应力波在透过力传感器和基座时产生反射,使得应力波在活塞中反复震荡导致应力不稳定[24];(2)活塞长径比比较小,爆轰波冲击活塞表面后产生了斜应力波导致应力不均匀[25]。
根据理论模型可知,对已知成分混合气爆轰波特性及固体材料参数,可以得到爆轰波冲击气固界面过程中,爆轰波到达界面时的压力峰值状态和反射稀疏波离去后的稳定状态,气固界面位置的相关参数,从而得到固体中产生的应力波的大小。这里的峰值状态和稳定状态分别对应理论模型过程的状态1 和状态4,这2 个状态的压力分别为p1和p4。
图9 为氢/氧混合气在不同初始压力下产生的稳定爆轰波冲击铝表面时峰值压力p1和稳定压力p4以及它们同爆压pcj的比值变化情况。可以看到,随着混合气初始压力变高,峰值压力p1和稳定压力p4呈线性增长,它们同爆压的比值均基本保持不变,p1/pcj为2.46,p4/pcj为0.37。
图9 混合气爆轰波冲击界面相关参数随初始压力的变化Fig. 9 Variation of detonation wave shock interface parameters with initial pressure
表2 所示为理论模型、数值模拟和实验得到的峰值压力p1和稳定压力p4。可以看到,峰值压力p1的理论和数值结果存在一定差距。这是由于数值模拟采用化学反应模型,模拟爆轰波中存在von Neumann 峰,与直接使用C-J 理论计算的理论模型相比峰值压力和速度偏大;另一方面,理论结果和实验值符合较好,这是因为固体中透射的von Neumann 峰值压力衰减非常快,固体部分仅考虑一维弹性过程的数值模拟无法表现应力波的衰减。稳定压力p4的理论、数值模拟和实验的结果基本相同。
表2 特征状态压力的理论、数值和实验结果Table 2 Theoretical, numerical and experimental results of characteristic state pressure
本文中基于C-J 爆轰理论和冲击波理论针对气相爆轰波冲击气固界面的过程建立了理论模型,同时,利用基元反应模型进行数值模拟,初步分析了气体中激波和固体中应力波的相关物理量的变化规律。并用实验加以验证。主要结论如下。
(1)气相爆轰波到达气固界面后,在固体中透射一弹性波并在气体中反射一道激波,爆轰波后的稀疏波与反射激波相交,削弱反射激波,最终形成稳定激波回传。气固界面在稀疏波和反射稀疏波的作用下,压力和速度逐渐下降,最终也形成稳定状态;随着混合气初始压力变高,峰值压力p1和稳定压力p4呈线性增长,它们同爆压的比值均基本保持不变,p1/pcj为2.46,p4/pcj为0.37。
(2)气体中反射激波在回传过程中与原稀疏波相交,随着相交过程进行,反射激波波后压力下降,波速变快,原稀疏波与反射稀疏波相交区域压力基本保持均匀;固体中的弹性波在开始时应力突变至峰值,波峰面后应力以指数形式衰减,最后到稳定状态保持不变。