稀薄流中基于碰撞分子特定次序的DSMC改进方法

张 利,史金光,韩 宇

(1.南京理工大学 能源与动力工程学院,江苏 南京 210094;2.海军驻沈阳弹药专业军事代表室,辽宁 沈阳 110045)

现代战争对弹箭突防能力提出了更高的要求,弹箭的射程越来越远,飞行高度越来越高。随着飞行高度的增加,大气密度随之降低,当空气密度降低到一定程度时,气体的间断粒子效应使得连续介质假设不再成立,无法采用传统的N-S方程求解弹箭的气动力特性,而需要采取稀薄气体动力学[1]的方法进行研究。Bird于1963年提出了采用直接模拟蒙特卡洛(direct simulation Monte Carlo,DSMC)方法解决高空飞行器过渡流区的气动特性问题[2-4]。DSMC方法的主要思想是:使用大量模拟分子代替真实气体分子,在计算机中存储模拟分子的位置坐标、速度和内能等信息,这些运动信息随着模拟分子的运动、与边界的作用以及分子间碰撞而改变,通过空间网格对模拟分子的诸量求和以得到流场的宏观量。DSMC方法已经得到了实验验证。超声速钝头体弹箭前方的绕流流场中,气体分子碰撞将导致各种化学反应,如氮气、氧气、一氧化氮、氮原子、氧原子等组分的离解、复合、置换等复杂的化学反应[5-10],这些气体分子间的化学反应会对弹箭的受力、受热分析产生不可忽视的影响。气体分子碰撞时可能发生的化学反应并不唯一,比如:双原子分子AB与双原子分子CD进行碰撞后,既可能是AB发生离解生成2个新的单原子分子,也可能是CD发生离解形成2个新的单原子分子,等等。但事实上,任何碰撞分子可能发生的化学反应都是唯一的,化学反应类型与碰撞对分子的次序是一一对应的关系。因此,本文在Bird的气体分子化学反应DSMC模拟程序的基础上,提出了一种化学反应类型取决于碰撞对中分子特定次序的DSMC方法,并利用FORTRAN90[11]编制了DSMC仿真程序,对过渡区域高超声速圆柱绕流进行模拟分析,验证了本文方法的可行性和合理性。

1 空气中的化学反应类型

高超声速稀薄气体流动中存在着多种不同形式的热化学非平衡现象。本文研究的是五组元(N2,O2,N,O,NO)气体分子间的离解、置换和复合3种化学反应类型34种化学反应过程:

N2+X2N+X
O2+X2O+X
NO+XN+O+X
N2+ONO+N
NO+OO2+N

(1)

式中:催化物X可以是N2,O2,N,O,NO。对于化学反应模拟,Bird提出了化学反应概率模型,当组分A的分子与组分B的分子发生碰撞时,以一定的概率发生化学反应。通常引入反应截面σR,而σR与总碰撞截面σT的比值σR/σT代表弹性碰撞导致化学反应的概率,并称为位阻因子。同时,碰撞中的总能量应超过分子的活化能。

①离解反应。空气中的离解反应有:

O2+X→O+O+X
NO+X→N+O+X
N2+X→N+N+X

(2)

式中:催化物X可以是N2,O2,N,O,NO,所以共有15种离解反应。

离解反应在DSMC模拟中的实现通常是利用Larsen-Bergnakke的振动松弛理论,当双原子分子振动级数激发到大于该分子离解能所对应的振动级数,则认为该双原子分子发生了离解,并标记该双原子分子以进行离解反应的计算。

②复合反应。空气中的复合反应主要有:

N+N+X→ N2+X
O+O+X→ O2+X
N+O+X→NO+X

(3)

复合反应在DSMC模拟中的实现通常是利用配分函数与细致平衡原理确定抽样概率,并根据积累概率法判断复合反应是否发生,当网格中的累计概率大于1时,即认为发生一次复合反应,标记该网格中的2个原子以进行复合反应的计算,并将该网格中的累计概率减去1。计算公式为

(4)

式中:a1,b1为复合反应常数,T为温度,n为数密度。

③置换反应。空气中的交换原子化学反应都是双原子分子和单原子分子间进行碰撞交换原子。空气中的置换反应有以下4种类型:

N2+O→NO+N
NO+N→N2+O
NO+O→O2+N
O2+N→NO+O

(5)

置换反应在DSMC模拟中的实现通常是利用配分函数与细致平衡原理确定置换反应抽样概率,并根据取舍法判断是否发生离解反应。若抽样概率大于随机数Rf,则发生置换反应,并标记该双原子分子和自由原子,若抽样概率小于随机数Rf,则进行下一碰撞对的判断。计算公式为

(6)

式中:a2,b2为置换反应常数。

2 基于碰撞分子特定次序的DSMC改进方法

在DSMC方法中,化学反应发生在分子对的碰撞过程中,Bird给出的方法是逐一判断碰撞对分子所有可能发生的化学反应。值得注意的是,分子碰撞时可能发生的化学反应并不唯一,没有必要逐一判断碰撞对所有可能的化学反应,且该过程效率低下。为此,本文提出限定碰撞对中分子的特定次序的DSMC改进方法,提高计算效率。

①判断碰撞对中第1个分子是否为双原子分子,若不是双原子分子,则转到②;若是双原子分子,根据抽样反应概率判断该分子是否发生离解反应,若发生了离解,则标记该碰撞对,若未发生离解反应,则进行下一碰撞对的判断。

