多次激波诱导正弦扰动预混火焰界面失稳的数值研究

陈霄; 董刚; 蒋华; 吴锦涛

doi:10.11883/1001-1455(2017)02-0229-08

多次激波诱导正弦扰动预混火焰界面失稳的数值研究

doi: 10.11883/1001-1455(2017)02-0229-08

南京理工大学瞬态物理重点实验室，江苏南京 210094

基金项目:

国家自然科学基金项目 11372140

详细信息

作者简介:
陈霄(1990-)，女，博士研究生

通讯作者:
董刚，dgvehicle@yahoo.com

中图分类号: O381
计量
- 文章访问数: 4376
- HTML全文浏览量: 1443
- PDF下载量: 332
- 被引次数: 0
出版历程
- 收稿日期: 2015-07-20
- 修回日期: 2015-11-20
- 刊出日期: 2017-03-25

Numerical studies of sinusoidally premixed flame interface instability induced by multiple shock waves

Key Laboratory of Transient Physics, Nanjing University of Science and Technology, Nanjing 210094, Jiangsu, China

摘要

摘要: 激波诱导火焰失稳是实际中常见的现象，为深入研究火焰失稳特性，采用三维单步化学反应的Navier-Stokes方程和9阶weighted essentially non-oscillatory (WENO)的高精度格式，对不同马赫数的入射激波及其反射激波多次诱导正弦型预混火焰界面失稳的现象进行了三维数值模拟，并对计算结果的可靠性进行了验证。研究结果显示，在激波的多次作用下，火焰界面的演变过程主要受Richtmyer-Meshkov (RM)不稳定特性和化学反应特性的双重影响，且随着入射激波强度的增强，上述2种特性均得到进一步强化。为构造体现反应性RM不稳定特性的参数，根据火焰界面混合区平均涡量和化学反应速率，提出了表征界面受不稳定性和化学反应影响的量纲一参数。通过分析发现，在同一入射激波强度下，该参数的对数形式随入射激波和反射激波的多次作用呈基本相同的线性增长趋势；对不同马赫数的入射激波，该参数对数形式的线性增长率也基本一致。这样的变化表明该量纲一参数能够反映反应性RM不稳定过程中火焰界面发展的内在规律。
- 激波 /
- 火焰失稳 /
- WENO /
- RM不稳定 /
- 化学反应
Abstract: In this work, to further study the features of the shock-wave induced flame instability, the two-dimensional Navier-Stokes (NS) equations with the single-step chemical reaction and the high resolution 9th-order weighted essentially non-oscillatory (WENO) scheme were adopted to simulate the instability of the sinusoidally premixed flame induced by incident shock waves with different Mach numbers and its reshock waves. The computational results were validated by the experimental results in the related literature. The computational results show that the evolutions of the flame are mainly influenced by both the Richtmyer-Meshkov (RM) instability and the chemical reaction. With the growth of the incident shock wave intensity, the interface instability and the chemical reaction are enhanced. To construct the parameter that can characterize the RM instability in the reactive flow, a dimensionless parameter 8338A131 that describes the interface RM instability and the chemical reactivity was proposed based on the average vorticity and the chemical reaction rate calculated in the mixing zone of the flame interface. The analysis of the parameter shows that, with the similar intensity of the incident shock wave, the logarithmic form of the parameter exhibits basically the same linear growth when an incident shock wave with a given Mach number and its reshocks successively pass through the flame interface. The linear growth rate of the logarithmic form of the parameter is also basically the same for different Mach numbers of the incident wave. Such variations of η suggest that the dimensionless parameter proposed in the present study can well characterize the intrinsic features of the flame interface development in the reactive RM instability process.
- shock wave /
- flame instability /
- WENO /
- Richtmyer-Meshkov instability /
- chemical reaction

HTML全文

