Structure and propagation mode of gaseous spinning detonation in rectangular tube
-
摘要: 为探索极限条件下矩形管道截面长宽比对于螺旋爆轰传播的影响,采用基于五阶WENO有限差分格式和两步总包反应模型的Euler方程,对三维气相螺旋爆轰波在矩形截面管道中的结构及其传播方式进行了数值研究。通过模拟不同管道截面尺寸下爆轰波的三波线运动轨迹、流场分布及高压印记结构,揭示了截面几何尺寸对气相临界爆轰波稳定传播的影响规律。结果表明:螺旋爆轰能在一定范围的小管道截面尺寸内通过横、竖两条三波线及其相互作用形成的斜三波线的运动来维持传播;随着管道截面尺寸长宽比的增加,螺旋爆轰在壁面上形成的高压印记逐渐由倾斜的条带结构变成局部点状分布结构,波阵面上的斜三波线的轨迹也由方管中沿着单一方向的圆周运动逐渐发展为具有转向机制的复杂运动轨迹;当长宽比进一步增加时,三维螺旋爆轰存在向二维结构的单头爆轰结构退化的趋势。Abstract: In order to explore the effect of the aspect ratio of rectangular tube on the propagation of the spinning detonation under the limiting detonation propagating conditions, the structure of the three-dimensional gas-phase spinning detonation wave and its propagation modes in rectangular cross-section tubes are numerically investigated based on Euler equations with a 5th-order WENO finite difference scheme and the two-step global reaction model. A linear stability theory of planar detonation wave based on the normal mode method is firstly adopted to examine the chemical reaction parameters for numerical simulations and then several cases with different aspect ratios in cross-section of rectangular tube are investigated to study the structure and propagation mode of three-dimensional gaseous spinning detonation waves. By recording motions of triple lines, flow-field distributions and high-pressure imprint of detonation wave under different sizes of tube cross-section, the effect of cross-sectional geometry on the stable propagation of gaseous detonation under the limiting detonation propagating condition is revealed. The results show that the spinning detonation propagation can be maintained within a certain range of small tube cross-section size, through the movements of horizontal and vertical triple lines and an oblique triple line that is produced by interaction between both horizontal and vertical triple lines. For a square tube with 1 of aspect ratio in cross-section, the high-pressure imprint of spinning detonation on the wall forms the helical strip pattern. With the increase of the aspect ratio of the cross-section size of the tube, the pattern of a high-pressure imprint formed by the spinning detonation on the channel wall varies from the strip structure to a dotted distribution structure, the trajectory of the oblique triple line on the wave front gradually develops from the circular motion in a single direction to a complex trajectory with varying direction. When the aspect ratio is further increased, there is a tendency for the three-dimensional spinning detonation wave to eventually degenerate into a two-dimensional single-head detonation wave structure.
-
Key words:
- spinning detonation /
- rectangular tube /
- triple line /
- two-step global reaction
-
三维爆轰波在管内传播时存在着几种复杂结构,包括单头螺旋爆轰、双头或多头胞格爆轰等。其中,螺旋爆轰属于极限条件下的爆轰,主要表现为横波数量较少,波后局部高压区沿管道轴向呈螺旋行进的传播模式。通常,螺旋爆轰都是在较小尺寸的管道或者在较低压力下测得的,相关理论、实验和数值模拟的研究较少,许多细节仍需进一步研究。
Campbell等[1]在1926年首次发现了螺旋爆轰现象的存在。Bone等[2]针对各种截面管道中的螺旋爆轰开展了实验研究,并且总结了螺旋爆轰在管道内壁上形成的高压螺旋结构的螺距与管道截面周长之间的比例关系,但研究并没有描述螺旋爆轰波的结构细节。Fay[3]使用声波耦合理论对圆管的螺旋爆轰进行了理论分析,得出在稳定螺旋爆轰模式下螺距与管直径的比值大约为3。Edwards等[4]采用理论分析与实验研究相结合的方式,进一步研究了方管内螺旋爆轰的波阵面结构。Huang 等[5]结合实验与理论的方法详细研究了圆管中螺旋爆轰的头部结构及头部周围流场物理量的变化,给出了稳定螺旋爆轰模式下流场结构与爆轰参数之间的关系。Dou等[6]基于高分辨率WENO计算格式,对方管内三维爆轰波进行了高精度数值模拟,结果表明,当管道截面尺寸较小时,爆轰会形成一种不稳定的螺旋模式进行传播。Dou等[7]分析了细方管下初始扰动在爆轰中的演化情况,证明了在一定的预混气体条件下螺旋结构是最稳定的爆轰传播模式。Tsuboi等[8-10]采用数值模拟的方法研究了圆管和方管内螺旋爆轰的传播过程,给出了爆轰波阵面结构的演化过程,并指出两种截面情况下的激波模式相似,且围绕壁面旋转的横向爆轰是维持螺旋爆轰稳定传播的重要因素。Nagao等[11]对比了方形截面管道内螺旋爆轰在不同当量比与不同压力下的混合物中的传播特性,指出螺旋爆轰在方管中的传播速度依赖于初始压力,并有周期性的变化趋势。Huang等[12]进一步开展了方管中三维爆轰结构从多头爆轰向螺旋爆轰模式转变的数值研究,利用三维可视化揭示了螺旋爆轰形成的重要特征。Sugiyama等[13-14]研究了圆管和方管内的螺旋爆轰模式,指出螺旋爆轰波阵面后面的横向爆轰波与产物区中的声波耦合是螺旋爆轰的主要传播机制。
近年来,Wang等[15-17]开发了一种新型高分辨的并行求解器,对三维长方管内的气相爆轰结构进行了细致研究,发现爆轰波的过驱度、不稳定性及初始扰动是影响螺旋爆轰形成的主要因素,并指出当管径截面尺寸足够小时,爆轰波由多头模态转变到单头模态以维持稳定传播。沈洋等[18-19]研究了矩形截面管道中的三维爆轰胞格结构,认为截面尺寸或长径比的变化最终会引起三维爆轰结构的显著差异。
由以上分析可知,目前关于圆管和方形截面管道中螺旋爆轰波传播特性的研究,对圆管和方管中波阵面的结构演化过程形成了基本共识。研究发现,不同燃料形成的螺旋爆轰在方管或圆管内会沿着单一方向旋转,且在壁面附近形成的横向爆轰波是维持旋转模式稳定传播的关键因素。虽然针对圆管和方管开展的研究较为充分,但针对矩形截面管道内螺旋爆轰波三维结构及其自持传播方式的研究报道却很少。矩形管道由于截面长宽比不同,会导致爆轰波螺旋传播时,在不同管道宽度方向的行进路线存在差别,影响了螺旋爆轰的传播和自持,因此,研究管道尺寸和几何形状对螺旋爆轰传播的影响规律十分必要。为此,本文中采用高精度数值模拟的方法,对不同尺寸的方形截面和不同长宽比下矩形截面管道内的预混气体在临界状态下的螺旋爆轰(单头模态)进行三维数值模拟,以进一步揭示不同管道截面尺寸下螺旋爆轰波的结构及其传播机制,加强对近极限条件下爆轰波自持传播特性的理解,为气相爆轰在实际系统,如爆轰推进系统中的应用提供理论依据。
1. 数值方法
1.1 控制方程
基于结构化网格有限体积框架,采用基于两步化学反应的三维非稳态反应性欧拉方程组研究矩形管内的临界爆轰现象,具体形式如下:
$$ \frac{{\partial {\boldsymbol{U}}}}{{\partial t}} + \frac{{\partial {\boldsymbol{F}}}}{{\partial x}} + \frac{{\partial {\boldsymbol{G}}}}{{\partial y}} + \frac{{\partial {\boldsymbol{H}}}}{{\partial {\textit{z}}}} = {\boldsymbol{S}} $$ (1) 式中:守恒量U,对流项F、G、H以及源项S分别为:
$$ {\boldsymbol{U}} = {\left( {\rho ,\;\rho u,\;\rho v,\;\rho w,\;E,\;\rho {Y_1},\;\rho {Y_2}} \right)^{\rm T}} $$ (2) $$ {\boldsymbol{F}} = {(\rho u,\;\rho {u^2} + p,\;\rho uv,\;\rho uw,\;(E + p)u,\;\rho u{Y_1},\;\rho u{Y_2})^{\text{T}}} $$ (3) $$ {\boldsymbol{G}} = {(\rho v,\;\rho uv,\;\rho {v^2} + p,\;\rho vw,\;(E + p)v,\;\rho v{Y_1},\;\rho v{Y_2})^{\rm T}} $$ (4) $$ {\boldsymbol{H}} = {(\rho w,\;\rho uw,\;\rho vw,\;\rho {w^2} + p,\;(E + p)w,\;\rho w{Y_1},\;\rho w{Y_2})^{\rm T}} $$ (5) $$ {\boldsymbol{S}}={\left(0,\;0,\;0,\;0,\;0,\;\rho {\dot{Y}}_{1},\;\rho \dot{{Y}_{2}}\right)}^{{\rm T}} $$ (6) $$ E = \frac{p}{{\gamma - 1}} + \frac{{\rho \left( {{u^2} + {v^2} + {w^2}} \right)}}{2} + \rho {q_2}{Y_2} $$ (7) 式中:ρ为密度,u、v和w为三个方向的速度,E为单位体积总能量,p为压力,γ为比热比,q2为反应2的放热,Y1和Y2分别为反应1和反应2的反应物质量分数,
$ \dot Y_1 $ 和$ \dot Y_2 $ 表示化学反应过程中反应物质量分数的消耗速率。化学反应采用包含诱导反应1和放热反应2的两步总包反应形式,反应速率表达式为[20]:$$ \left\{ \begin{gathered} \rho \dot Y_1= H\left( {{Y_1}} \right){A_1}{\rho ^{2.2}}\exp \left( {{E_{\mathrm{a}}}/T} \right) \\ \rho \dot Y_2 = \left[ {1 - H\left( {{Y_1}} \right)} \right]{A_2}{\rho ^2}Y_2^2 \\ H\left( {{Y_1}} \right) = \left\{ {\begin{array}{*{20}{l}} 0& {Y_1} = 1 \\ 1& 0 {\text{<}} {Y_1} {\text{<}} 1 \end{array}} \right. \end{gathered} \right. $$ (8) 式中:A1和A2分别为反应1和反应2的指前因子,H(Y1)为混合物函数,Ea为活化能。活化能Ea表征了气体的不稳定程度,其值越大说明气体越不稳定。指前因子A1和A2可控制反应区的宽度和反应时间(时空尺度)。本文中采用的两步化学反应能够合理模拟反应系统的诱导过程和反应放热过程,在爆轰波数值模拟研究中使用广泛[20-21]。化学反应使用的具体参数见表1,其中p0和T0分别为环境压力和温度,R为通用气体常数。
表 1 初始条件与化学反应参数Table 1. Initial conditions and chemical reaction parametersp0/MPa T0/K A1/[s∙(kg∙m−3)2.2] A2/[s∙(kg∙m−3)2] Ea/RT0 γ q2/RT0 0.1 300 2×107 7×105 21 1.3 24.4 为考察本文选取的化学反应参数的特性,利用基于正则方法的线性稳定性理论[22-23],考察了两步化学反应下平面爆轰波的二维线性稳定性。图1给出了表1参数条件下,扰动增长率Re(α)随扰动波数κ的变化规律。这里Re(α)表示小扰动的实部,正值代表爆轰波状态不稳定(扰动增长),负值则代表爆轰波状态稳定(扰动衰减)。图1的结果表明,当波数κ=0时,扰动是衰减的,这意味着反应系统的扰动从时间上是衰减的,而在κ=0.1~2.3的范围内扰动存在一定的增长,其最大增长率大约位于κ=1的位置。这说明在本文所使用的参数下,反应系统形成的爆轰波,其内在扰动会随时间衰减,而在横向空间上会增长,这样的结果表明,本文采用的反应参数对应于反应活性相对稳定的燃料,因而可以去除随时间的扰动所带来的干扰,从而在后继研究中关注空间扰动对于不同管道长宽比下螺旋爆轰波传播的影响。
1.2 数值计算条件
对控制方程式(1),采用五阶WENO格式对对流项进行求解,采用三阶Runge-Kutta法对时间推进和化学反应源项进行处理。为保证计算稳定性,设置库朗数CFL(Courant-Friedrichs-Lewy)为0.3。
采用激波坐标系进行管内三维爆轰波的数值模拟,计算模型如图2所示。图中管道进口(右边界)为均匀来流边界,来流压力和温度分别取环境压力和温度(见表1),来流速度为理论CJ (Chapman-Jouguet )爆轰速度,uCJ=1887 m/s(
$ {u_{{\text{CJ}}}} = {M_0}\sqrt {\gamma {\text{R}}{T_0}} $ ),来流马赫数$ {M_0} = \sqrt {1 + {K_{{\text{CJ}}}}/2 + \sqrt {{K_{{\text{CJ}}}} + K_{{\text{CJ}}}^2/2} } $ ,其中$ {K_{{\text{CJ}}}} = 2\left( {{\gamma ^2} - 1} \right){q_2}/\gamma R{T_0}$ 。管道边界条件均选用刚性反射壁面,而管道出口(左边界)采用零梯度边界条件。初始时刻在计算域中爆轰波阵面左侧区域赋予一维CJ爆速下的ZND(Zeldovich-von Neumann-Doring)结果,并在来流中给定一个高温高压的长方体区域作为初始扰动,以快速促进爆轰波阵面失稳,如图2所示。为确定用于形成螺旋爆轰的计算域的临界尺寸,首先采用Vasil’ev根据实验观察得到的临界管径判定方法[24]来确定该尺寸,即发生螺旋爆轰的最小临界尺寸约为λ/4,其中λ代表爆轰波的胞格宽度。为获得本文初始条件下预混气体的胞格尺寸,首先计算了二维条件下的几种不同管道宽度下的正爆轰波传播的胞格结构,图3给出了三种管道宽度下计算得到的爆轰波稳定传播时的胞格结构。
从图3可知,当管道宽度为2 mm时,爆轰波呈现出一个完整的胞格结构,其胞格尺寸(宽度)λ为2 mm;当管道宽度为4 mm时,沿横向的胞格数量为1.5个,当管道宽度为8 mm时,沿横向的胞格数量为3个,即胞格尺寸λ约为2.6 mm。上述结果说明,本文条件下得到的爆轰波的胞格尺寸在2~2.6 mm之间。根据文献[24]的判断准则,螺旋爆轰的临界尺寸在0.5~0.65 mm之间,因此,本文中设计了如表2所示的几种不同截面尺寸和不同截面长宽比条件下的算例,其中,算例1~3代表方形截面管道,算例1为基准算例,从算例1到算例3其尺寸逐渐变小,直至临界尺寸的范围(0.6 mm);算例1和算例4~6代表了管道截面长宽比的影响,在保证管道截面中长(1.0 mm)不变的情况下,其宽度逐渐由1 mm缩短至0.4 mm。其中,算例1采用的方管尺寸(1.0 mm)与文献[10]中的是相同的。
表 2 计算工况Table 2. Computational cases工况 管道长/mm 管道宽/mm 1 1.0 1.0 2 0.8 0.8 3 0.6 0.6 4 1.0 0.8 5 1.0 0.6 6 1.0 0.4 2. 结果与讨论
2.1 网格无关性验证
为了进一步验证网格分辨率对计算结果的影响,首先针对表2中的算例1(1.0 mm×1.0 mm的方形截面管)开展不同网格尺寸下螺旋爆轰波传播过程的计算。
图4(a)为网格尺寸为10 μm时螺旋爆轰波在展开的方管壁面上留下的倾斜的高压螺旋印记图。图4(b)为cdcʹdʹ处管壁上线段AB间的压力曲线。由图4(a)与图3(b)可见,壁面上的高压螺旋条带印记呈现出压力大小交错的周期性变化规律,将图4(a)高压条带间的宽度(图4(b)中各压力峰之间的水平距离)定义为螺距P,进而可以测量不同网格分辨率下方管螺旋爆轰高压印记条带的平均螺距
$ \bar P $ 与方管外接圆直径长度D的比值,结果如表3所示。Bone等[2]针对不同截面管道中的螺旋爆轰实验研究指出,螺旋爆轰在管道内壁上形成的高压螺旋结构的螺距与管道外接圆直径之间的比值约为3,与化学反应所用的燃料无关,因此,可以与本文采用的两步反应的结果进行对比。近年来,Wang等[15]通过数值模拟计算得到了方管下的螺距尺寸约为3.128,这与本文的结果相近。因此,表3的结果表明,在本文给定的3种网格尺寸下,螺距尺寸均未显示出明显差距,且
$ \bar P/D $ 的计算结果与其他学者的数值模拟结果基本一致。为尽可能精确刻画三维爆轰波的结构,采用10 μm的网格尺寸开展研究。由于计算采用了基于激波坐标系下的计算区域和并行计算手段,有效提高了计算效率,因此采用10 μm的网格尺寸的计算花费可以控制在合理范围内。在该网格尺寸下,本文条件的爆轰波ZND结构中的每半反应区宽度(爆轰前导激波与最大释热速率位置之间的距离)包含了16个网格,因而足以分辨螺旋爆轰波的主要结构特征。2.2 流场结构分析
为研究管道尺寸对螺旋爆轰波结构的影响,首先针对算例1研究螺旋爆轰波传播过程中的结构特征。
由于爆轰波的三维不稳定性,爆轰波阵面上存在一系列三波线(triple lines,TL)。从管道截面方向看,若三波线与水平壁面基本平行,则称为水平三波线(horizontal triple line,HTL);若三波线与垂直壁面基本平行,则称为垂直三波线(vertical triple line,VTL)。在HTL和VTL的交汇处,由于马赫反射会形成斜三波线(oblique triple line,OTL)。图5展示了1.0 mm×1.0 mm的方形截面的管道中,六个典型时刻下从来流方向观察到的爆轰波阵面形态,图中蓝色为p=1 MPa的压力等值面。此外,三波线也会将爆轰波阵面分成几个部分:入射激波阵面I(incident shock)代表其周围的三波线均向其内部运动,马赫阵面M(Mach)表示其周围的三波线均在远离该阵面。当马赫面周围的三波线正在壁面上发生反射时,则称为强马赫阵面(strong Mach,SM),其他情况下均为弱马赫面(weak Mach,WM)。当波阵面周围既有向内又有向外传播的三波线,则这些波阵面表现为过渡面(Mach to incident,MI)。对于螺旋爆轰,如果过渡面周围的三波线传播方向均与螺旋轨迹方向相同时称为弱过渡面(weak Mach to incident,WMI),否则称为强过渡面(strong Mach to incident,SMI)。上述阵面类型涵盖了所有可能的三维爆轰波阵面情况,且均在图5中进行了标记。
当t=8.70 μs时(图5(a)),向上运动的HTL被向左运动的VTL分为HTL1和HTL2两部分,其中HTL1由于靠近上壁面而不可见;同样向左运动的VTL也被HTL分为VTL1和VTL2两部分,其中VTL1较短,未在图中标出。HTL和VTL交汇形成的OTL正靠近上壁面。此时刻波阵面大致分为三个区域,在左上区域VTL1与HTL2均在向其内部运动,所以此区域为入射激波阵面I;在左下方区域,HTL2远离该区域,VTL2指向该区域,但二者运动方向均与螺旋轨迹方向相反,因此该区域为SMI;而在右侧区域,VTL2远离该区域,HTL1也在远离该区域且正在与壁面发生碰撞,因此该区域为强马赫阵面SM。在t=8.94 μs时刻,爆轰波阵面被经上壁面反射而向下传播的HTL1、HTL2和持续向左传播的VTL1、VTL2分割为四个区域:根据前面的命名规则,分别为WM、WMI、SMI和I(见图5(b))。此外,HTL和VTL由于马赫反射在交接处形成的OTL清晰可见,此时OTL按逆时针方向经上壁面反射后向左下方传播。到t=9.18 μs时,三波线未发生碰撞导致运动方向改变,因此波阵面结构与8.94 μs时相似,HTL和VTL的运动方向也与8.94 μs时相同,见图5(c)。当t=9.40 μs时,VTL与cd壁面碰撞,导致HTL2消失和VTL1和VTL2向右反射,此时波阵面大致被分为SMI和I两个区域(图5(d))。到t=9.62 μs时,VTL碰撞后反射向右传播,OTL反射后向右下方传播,其中VTL2较短,未在图中标出,因为从上一时刻到9.62 μs时存在三波线与壁面碰撞后改变运动方向,所以此时波阵的结构也因三波线运动的轨迹发生了改变,在左上方区域,HTL2和VTL1均在远离该区域,且无三波线与壁面碰撞,所以该区域为WM;在左下方区域,HTL2向其内部运动,而VTL2在远离该区域,且二者运动方向与螺旋轨迹方向相同,所以该区域为WMI,而其他两个阵面结构(SMI和I)与前一个时刻的相同,见图5(e)。t=9.86 μs时,三波线未发生碰撞导致运动方向改变,所以 HTL和VTL的运动方向也与9.62 μs时相同,但由于此时壁面处正在发生碰撞,所以波阵面结构与上一时刻有所差异,此时左侧区域为SM,而右侧区域仍为SMI,见图5(f)。由于图4给出了约1/2周期的三波线运动轨迹图,观察其OTL的轨迹可以发现,其在方形管道中的完整轨迹应为逆时针旋转的菱形结构,这一结构与Wang等[17]的计算结果是相同的。总体来说,三波线(尤其是斜三波线)在壁面上的碰撞是维持螺旋爆轰波稳定传播的关键,爆轰波阵面通过水平和垂直三波线在壁面的往复碰撞以及斜三波线的转向碰撞,维持稳定的传播。
图6中的斜三波线(OTL)的运动轨迹与图3中的高压印记密切相关。为说明这种相关性,图6进一步给出了图5中6个时刻下螺旋爆轰波的三维结构。图中,绿色等值面代表波阵面,蓝色和黄色等值面代表反应面,其中蓝色一侧面向未燃气体,而黄色一侧面向已燃气体。红色的等值面代表此瞬时流场中的高压区,其值约为30 MPa,这一区域实际上是由OTL后方的横波结构与壁面碰撞诱导的局部爆轰。由图6可见,在8.70 μs时,OTL接近上壁面,而此时高压区也与上壁面接触,故在图4(a)所示的壁面上呈现出明显的高压图案。在8.94 μs时,OTL与baaʹbʹ壁面碰撞后旋转方向改变,因此OTL后方的横波也开始向addʹaʹ壁面运动,此时因高压区域开始远离baaʹbʹ壁面,因此在baaʹbʹ留下的高压印记的面积逐渐减小,其压力大小也有明显下降。到9.18 μs时,由于期间不存在三波线与壁面的碰撞,所以斜三波线后方的高压区持续减弱。9.40 μs时,OTL与addʹaʹ壁面碰撞,可以看到此时OTL后方高压区面积增大且OTL运动方向改变。9.62 μs时,OTL按照固定的螺旋方向往addʹaʹ壁面运动,因为在此阶段内不存在OTL与壁面的碰撞,所以此时OTL后方的高压区面积开始减小。9.86 μs时,OTL与addʹaʹ壁面碰撞,高压区面积增大。由图6的结果可以看出,波阵面上斜三波点的位置对应于壁面上局部爆轰区的位置,因此,斜三波线的运动轨迹实际上代表了壁面上高压印记的轨迹。
为进一步澄清导致壁面上高压印记的局部高压区的结构,图7给出了算例1在t=8.94 μs时刻管道壁面上的温度和压力分布图。由图7(a)可以看出,在add′a′壁面上爆轰波阵面呈现出复杂的马赫反射结构。图7(b)表明,三波线VTL1在壁面add′a′上留下的三波点A连接了分别由马赫面WM和WMI在壁面上留下的波线以及延展至产物区的横波。图7(c)为该结构的相应的示意图,图中AB、BD为反射激波,BC为横激波,CE为入射激波后方,前导激波由入射激波和马赫杆组成,虚线表示反应面,带方向的虚线表示流体运动方向。由于BC段的横波强度较大进而引发了横向爆轰,因此BC段也称为横向爆轰波。正是高强度的横向爆轰成为留在壁面上的高压印记。因此,图7的结果表明,在该时刻局部爆轰区为图7(c)所示的“强横波结构”,这一结构在文献[5, 25-26]中得到了印证。
总之,图5中的斜三波线会在四个壁面的中点附近碰撞,其对应于壁面上横向爆轰波留下的高压印记。因此,通过综合分析图5的三波线运动轨迹与图4高压印记条带在展开壁面上的形态就能揭示管道中三维螺旋爆轰波的运动轨迹。为更加直观地理解不同管道截面尺寸下螺旋爆轰波的传播方式,后续将采用类似于图4~5的计算结果来分析矩形管内三维爆轰波的传播方式及螺旋结构。
2.3 管道截面尺寸和长宽比的影响
图8首先给出了不同管道宽度下的方管中,螺旋爆轰在展开壁面上的高压印记图。可以发现,在1.0 mm×1.0 mm到0.6 mm×0.6 mm的三个不同截面尺寸的方管中,展开的壁面上均能观察到斜条带形的高压印记,这表明螺旋爆轰可沿管道轴向方向自持传播,但是不同尺寸下螺旋爆轰波的局部爆轰(高压)区沿管壁的旋转方向不同。图7的结果还表明,随着管道截面尺寸的减小,螺旋爆轰在壁面上留下的最大压力也在减小,其中,1.0 mm×1.0 mm的方管壁面的压力峰值为52 MPa,而0.6 mm×0.6 mm的方管壁面压力峰值则减小到16 MPa,这主要是由于管道尺寸减小后,没有足够的空间支持气体在短时间内完全燃烧放热,因而导致压强峰值下降。
为进一步研究矩形管道中截面长宽比的影响,以算例1的截面尺寸为基准,通过改变其截面的长宽比来研究在不同长宽比下矩形管中临界爆轰波的传播特性。
图9给出了算例1、4~6在展开壁面上的高压印记图的变化。可以发现,随着长宽比的增加(算例4~5),壁面上显示的螺旋结构开始变得不规则,原来较为明显的倾斜的条带结构(图9(a))逐渐变成存在局部高压区的点状结构,其条带走向也变得不再明显(图9(b)~(c)),当长宽比较大时(图9(d)),高压印记在展开壁面上形成了规则的六边形,其中在窄壁面上存在高压区并构成了六边形的短边,而在宽壁面上形成了较细的高压印记并构成了六边形的长边。
2.2节的讨论表明,矩形管道展开壁面上的高压印记与波阵面上斜三波线(OTL)的运动轨迹密切相关,为此进一步考察了不同截面长宽比下爆轰波阵面上OTL的运动情况。其中,长宽1.0 mm×1.0 mm的算例1的三波线在半个运动周期内的轨迹变化已在图5中给出,而图10~12分别给出了长宽为1.0 mm×0.8 mm、1.0 mm×0.6 mm和1.0 mm×0.4 mm条件下,螺旋爆轰在一个周期内的三波线运动过程。
图10的结果表明,在爆轰波阵面上,OTL在一个周期内形成了从左上到右下的对角线方向行进的重复性运动轨迹,这与方形管中得到的菱形轨迹(图4)有所不同。此外,图10中,从7.19 μs到7.43 μs ,OTL在管道左上角发生连续碰撞,因而会在壁面留下较大区域的高压印记,同时,从8.14 μs到8.39 μs,在管道右下角也会经历一次类似的连续碰撞并留下较大的高压印记。因此,螺旋爆轰在一个周期内会在壁面上留下两处面积较大的高压区域,这与图8(b)中高压印记的点状分布是一致的。
图11给出了算例5的长宽(1.0 mm×0.6 mm)条件下,螺旋爆轰在一个完整周期内的三波线运动过程。图中结果表明,在壁面宽度更窄的情况下,OTL在一个周期内运动形成了“α”形的运动轨迹。这种运动轨迹与算例4中在对角线附近的碰撞不同,算例5的OTL在一个周期内的运动是顺时针与逆时针相结合的周期性运动,在每个方向上运动时都会在壁面上产生3次碰撞,因而会在一个周期内留下对称的高压印记,主要表现为因与ab壁面接触的位置不同而分别形成靠近cd或ab管壁的高压区,反映在三维高压印记中则如图9(c)所示。
图12给出了算例6的长宽(1.0 mm×0.4 mm)条件下,螺旋爆轰在一个完整周期内的三波线运动过程。图中结果表明,在进一步狭窄的矩形管内,OTL在一个周期内形成了“∞”形的运动轨迹,该轨迹类似于图5中两个相差90°的菱形相叠加,因此相较于算例5具有更好的对称结构。由于OTL在一个周期内在四个顶角处各有一次连续碰撞,因而相较于算例4的高压印记,算例6在一个周期内存在四处明显的高压区,由于此时管道存在窄壁面,靠近窄壁面处的两次连续碰撞较为接近,因此可以视为一次范围较大的连续碰撞过程,因此在一个周期内在窄壁面处存在两次高压印记,这与图9(d)显示的结果是一致的。
通过对不同截面长宽比的数值模拟结果的分析,可以发现,随着长宽比的增加,波阵面上斜三波线(OTL)的轨迹会逐渐改变,对方形截面管道(算例1),OTL的轨迹是依次与方管壁面的中心处碰撞并沿单一方向做周期性运动;当长宽比发生变化时(算例4),OTL在壁面的碰撞位置开始偏离壁面中心而向壁面角端移动,这就使得一个周期内的4次碰撞中每两次的位置离得更近,高压区发生部分重叠,因此逐渐形成了图9(b)中两处面积较大的高压区;当长宽比进一步增加时(算例5~6),OTL在壁面的角端碰撞时,由于过于靠近壁面角部,因此会导致碰撞反射后出现运动轨迹的转向。从图11~12中可以看出,在算例5~6的长宽比下,OTL在运动过程中,会越过管壁右侧的cʹdʹ壁面直接与上下壁面发生碰撞,进而导致OTL的运动方向发生改变,产生更为复杂的“α”形和“∞”形运动轨迹。
由于OTL是HTL和VTL经马赫反射得到的,因此OTL的轨迹变化实际上是由HTL和VTL的运动所决定的。为此,进一步计算了不同长宽比下HTL和VTL的运动速度,图13给出了不同算例在一个OTL运动周期内HTL和VTL的平均运动速度,为便于比较,图中还给出了CJ爆速DCJ和爆轰波诱导区声速ci的值。可以看出,在方管条件(算例1)下,HTL 和VTL的速度都与诱导区声速ci十分接近,但远小于DCJ,表明此种情况下HTL/VTL运动与声波横向振动之间存在良好的耦合关系,符合Fay等提出的声波耦合理论[3]。对于矩形截面管道(算例 4~6),HTL/VTL运动速度之间存在差别,其中沿着管道较窄方向运动的HTL的速度比沿着管道较宽方向运动的VTL的速度要低,这是螺旋爆轰稳定传播时V/HTL根据管道截面尺寸自行调整的结果。由于管道较宽一侧的尺寸决定了螺旋爆轰的稳定传播,因此在矩形管内沿着管道较宽方向运动的VTL会起到主导作用。图中结果表明,尽管VTL的速度略大于声速,但差别并不大,表明矩形截面管道的横波的传播也基本符合声波耦合理论。根据上面的分析可知,当VTL/HTL以接近声波的速度传播时,由于长宽比的不同,会导致VTL和HTL存在不同的相位差,因此两者相互作用产生的OTL的轨迹就会随着长宽比的不同而不同。
图14给出了根据本文计算的不同管道截面长宽比下的OTL的运行轨迹所绘制的示意图,可以看出,其轨迹逐渐从方形管(1.0 mm×1.0 mm)中的菱形轨迹转变为窄矩形管(1.0 mm×0.4 mm)中“∞”形轨迹,这意味着随着管截面长宽比的变化,爆轰波阵面逐渐从三维螺旋结构向二维单头结构过渡。因此本文认为,只要矩形管道截面上较长一侧的长度大于产生临界爆轰的极限尺寸,则爆轰波仍会稳定传播;但当截面长宽均小于产生临界爆轰的极限尺寸时,爆轰波将会熄爆。
3. 结 论
利用两步诱导-放热反应对矩形截面管道内的螺旋爆轰波的结构和传播模式进行了高精度三维数值研究,着重考察了不同管道截面尺寸和形状(长宽比)对螺旋爆轰波结构的影响。得到以下主要结论。
(1)在本文的研究条件下,螺旋爆轰波阵面上的斜三波线在阵面后下游方向存在强横波结构,其中的横向爆轰波跟随斜三波线在管内做沿轴向旋转的周期运动,并在壁面上留下高压印记,因此横向爆轰是维持爆轰波稳定传播的重要因素。
(2)对方形截面管道,螺旋爆轰均能在壁面形成稳定的高压条带印记,随着方形截面尺寸的减小,壁面上高压印记的最大压力也逐渐减小。
(3)当管道截面长宽比较小时,即接近正方形结构时,其斜三波线的轨迹呈菱形变化;随着长宽比的增加,即接近长方形时,斜三波点轨迹呈现出“α”形的变化规律;当长宽比进一步增加时,三波点轨迹呈现出“∞”形的变化规律。即长宽比的增加会使斜三波线轨迹逐渐趋于二维往复运动的模式。
-
表 1 初始条件与化学反应参数
Table 1. Initial conditions and chemical reaction parameters
p0/MPa T0/K A1/[s∙(kg∙m−3)2.2] A2/[s∙(kg∙m−3)2] Ea/RT0 γ q2/RT0 0.1 300 2×107 7×105 21 1.3 24.4 表 2 计算工况
Table 2. Computational cases
工况 管道长/mm 管道宽/mm 1 1.0 1.0 2 0.8 0.8 3 0.6 0.6 4 1.0 0.8 5 1.0 0.6 6 1.0 0.4 -
[1] CAMPBELL C, WOODHEAD D W. CCCCI.—the ignition of gases by an explosion-wave. part I. carbon monoxide and hydrogen mixtures [J]. Journal of the Chemical Society (Resumed), 1926, 129: 3010–3021. DOI: 10.1039/JR9262903010. [2] BONE W A, FRASER R P, WHEELER W H. A photographic investigation of flame movements in gaseous explosions Ⅶ—the phenomenon of spin in detonation [J]. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences, 1935, 235(746): 29–68. DOI: 10.1038/136974a0. [3] FAY J A. A mechanical theory of spinning detonation [J]. The Journal of Chemical Physics, 1952, 20(6): 942–950. DOI: 10.1063/1.1700655. [4] EDWARDS D H, PARRY D J, JONES A T. The structure of the wave front in spinning detonation [J]. Journal of Fluid Mechanics, 1966, 26(2): 321–336. DOI: 10.1017/S0022112066001265. [5] HUANG Z W, LEFEBVRE M H, VAN TIGGELEN P J. Experiments on spinning detonations with detailed analysis of the shock structure [J]. Shock Waves, 2000, 10(2): 119–125. DOI: 10.1007/s001930050185. [6] DOU H S, TSAI H M, KHOO B C, et al. Simulations of detonation wave propagation in rectangular ducts using a three-dimensional WENO scheme [J]. Combustion and Flame, 2008, 154(4): 644–659. DOI: 10.1016/j.combustflame.2008.06.013. [7] DOU H S, KHOO B C. Effect of initial disturbance on the detonation front structure of a narrow duct [J]. Shock Waves, 2010, 20(2): 163–173. DOI: 10.1007/s00193-009-0240-8. [8] TSUBOI N, ASAHARA M, ETO K, et al. Numerical simulation of spinning detonation in square tube [J]. Shock Waves, 2008, 18(4): 329–344. DOI: 10.1007/s00193-008-0153-y. [9] TSUBOI N, ETO K, HAYASHI A K. Detailed structure of spinning detonation in a circular tube [J]. Combustion and Flame, 2007, 149(1/2): 144–161. DOI: 10.1016/j.combustflame.2006.12.004. [10] TSUBOI N, HAYASHI A K. Numerical study on spinning detonations [J]. Proceedings of the Combustion Institute, 2007, 31(2): 2389–2396. DOI: 10.1016/j.proci.2006.07.262. [11] NAGAO T, ASAHARA M, HAYASHI A K, et al. Numerical analysis of spinning detonation dependency on initial pressure using AUSMDV scheme [C] // 51st AIAA Aerospace Sciences Meeting Including the New Horizons Forum and Aerospace Exposition. Grapevine: AIAA, 2013. DOI: 10.2514/6.2013-1177. [12] HUANG Y, JI H, LIEN F, et al. Numerical study of three-dimensional detonation structure transformations in a narrow square tube: from rectangular and diagonal modes into spinning modes [J]. Shock Waves, 2014, 24(4): 375–392. DOI: 10.1007/s00193-014-0499-2. [13] SUGIYAMA Y, MATSUO A. Numerical study of acoustic coupling in spinning detonation propagating in a circular tube [J]. Combustion and Flame, 2013, 160(11): 2457–2470. DOI: 10.1016/j.combustflame.2013.06.003. [14] SUGIYAMA Y, MATSUO A. Numerical analysis on acoustic coupling of spinning detonation in a square tube [J]. Journal of Thermal Science and Technology, 2016, 11(1): JTST0010. DOI: 10.1299/jtst.2016jtst0010. [15] WANG C, SHU C W, HAN W H, et al. High resolution WENO simulation of 3D detonation waves [J]. Combustion and Flame, 2013, 160(2): 447–462. DOI: 10.1016/j.combustflame.2012.10.002. [16] 王成, 韩文虎, 宁建国. 三维气相爆轰动态并行计算程序设计与开发 [J]. 计算力学学报, 2012, 29(6): 948–953. DOI: 10.7511/jslx201206022.WANG C, HAN W H, NING J G. Design and development of dynamic parallel computing code for three-dimensional gaseous detonation [J]. Chinese Journal of Computational Mechanics, 2012, 29(6): 948–953. DOI: 10.7511/jslx201206022. [17] WANG C, LI P, GAO Z, et al. Three-dimensional detonation simulations with the mapped WENO-Z finite difference scheme [J]. Computers & Fluids, 2016, 139: 105–111. DOI: 10.1016/j.compfluid.2016.04.016. [18] 沈洋, 申华, 刘凯欣. 矩形通道中三维气相爆轰的三波线结构分析[C]//北京力学会第21届学术年会暨北京振动工程学会第22届学术年会论文集. 北京: 北京力学会, 北京振动工程学会, 2015. DOI: 10.11883/1001-1455(2016)05-0577-06. [19] SHEN Y, SHEN H, LIU K X, et al. Three-dimensional detonation cellular structures in rectangular ducts using an improved CESE scheme [J]. Chinese Physics B, 2016, 25(11): 114702. DOI: 10.1088/1674-1056/25/11/114702. [20] XIAO Q, SOW A, MAXWELL B M, et al. Effect of boundary layer losses on 2D detonation cellular structures [J]. Proceedings of the Combustion Institute, 2021, 38(3): 3641–3649. DOI: 10.1016/j.proci.2020.07.068. [21] SHORT M, SHARPE G J. Pulsating instability of detonations with a two-step chain-branching reaction model: theory and numerics [J]. Combustion Theory and Modelling, 2003, 7(2): 401–416. DOI: 10.1088/1364-7830/7/2/311. [22] SHORT M. Multidimensional linear stability of a detonation wave at high activation energy [J]. SIAM Journal on Applied Mathematics, 1997, 57(2): 307–326. DOI: 10.1137/s0036139995288101. [23] SHORT M, DOLD J W. Linear stability of a detonation wave with a model three-step chain-branching reaction [J]. Mathematical and Computer Modelling, 1996, 24(8): 115–123. DOI: 10.1016/0895-7177(96)00144-6. [24] VASIL’EV A A. Near-limiting regimes of gaseous detonation [J]. Combustion, Explosion and Shock Waves, 1987, 23(3): 358–364. DOI: 10.1007/BF00748799. [25] VOITSEKHOVSKII B V, MITROFANOV V V, TOPCHIYAN M E. Structure of the detonation front in gases(survey) [J]. Combustion, Explosion and Shock Waves, 1969, 5(3): 267–273. DOI: 10.1007/BF00748606. [26] FICKETT W, DAVIS W C. Detonation [M]. Berkeley: University of California Press, 1979. -