旋转爆轰的三维数值模拟

2010年2月第31卷 第1期
推 进 技 术
J OURNAL OF PRO PUL SI ON TECHNOLOGY
Feb 12010
V ol 131 No 11
旋转爆轰的三维数值模拟
*
归明月,范宝春,张旭东,潘振华,姜孝海
(南京理工大学瞬态物理实验室,江苏南京210094)
摘 要:旋转爆轰利用爆轰波在燃烧室的供气端传播,而爆轰产物从另一端排出。

为了揭示其宏观特性,采用带化学反应的三维Euler 方程,在贴体坐标系中,利用带限制函数的波传播算法,对圆环形燃烧室内的旋转爆轰进行了数值研究,讨论了爆轰传播过程中的波系结构。

数值结果表明,即使仅在燃烧室内壁面的部分区域充入可燃气,爆轰波也能够以旋转方式在该区域稳定传播。

关键词:旋转爆轰;数值仿真;激波
中图分类号:O 382;V 235122 文献标识码:A 文章编号:
1001-4055(2010)01-0082-05
* 收稿日期:
2008-12-03;修订日期:2009-06-04。

基金项目:国家自然科学基金(10872096)。

作者简介:归明月(1977)),男,博士,研究领域为燃烧、爆轰等方面。

E -m ai:l gu i m i ngyue @yahoo .co
Three -di m ensional si m ul ation of conti nuous spin det onation
GU IM i n g -yue ,FAN Bao -chun ,Z HANG Xu -dong ,PAN Zhen -hua ,JI A NG X iao -hai
(Lab .o f T ransient Phy si cs ,N anji ng U n i v .of Science and T echno logy ,N anji ng 210094,Ch i na)
Abstrac t :The de tonation w ave propagates a t the c l o sed one end of t he annua l co m bustor and t he de t onati on product is e -jected a t the opened one end i n continuous spi n detonation .In order to study the m acro -charac ter ,based on the three -d i m en -s i ona l Euler equati ons coup li ng w ith che m ical reacti ons ,the phenom ena o f the conti nuous sp i n de t onati on i n an annular com -bustor were nu m er ica lly stud ied by using wave propaga ti on a l go rith m w ith li m ite r i n the body -fitted coord i nate system .The w ave structure i n t he fl ow fi e l d w as discussed acco rding to the calculated resu lts .Even the combustor are pa rtia lly filled w it h the combusti b l e gas ,t he rotating detona ti on w ave can propaga te stead il y in the layer o f prem ixed gas .
K ey word s : C onti nuous sp i n detonation ;N u m er ica l si m ulati on ;Shock w av e
1 引 言
与普通燃烧不同,爆轰是一种由激波诱导的,在
高温高压下进行的燃烧,具有数千米/秒的传播速度,故具有更高的燃烧效率、更少的污染和更宽广的飞行马赫数的适用范围。

由于现代飞行,特别是近空间飞行的需要,以爆轰为能量释放方式的爆轰发动机引起了人们的关注。

由于爆轰具有/稍纵即逝0的特点,故爆轰发动机的首要问题是如何使爆轰波长期停留在燃烧室中。

目前,有三种方式可以达到这一目的。

第一种是让爆轰以脉冲的方式和极高的频率,在燃烧室中重复出现,即在燃烧室内实现周期性的充气、点火、爆轰、排气和再充气再点火的循环过程,称此为脉冲爆轰发动机(Pu lse Detonati o n Eng i n e ,简称PDE )
[1]。

第二种是
让爆轰驻定在燃烧室内,即使爆轰相对于燃烧室处于
静止状态,称此为驻定爆轰发动机(Ob li q ue Detona -tion W ave Eng ine ,简称OD W E )。

第三种是利用爆轰波在燃烧室的供气端连续旋转,而爆轰产物从另一开口端被高速甩出进而产生推力,称此为旋转爆轰发动机(Conti n uous Sp i n Detonation Eng i n e ,简称CSDE )。

其工作原理是仅在燃烧室的供气端充入可燃气,而其它部分则是空气,当爆轰波在燃烧室内的供气端旋转时,爆轰产物在离心力和膨胀的作用下,离开爆轰波阵面,同时,新鲜的未燃气体及时充入,以维持爆轰波的持续旋转(见图1)。

旋转爆轰发动机非常像无数多个小PDE 以极高的频率(数k H z)运行,其燃烧室内的平均压力要高于单管PDE 。

早在1959年,V oitsekhovsk ii [2]
就提出,爆轰波可
以围绕圆筒持续旋转。

此后,M i k hailov [3]
分析了此类
第31卷 第1期旋转爆轰的三维数值模拟
F i g .1 S i m p lif i ed d i agra m of a CSDE
co m busti on cha mber
爆轰应该具有的波结构,但最终他们都未能实现这样的爆轰。