②第1个分子为单原子分子,判断第2个分子是否为双原子分子,若不是双原子分子,则转到③;若是双原子分子,则根据抽样反应概率判断第2个分子是否发生置换反应,如果发生了置换,则标记该碰撞对,如果未发生置换,则进行下一碰撞对的判断。

③在上述步骤中判断出碰撞对均为单原子分子,则利用抽样概率判断该碰撞对是否发生复合反应,如果发生了单原子复合,则标记该碰撞分子对,如果未发生置换反应,则进行下一碰撞对的判断。

仿真过程如图1所示。图中,Rf为随机数,Pcol,Prec,Pdis,Pexch分别表示碰撞概率、复合反应概率、离解反应概率和置换反应概率。

图1 化学反应的DSMC流程图

3 数值算例

3.1 计算模型与计算条件

现代战争对弹箭突防能力提出了更高的要求,弹箭的射程越来越远,飞行高度越来越高。随着弹箭飞行高度的增加,流场的稀薄效应随之增强,弹箭的受力、受热分析无法采用连续流理论,需要采取稀薄气体动力学的方法。本节在直角坐标网格框架下,在间断能级的Larsen-Borgnakke碰撞模型和Bird的化学反应概率模型基础上,利用DSMC方法模拟了弹箭飞行速度为7.5 km/s,飞行高度h分别为80 km,85 km,90 km时的圆柱绕流,为后续弹箭在稀薄流区的受力、受热分析提供理论基础。本文来流条件参考了文献[12],在80 km和85 km处,来流组分摩尔分数为x(N2)=0.763,x(O2)=0.237,在90 km处,x(N2)=0.788,x(O2)=0.209,x(O)=0.003,物面温度均为1 000 K。

3.2 计算结果及分析

将壁面处的热流密度模拟结果无量纲化,得到热流系数Ch的表达式:

(7)

式中:q为热流密度,ρ为来流密度,u∞为来流速度。

图2是本文计算的沿圆柱表面热流系数Ch和沿驻点线流场温度T的分布曲线。图2(a)～2(c)中原点处为驻点,横坐标是圆柱表面不同位置θ,纵坐标是热流系数Ch;图2(d)中横坐标是距壁面的距离d。表1为本文数值计算结果和文献[12]的比较,结果表明,在80 km,85 km和90 km处,热流系数峰值的误差分别1.1%,2.0%和1.9%,结果吻合得较好,并且热流系数峰值随高度增加而增大,其值从80 km处的0.380增大到90 km处的0.632。从热流系数曲线还可以看出,在驻点附近热流密度有最大值,过了驻点后热流密度逐渐减小,在尾部热流密度非常小,接近于0。图3是壁面热流系数Ch和摩阻系数Cf随高度变化曲线,可以看出热流密度和摩阻系数随高度变化比较敏感。

表1 不同高度圆柱表面热流密度和热流系数峰值

图4(a)是流场平动温度、转动温度和振动温度沿驻点线的分布,可以看出壁面附近存在明显的热力学非平衡现象。由此不难得出,这一区域存在脱体的弓形激波,激波内流场振动温度明显低于平动温度和转动温度,且振动温度在此两者之后开始上升,这是由于振动松弛速率要远小于平动与转动松弛的速率,气体内分子的振动平衡需要更多次的分子碰撞,振动达到平衡就会比平动和转动要迟一些。

图2 计算结果比较

图3 稀薄程度对壁面热流系数和摩阻系数的影响

图4 稀薄程度对驻点线物理量的影响

图4(b)给出了80 km高空处沿驻点线各气体组分的摩尔分数xB(B分别代表N2,O2,NO,N)分布曲线。从图中可以看出,靠近壁面过程中,由于气体分子密度逐渐增加,分子碰撞频率逐渐增大,离解反应、复合反应和置换反应发生的概率增大,使得氮气和氧气摩尔分数逐渐降低,氮原子和氧原子逐渐增加,靠近壁面时冷壁效应导致氮气和氧气摩尔分数略有升高。值得注意的是,来流气体分子在流场温度达到峰值前就发生了离解反应,表明在高空稀薄流条件下,高超声速弹箭壁面前方激波层附近存在明显的化学反应,这对弹头部的受力、受热有着很大影响。图4(c)是氮原子摩尔分数随高度的分布情况,在90 km高空处,基本上没有氮原子的生成,流场化学反应变得很微弱,可以得出弹箭表面的化学反应随高度增加而减弱,飞行高度达到自由分子流区时可以不再考虑化学反应。

4 结论

本文详细讨论了包括离解反应、复合反应和置换反应在内的3种化学反应类型,基于Bird建议的化学反应抽样概率模型,给出了这3种化学反应在DSMC方法中的实现过程。在Bird给出的模拟二维圆柱绕流的化学反应DSMC程序的基础上,提出了化学反应类型取决于分子碰撞对特定次序的DSMC改进方法,并以高超声速圆柱绕流为例,在80 km,85 km,90 km高空时模拟了五组元混合气体的32种化学反应,在80 km,85 km和90 km时的热流系数峰值的误差分别为1.1%,2.0%和1.9%,与理论数据吻合得较好,表明本文提出的改进方法是合理可行的,且效率较高,为后续高超声速弹箭在高空稀薄流区的受力、受热分析提供了理论基础。