高台阶抛掷爆破技术主要应用在大型露天矿生产中, 能大幅度降低剥离成本, 提高露天矿的生产水平。对高台阶抛掷爆破技术有了大量研究, 在爆堆形态预测方面应用较成熟的理论为弹道理论模型和WEIBULL模型^[1-2], 然而高台阶抛掷爆破技术成功与否的关键在于岩石抛掷过程中的最大起始速度。由于爆破过程的复杂性, 对爆破现象进行详细的观测比较困难, 利用数值模拟技术, 可以再现爆破过程。软件AUTODYN在求解碰撞冲击、爆炸、冲压成型等瞬态问题方面有独到的优势, 在相关领域已经得到了广泛的应用^[3-4]。本文中, 选用软件AUTODYN, 算法上采用光滑粒子流体动力学方法(smoothed particle hydrodynamics, SPH)和有限元法(finite element method, FEM), 对高台阶抛掷爆破进行模拟, 同时通过现场高速摄影设备及其技术方法, 分析爆炸作用下的台阶坡面岩石速度场分布特征, 并对数值计算过程中所采用的RHT(Riedel-Hiermaier-Thoma)材料本构模型的实用性及其参数的可靠性进行验证。

1. 计算模型及材料参数

1.1 模型构建

以哈尔乌素露天矿高台阶抛掷爆破工程为背景, 构建二分之一等比例的三维高台阶抛掷爆破模型。模型中高台阶坡面及地表面部分定义为自由边界, 其他与岩层接触部位定义为无反射边界; 模型的主体部分采用粒子单元进行网格划分, 考虑粒子单元无法定义边界条件, 在需要定义无反射边界条件的区域构建有限元单元, 采用SPH-FEM耦合的方法进行边界约束, 再对有限元边界定义无反射边界条件。台阶高度35 m, 坡面延伸方向长20 m, 坡面倾角65°, 坡面背侧方向上坡顶及坡底长度分别为17.00和33.29 m; 炮孔布置方面, 模型中构建一个半炮孔, 炮孔采用与台阶坡面平行的倾斜炮孔, 直径310 mm, 孔距11 m, 炮孔至坡面的水平距离为9 m; 炮孔填塞长度7 m; 起爆方式采用孔底起爆。在与对称面交叉部位的台阶坡面上, 由上向下依次在距离台阶顶端7、22、30 m处选定3个高斯记录点。有限元及其他部分布置情况、采用软件AUTODYN构建的二分之一三维模型见图 1。

图 1 高台阶抛掷爆破的三维模型

Figure 1. Three dimensional model of high bench cast blasting

下载: 全尺寸图片幻灯片

1.2 本构模型及参数

RHT材料本构模型综合考虑了在材料破坏过程中所具有的应变硬化、应变率敏感性、压力依赖性和压缩损伤软化等特性, 考虑了压静水和拉静水区应变速率敏感性的差异性, 引入了偏应力张量第三不变量对破坏面形状的影响, 采用了不同的动力放大系数。另外, RHT模型还引入了拉伸损伤, 拉伸和压缩损伤均取决于等效塑性应变, 与材料塑性体积变化无关。RHT模型可应用到弹体贯穿、侵彻、爆炸等数值模拟研究中。该模型引入了最大失效面、弹性极限面和残余失效面3个控制破坏面^[5-7]。失效方程为:

$\sigma_{\mathrm{eq}}^{*}(p, \theta, {\dot \varepsilon })=Y_{\mathrm{TXC}}^{*}(p) R_{3}(\theta) F_{\mathrm{rate}}({\dot \varepsilon })$

(1)

式中: 表示压缩子午线上的等效应力强度, 为应变率强化因子, 分别为:

$Y_{\mathrm{TXC}}^{*}(p)=A\left(p^{*}-p_{\mathrm{spall}}^{*} F_{\mathrm{rate}}({\dot \varepsilon })\right)^{N}$

(2)