直至20世纪80年代,Bykovsk ii
[4,5]
等利用
气态、液态推进剂为燃料,终于在实验室实现了爆轰
波的持续旋转,并对该现象进行了初步的研究。

随后,该课题组
[6]
通过改变燃烧室的形状和尺寸、燃料
种类,对影响旋转爆轰的物理、化学和形状参数进行了较深入的探讨。

近年来,随着爆轰发动机研究的深
入,基于此类爆轰的旋转爆轰发动机的研究引起人们的关注,俄国
[2~7]
、法国
[8~10]
、波兰和日本都相继开
展了这方面的研究,并取得一定程度的成功。

国内在
这方面的研究较少,作者[11]
对该现象进行了初步的数值研究,分析了旋转爆轰流场的动力学特征,验证了旋转爆轰的可行性。

关于旋转爆轰发动机的研究大多集中在实验方面[2~6,8,9],但它们对流场中爆轰波的结构探讨较少,而关于该波系结构的数值计算方面的研究不多
[7,10]
,
且使用的都是二维模型。

由于爆轰的本质是三维的,
本文基于带化学反应的三维Euler 方程,以氢气/空气预混气为研究对象,对该现象进行了初步的数值研究。

根据计算结果分析了旋转爆轰的波系结构。

2 物理模型和计算方法
2.1 基本方程
曲线坐标系下,三维可压缩化学反应流的Eu ler 方程
5Q 5t +5F 5N +5G 5G +5H
5F
=S (1)式中Q 表示守恒变量,F ,G 和H 分别表示N ,G 和F 方向的对流通量,S 表示化学反应源项,这些矢量可以表示为
Q =
Q 1s Q K
Q u Q v Q w
E
,F =Q 1
U s Q K U Q uU +p N x Q vU +
p N y Q w U +p N z U p +E ,G =Q 1V s
Q K V Q uV +p G x Q vV +p G y Q w V +p G z V p +E ,H =
Q 1W s Q K W
Q uW +p F x Q vW +p F y Q wW +p F z W p +E
,S =ÛX 1s ÛX K 0000
(2)
式中u,v 和w 分别表示笛卡儿坐标系(x,y 和z)下的速度分量,N x i +N y j +N z k,G x i +G y j +G z k 和F x i +F y j +F z k 分别表示曲线坐标N
,G 和F 方向的单位矢量,U,V 和W 分别表示曲线坐标系(N ,G 和F )下的速度分量,即
U =u N x +v N y +w N z ,V =u G x +v G y +w G z ,
W =u F x +v F y +w F z
(3)
混合物的密度:Q =
E
K
k=1Q k
,其中Q k =Q
Y k ,Y k 是组分k 的质量分数,u,v 和w 分别是x,y 和z 方向的分速度。

p 为体系的压力,p =
E K
k=1Q k
M k
RT ,E 为体系单位体积的总内能,E =Q Q
T
0c v d T +12
Q (u 2+v 2
)+
E
K
k=1
Q k h 0
k
(0),c v =
E K
k=1
X
k c vk
,c vk 为组分k 的等容比
热,M k ,X k ,h 0
k (0)分别为组分k 的摩尔质量、摩尔分数、标准生成焓。

T 为体系的温度。

ÛX k 为组分k 的净生成速率ÛX k =
E
I
i=1
C d ki -C c ki
k f i F K
k=1
X k
C c
k i
-k bi F K
k=1
X k
C d
ki
(4)式中C c ki ,C d ki 分别表示第i 个基元反应中组分k 的正、逆反应计量系数;[X k ]为组分k 的摩尔浓度;k fi ,k bi 分别表示第i 个基元反应的正、逆反应速率常数,它们遵循Arr hen i u s 定律
k f i =A f i T B
fi exp -E f i
RT
(5)
式中A f i 表示第i 个正反应的指前因子;B fi 表示第i 个
83
推 进 技 术2010年
正反应的温度指数;E fi 表示第i 个正反应的活化能。

2.2 计算方法
本文采用分裂格式,将控制方程中的对流项和化学反应项分开处理。

守恒变量Q 的求解可表示为
Q n+1
=L S L H L G L F Q n
(6)
式中L 为操作算符,下标F,G,H 和S 代表的意义与文中前面一致。

由于计算流场中存在激波这样的强间断,因而流场的计算要能反映间断(波)的传播特性,本文的对流项采用二阶精度的波传播算法[12]。

而在化学反应过程中不同反应物质的特征时间相差很大,因此会带来/反应刚性0的计算问题,故在该过程中采用基于隐式Gear 格式的LSODE 程序计算。

时间导数项采用二步Runge -Kutta 法
[13]
Q (1)ijk =Q n ijk +$t #(L S L H L G L F U n )Q n+1
ijk =
12
[Q n ijk +Q (1)ijk +$t #(L S L H L G L F U (1)
)](7)
计算时,对流项采用无量纲量(无量纲的参考值为:压力p 0=01101325M Pa ,温度T 0=298.15K,长
度L 0=1m ),而化学反应项采用有量纲计算。

本文所采用的计算格式已经在文献[14]中得到验证。

2.3 算 例
旋转爆轰燃烧室的示意图,见图2(a)所示,为管状的圆环形燃烧室。

燃烧室内壁半径为1.4L 0,外壁为1.8L 0,高度为2.0L 0。

燃烧室中只有靠近供气端的内环区域充满可燃气,见图2(a)中的深色区域,其内径为1.4L 0,外径为1.6L 0,高度为0.1L 0。

数值计算在计算空间内进行,如图2(b)所示,点火区的点火温度为T ign ition =5.0T 0,压力为p ign ition =10.0p 0。

燃烧室内其它区域的温度和压力分别为T 0和p 0。

燃烧室左右两个边界均为周期边界,前后两个边界均为绝热固壁边界,顶端(供气端)也采用绝热固壁边界,低端采用出口边界。

采用的网格数为N @G @F =25@400@120。

对于氢气-空气混合物的反应,本文采用含19个基元反应和9组分的化学反应机理
[15]。

3 结果与讨论
图3(a)为某时刻爆轰燃烧室压力分布的三维结构图,爆轰波呈顺时针旋转,深色区域表示爆轰波阵面,它仅出现在充满可燃物的靠近供气端的内环区
域。

由于爆轰波通过可燃气与空气的交界面,(包括
Fig .2 S chematic of co m putati on al do m ain
外侧面和下侧面),在空气中形成透射激波,下侧的透射激波传向管外,外侧的透射激波在燃烧室内、外壁来回反射,形成激波系列。

为了进一步分析流场,
本文选择三个截面:进气端壁面、内环壁面和外环壁面,分别用A,B 和C 表示,图3(b)为三个截面上的压力等值线分布图,图中虚线表示初始时刻,可燃气与空气的分界面。

在进气端(A 面),分界面1(Sur -face 1)内侧充满可燃气,外侧是空气。

爆轰波仅在可燃气区域中传播,阵面基本为平面。

爆轰波在外侧的空气中形成透射激波,该透射激波在内外壁面来回反射,形成的反射激波在爆轰产物和空气中逐渐衰减成声波。

在内壁面(B 面),分界面2(Surface 2)上侧充满可燃气,下侧充满空气。

可燃气区域的情况与A 面是一样的,爆轰波在该区域稳定传播,阵面基本为平面。

爆轰波在下侧空气中形成的透射激波,向管外传播,此外,从外壁反射的激波,会在内壁形成类似的第二道激波。

在外壁面(C 面),附近的波前介质完全是空气,由于透射激波在A 面的外壁面边界形成反射,同时也会在外壁面上形成透射激波,其在空气中传播,逐渐衰减。

图4,图5和图6分别表示A,B 和C 三个平面上,爆轰产物(H 2O 的质量分数)、密度和温度的等值线分布。

在A 平面,爆轰产物在透射激波和离心力的作用下,向外侧膨胀,占据了原来为空气的外侧空间。

爆轰产物的密度和温度明显高于未燃物和空气。

84
第31卷第1期旋转爆轰的三维数值模拟
F i g.3D istribution of the pressure at t=426.5L
s
Fig.4D istr i bu ti on of the H
2
O mass
fraction at t=426.5L
s
F i g.5D istribution of the
density at t=426.5L
s
Fig.6D istri bu tion of th e
te m perature at t=426.5L
s
Fig.7D istribu tion of th e OH
mass frac ti on at t=426.5L s
同时,由于透射激波的压缩作用,其附近的密度和温
度也有所增加。

在B平面,同样在透射激波作用下,
爆轰产物向出口端膨胀,加速爆轰阵面后的稀疏过
程,爆轰产物存在区域的密度下降,温度升高。

透射
和反射激波的波后区域,密度和温度有所提高。

在C
平面,爆轰产物没有传播到该面,在该面上密度和温
度的变化皆由透射和反射激波所引起的。

图7表示流场中的OH质量分数,它反映了化学
反应阵面的形状。

从图中可以看到,爆轰波阵面上,
OH质量分数较高,而其它区域则没有。

可燃物与空
气间接触间断的存在,以及由此导致的侧向稀疏波使
反应区的外缘发生了弯曲。

图8为不同时刻流场中的压力分布。

图8(a)为
反应前的初始流场,火团在高温高压的作用下,
迅速
F i g.8Pressure d istr i bu ti on at d ifferen t ti m e
85
推进技术2010年
形成爆轰,如图8(b)所示,同时,爆轰波透过未燃物与空气的界面,在空气中形成透射激波,该激波会在外壁边界反射,见图8(c)。