${F_{{\rm{rate }}}}(\dot \varepsilon ) = \left\{ {\begin{array}{*{20}{l}} {{{\left( {\frac{{\dot \varepsilon }}{{{{\dot \varepsilon }_0}}}} \right)}^\alpha }}&{p > \frac{{{f_c}}}{3},\quad {{\dot \varepsilon }_0} = 3 \times {{10}^{ - 5}}{{\rm{s}}^{ - 1}}}\\ {{{\left( {\frac{{\dot \varepsilon }}{{{{\dot \varepsilon }_0}}}} \right)}^\alpha }}&{p > \frac{{{f_{\rm{t}}}}}{3},\quad {{\dot \varepsilon }_0} = 3 \times {{10}^{ - 6}}{{\rm{s}}^{ - 1}}} \end{array}} \right.$

(3)

R₃(θ)是相似角θ以及偏平面上失效曲线在拉、压子午线处的偏应力的比Q₂的函数, 可以描述失效强度与偏应力张量第二不变量J₂、第三不变量J₃及净水压力有关的实验结果:

$R_{3}(\theta)=\frac{2\left(1-Q_{2}^{2}\right) \cos \theta+\left(2 Q_{2}-1\right)\left(4\left(1-Q_{2}^{2}\right) \cos ^{2} \theta+5 Q_{2}^{2}-4 Q_{2}\right)^{\frac{1}{2}}}{4\left(1-Q_{2}^{2}\right) \cos ^{2} \theta+\left(1-Q_{2}^{2}\right)^{2}}$

(4)

$\begin{array}{r} \theta=\frac{1}{3} \arccos \left(\frac{3 \sqrt{3} J_{3}}{2 J_{2}^{\frac{1}{2}}}\right) \quad 0 \leqslant \theta \leqslant \frac{\pi}{3} \\ Q_{2}=\frac{r_{1}}{r_{\mathrm{c}}}=Q_{0}+B_{Q} P^{*} \quad 0.51 \leqslant Q_{2} \leqslant 1.00 \end{array}$

式中:量纲一层裂强度为层裂强度; r_t是拉子午线处的偏应力, r_c是压子午线处的偏应力, B_Q是压力影响系数, A、N、α、δ和Q₀是材料常数。

计算中采用的岩石RHT本构模型失效强度参数分别为: 其他参数采用文献[5]中的低强度混凝土参数。

炸药采用铵油炸药, 结合JWL状态方程, 模拟炸药爆炸过程中的压力与体积的关系:

$p=A\left(1-\frac{\omega \eta}{R_{1}}\right) \mathrm{e}^{\frac{R_{1}}{\eta}}+B\left(1-\frac{\omega \eta}{R_{2}}\right) \mathrm{e}^{\frac{R_{1}}{\eta}}+\omega \rho e$

(5)

式中:A、B、R₁、R₂、ω为实验拟合参数; e为比内能; η=ρ/ρ₀, ρ₀为初始密度。采用的炸药材料参数见文献[8]。

2. 数值计算结果

运用软件AUTODYN, 按上述模型尺寸及本构参数构建模型并进行计算, 获得了高台阶坡面体中岩石在抛掷过程中的形态和速度变化规律。在模型的对称面部位进行切片处理, 提取坡面上设置的3个高斯点的速度时程曲线, 如图 2所示; 提取6个时间点的岩石抛掷形态及速度云图, 如图 3所示。

图 2 高斯点速度

Figure 2. The velocity curves of Gauss points

下载: 全尺寸图片幻灯片

图 3 起爆后不同时间岩石运动形态

Figure 3. Variety of damage and deformation of rock when different sequence

下载: 全尺寸图片幻灯片

对计算结果进行分析研究后得到如下结论:

(1) 起爆点设置在孔底部位, 起爆方式采用反向起爆, 整个炮孔内的炸药由起爆点开始起爆至爆轰波传递完成耗时约10 ms, 和装药长度与炸药爆轰速度的比保持一致。

(2) 起爆后, 台阶坡面岩石由静止到产生明显位移时历时约50 ms, 大约在100 ms时抛掷速度达到最大, 最大初速度在22~30 m/s之间, 3个高斯点的最大初速度为22.5、29.6、24.1 m/s。

(3) 台阶坡面岩石运动初期, 各个测点的运动速度差别不大, 随着时间的推移, 在速度增大过程中, 台阶中部岩石移动速度逐渐高出上下两端, 致使岩体整体运动状态呈现中间鼓起、两端略低的形态。

(4) 台阶坡面岩石最大水平抛掷距离约60 m, 整个抛掷过程历时约3.95 s。

3. 高速摄影实验

3.1 实验

哈尔乌素露天矿爆破爆区范围80 m×670 m, 高台阶抛掷爆破现场实验所用的高速摄影设备为Mega Speed Corp公司MS系列CCD和CMOS高速摄影机。为保证人员和仪器的安全, 高速摄影机架设点布置在东观礼台附近(见图 4(a)), 采用手持式GPS获取测点和摄影设备布设点的三维坐标, 确定观测距离为1 073 m, 摄影主光轴倾角和斜角分别为7.01°和17.37°。地面标志点作为高速摄影实验的主要参照物, 为保证在爆破过程中不产生大地位移, 在高速摄影设备一侧距爆区边界50 m处, 在地面树立旗子, 旗子大小120 cm×180 cm(见图 4(b))。

图 4 摄影机架设及地面标志点

Figure 4. Erect high-speed cameras and ground mark

下载: 全尺寸图片幻灯片

3.2 结果

通过现场高速摄影实验, 采用图像分析软件AVI Player, 通过Grid Setings功能在图像上设置网格, 图像分析中采用框标坐标系, 以左右两边的中点连线为x轴, 上下两边的中点连线为y轴, 以两联系交点为原点。读取观测点的坐标, 每个观测点及其坐标测试结果见表 1。

表 1 观测点坐标

Table 1. Imagecoordinates of observation points

t/ms	l=7 m		l=22 m		l=30 m
t/ms	x/m	y/m	x/m	y/m	x/m	y/m
42	37.65	17.57	37.82	13.40	39.36	10.92
95	37.95	17.85	38.15	13.87	39.65	11.25
203	38.25	18.09	38.37	14.05	39.88	11.46
310	38.60	18.37	38.80	14.20	40.01	11.67
428	38.92	18.60	39.12	14.50	40.50	11.91
541	39.21	18.75	39.38	14.73	40.75	12.20
655	39.39	19.05	39.63	15.05	41.17	12.35
778	39.90	19.37	39.95	15.31	41.60	12.55
889	40.55	19.65	40.38	15.43	41.79	12.75
1 009	40.90	19.87	40.80	15.62	42.25	12.93

下载: 导出CSV

| 显示表格

运用高速摄影分析原理及方法^[9], 对高速摄影录像分析后得出如下结论:

(1) 根据高速摄影过程, 将台阶坡面岩石抛掷运动过程分为静止、加速运动和自由落体运动3个阶段。

(2) 根据测点、地面标志点及拍摄点的空间坐标及记录时间的关系, 获得了台阶坡面上3个记录点的位移及速度变化规律, 如图 5所示。台阶坡面岩石抛掷过程中, 由静止到产生明显移动所需时间在33~42 ms之间, 达到最大抛掷速度的时间在93~105 ms之间, 最大初速度在18~28 m/s之间。

图 5 高速摄影速度曲线

Figure 5. The velocity curves from high-speed photography

下载: 全尺寸图片幻灯片

(3) 岩石抛掷速度达到最大值后, 呈现带有起伏波动的下降趋势, 这主要由岩体破碎后岩块间的相互碰撞作用引起的。

4. 结论

(1) 构建了高台阶抛掷爆破的三维数值计算模型, 运用软件AUTODYN对岩石抛掷过程进行了数值模拟, 模拟结果与高速摄影实验结果非常接近, 成功获得了高台阶抛掷爆破岩石抛掷速度变化规律。

(2) 采用高速摄影技术, 对哈尔乌素露天矿高台阶抛掷爆破工程进行现场高速摄影实验, 得到了爆炸瞬间完整的岩石鼓包过程, 将台阶坡面岩石抛掷运动过程分为静止、加速运动和自由落体运动3个阶段。岩石抛掷过程中, 台阶中部岩石移动速度逐渐高出上下两端, 使岩体整体运动状态呈现中间鼓起、两端略低的形态。

(3) 通过对标记点的研究分析, 获得了岩石鼓包过程中的最大抛掷速度, 与数值计算结果一致, 验证了数值计算结果的正确性。说明RHT本构模型在高台阶抛掷爆破数值模拟中具有实用性。

图 1 实验结果(上图)^[15]与本文模拟结果(下图)对比

Figure 1. Comparison of experimental results (upper pictures)^[15] with computational results (lower pictures)

下载: 全尺寸图片幻灯片

图 2 二维计算模型

Figure 2. 2D numerical model

下载: 全尺寸图片幻灯片

图 3 不同时刻反应物浓度

Figure 3. Distributions of reactant concentration at different times

下载: 全尺寸图片幻灯片

图 4 不同时刻的涡量和化学反应放热率分布(Ma=2.2)

Figure 4. Distributions of vorticity and heat release rate of chemical reaction at different times (Ma=2.2)

下载: 全尺寸图片幻灯片

图 5 平均涡量幅值随时间的变化

Figure 5. Time histories of mean vorticity magnitude

下载: 全尺寸图片幻灯片

图 6 平均化学反应速率随时间的变化

Figure 6. Time histories of the mean rate of chemical reaction

下载: 全尺寸图片幻灯片

图 7 η随时间的变化

Figure 7. Time histories of η

下载: 全尺寸图片幻灯片

表 1 计算区域尺寸

Table 1. Size of computational domain

Ma l₁/m l₂/m l₃/m λ₀/m

1.2 0.4 0.2 0.01 0.02

1.7 0.1 0.5 0.01 0.02

2.2 0.01 0.5 0.01 0.02

下载: 导出CSV

表 2 η的拟合参数

Table 2. Fitting perameters of η

Ma A₁ A₂ R²

1.2 1.268 -4.020 0.9150

1.7 1.465 -3.256 0.9179

2.2 1.323 -2.703 0.9082

下载: 导出CSV

参考文献(17)

[1]	Marble F E, Sonneborn G, Pun C S J, et al. Physic conditions in circumstellar gas surrounding SN 1987A 12 years after outburst[J]. The Astrophysical Journal, 2000, 545:390-398. doi: 10.1086/apj.2000.545.issue-1
[2]	Oran E S, Gamezo V N. Origins of the deflagration-to-detonation transition in gas-phrase combustion[J]. Combustion Flame, 2007, 148(1-2):4-47. doi: 10.1016/j.combustflame.2006.07.010
[3]	Yang J, Kubota T, Zukoski E E. Applications of shock-induced mixing to supersonic combustion[J]. AIAA Journal, 1993, 31(5):854-862. doi: 10.2514/3.11696
[4]	Markstein G H. A shock-tube study of flame front-pressure wave interaction[C]//6th Symposium (International) on Combustion. Pittsburgh, USA: The Combustion Institute, 1957: 387-398.
[5]	Ton V T, Karagozian A R, Marble F F, et al. Numerical simulations of high speed chemically reactive flow[J]. Theoretical and Computational Fluid Dynamics, 1994, 6:161-179. doi: 10.1007/BF00312347
[6]	Ju Y, Shimano A, Inoue O. Vorticity generation and flame distortion induced by shock flame interaction[C]//27th Symposium (International) on Combustion. Pittsburgh. USA: The Combustion Institute, 1998: 735-741.
[7]	朱跃进, 董刚, 刘怡昕, 等.激波诱导火焰变形、混合和燃烧的数值研究[J].爆炸与冲击, 2013, 33(4):430-437. doi: 10.3969/j.issn.1001-1455.2013.04.016 Zhu Yuejin, Dong Gang, Liu Yixin, et al. A numerical study on shock induced distortion, mixing and combustion of flame[J]. Explosion and Shock Waves, 2013, 33(4):430-437. doi: 10.3969/j.issn.1001-1455.2013.04.016
[8]	Zhu Yuejin, Dong Gang, Liu Yixin, et al. Three-dimensional numerical simulations of spherical flame evolutions in shock reshock accelerate flows[J]. Combustion Science and Technology, 2013, 185(10):1415-1440. doi: 10.1080/00102202.2013.798656
[9]	Khokhlov A M, Oran E S, Chtchelkanova A Y. Interaction of a shock with a sinusoidally perturbed flame[J]. Combustion and Flame, 1999, 117:99-116. doi: 10.1016/S0010-2180(98)00090-X
[10]	Kilchyk V, Nalim R, Merkle C. Laminar premixed flame fuel consumption rate modulation by shocks and expansion waves[J]. Combustion and Flame, 2011, 158:1140-1148. doi: 10.1016/j.combustflame.2010.10.026
[11]	Massa L, Jha P. Linear analysis of the Richtmyer-Meshkov instability in shock-flame interactions[J]. Physics of Fluids, 2012, 24:056101. doi: 10.1063/1.4719153
[12]	Leinov E, Malamud G, Elbaz Y, et al. Experimental and numerical investigation of the Richtmyer-Meshkov instability under re-shock conditions[J]. Fluid Mechanics, 2009, 626:449-475. doi: 10.1017/S0022112009005904
[13]	Ukai S, Balakrishnan K, Menon S. Growth rate predictions of single- multi-mode Richtmyer-Meshkov instability with reshock[J]. Shock Wave, 2011, 21:533-546. doi: 10.1007/s00193-011-0332-0
[14]	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:405-452. doi: 10.1006/jcph.2000.6443
[15]	Thomas G O, Bambrey R, Brown C. Experimental observations of flame acceleration and transition to detonation following shock-flame interaction[J]. Combustion Theory and Modelling, 2001, 5(4):573-594. doi: 10.1088/1364-7830/5/4/304
[16]	蒋华, 董刚, 陈霄, 小扰动界面形态对RM不稳定性影响的数值分析[J].力学学报, 2014, 46(4):544-552. doi: 10.7638/kqdlxxb-2013.0008 Jiang Hua, Dong Gang, Chen Xiao. Numerrical study on the effects of small amplitude initial perturbations on RM instability[J]. Chinese Jounal of Theoretical and Applied Mechanics, 2014, 46(4):544-552. doi: 10.7638/kqdlxxb-2013.0008
[17]	Latini M, Schilling O, Don W S. Effects of WENO flux reconstruction order and spatial resolution on reshocked two-dimensional Richtmyer-Meshkov instability[J]. Journal of Computational Physics, 2007, 221:805-836. doi: 10.1016/j.jcp.2006.06.051

施引文献

资源附件(0)

访问统计

图(7) / 表(2)

计量

文章访问数: 4376
HTML全文浏览量: 1443
PDF下载量: 332
被引次数: 0

1. 计算模型及材料参数
1.1 模型构建
1.2 本构模型及参数
2. 数值计算结果
3. 高速摄影实验
3.1 实验
3.2 结果
4. 结论

Ma	l₁/m	l₂/m	l₃/m	λ₀/m
1.2	0.4	0.2	0.01	0.02
1.7	0.1	0.5	0.01	0.02
2.2	0.01	0.5	0.01	0.02

Ma	A₁	A₂	R²
1.2	1.268	-4.020	0.9150
1.7	1.465	-3.256	0.9179
2.2	1.323	-2.703	0.9082

多次激波诱导正弦扰动预混火焰界面失稳的数值研究

doi: 10.11883/1001-1455(2017)02-0229-08

作者简介: 陈霄(1990-)，女，博士研究生

通讯作者: 董刚，dgvehicle@yahoo.com

计量

出版历程