在进气端形成的反射激波会在内壁边界再次反射,形成二次反射激波,同时会在内壁面形成第二道透射激波,见图8(d)。

该二次反射激波又会在外壁面边界反射,形成三次反射激波,同时在外壁面形成第二道透射激波,见图8(e)。

随后的过程不断重复前面的反射,但强度不断衰减,直致成为声波,见图8(f)~(h)。

即使可燃气仅存在于有限区域内,仍可在其中形成稳定传播的爆轰波。

4结论
基于三维多组分的Euler方程,本文对旋转爆轰进行了数值模拟。

结果表明,即使仅在爆轰燃烧室供气端的内壁面的部分区域充入可燃气,由于能量的支持,爆轰波也可以呈旋转状态在该层稳定传播,但此时,波后是复杂的三维流场,存在复杂的激波列,同时可燃物与空气接触间断引起的侧向稀疏波使爆轰波内的反应区的外缘发生了弯曲。

这些为旋转爆轰流场的深入研究以及旋转爆轰发动机的实验研究提供了有意义的参考。

参考文献:
[1]R oy G D,F ro l ov S M,Bor isov A A,et a.l Pulse detona-
tion propulsion:cha llenges,current stat ues and f uture per-
spective[J].Progress in Energy and Co m busti on Science,
2004,30(6):545~672.
[2]V o i tsekhovsk ii B V.Stationary detona tion[J].D ok lady
A kade m iiN auk UzSSR,1959,129(6):1254~1256.
[3]M i kha ilov V V,T opch i anM E.To the stud i es of conti nu-
ous de t onati on in an annu l ar channe l[J].Com bustion,
Exp losion,and Shock W ave,1965,1(4):12~14.
[4]BykovskiiF A,M itro f anov V V.D etonation combustion of
a gas m i x t ure i n a cy li ndrical chamber[J].Com bustion,
Exp losion,and Shock W ave,1980,16(5).
[5]Bykov s k ii F A,M i tro fanov V V.A conti nuous spi n deto-
na ti on i n li quid fuel sprays[C].M osco w:Control of D eto-
nation P rocesses,E le x-KM P ub l,2000:209~211.
[6]Bykovsk ii F A,Zhdan S A,V ede rnikov E F.Conti nuous
sp i n detonati ons[J].J ournal of P ropulsi on and Pow er,
2006,22(6):1204~1216.
[7]Zhdan S A,Bykovsk ii F A,V edern i kov E F.M athe m at-i
cal modeli ng of a rotati ng detona ti on wave in a hydrog en-
oxygen m i x ture[J].Co m bustion,Exp lo si on,and Shock
W ave,2007,43(4):449~459.
[8]L entsch A,Bec R,Serre L,et a.l O ve rv i ew of current
F rench acti v iti es on PDR E and conti nuous de t onati on
wave rocket eng i nes[R].A IAA2005-3232.
[9]D an i au E,F ale m pin F,Bykov s k ii F A,et a.l Conti nuous
de t onati on wave propu l s i on syste m s:first step to w ard op-
erati onal eng i nes[R].ISABE2005-1302.
[10]D an i au E,F ale m p i n F,G eti n N,et a.l D esi gn o f a con-
ti nuous de tonation w ave eng i ne for space applicati on
[R].A IAA2006-4794.
[11]姜孝海,范宝春,董刚,等.旋转爆轰流场的数值模
拟[J].推进技术,2007,28(4).(JI ANG X iao-ha,i
FAN Bao-chun,DONG G ang,et a.l N u m er ica l i nvestiga-
ti on on t he flo w field o f ro tati ng detonati on w ave[J].
Journal of Prop ulsion T echnology,2007,28(4).) [12]Leveque R J.W ave propagation a l go rith m s for m ultid-i
m ens i ona l hyperbo lic syste m s[J]putational
Phy sics,1997,131(2):327~353.
[13]Shu C W,O she r S J.Effic i ent i m p l ementation of essen-
tiall y non-osc ill a tory shock-capt ur i ng sche m es,II[J].J.
Computational Phy sics,1989,83(1):32~78.
[14]M ingyue G u,i Baochun Fan,G ang D ong,et a.l Interac-
ti on of a refl ected sho ck from a concave w a ll w ith a fla m e
disto rted by an i nc i dent sho ck[J].Shock W aves,2009,
18(6):487~494.
[15]归明月,范宝春,于陆军,等.聚心火焰与诱导激波相
互作用及爆燃转爆轰过程[J].推进技术,2007,28
(3).(GU IM i ng-yue,FAN Bao-chun,YU L u-j un,et a.l
In terac tion o f i m plosion fla m e and induced shock wave
and DDT[J].Journal of P ropulsion T echnology,2007,28
(3).)
(编辑:梅瑛)
86。