CN104050334B - 一种火箭羽流仿真方法 - Google Patents

Abstract

本发明公开了一种火箭羽流仿真方法，其步骤如下：一、利用Gambit软件生成非结构化网格，并设置边界条件；二、将网格导入到Fluent中，利用Fluent软件求解羽流流场物理模型：湍流模型、离散相模型和燃烧模型；三、在Fluent中根据模型需要或实际情况设置边界条件和迭代初始值，模拟仿真获得流场数据。本发明基于火箭羽流流场的形成机理研究，根据燃烧理论、流体力学、气体动力学，考虑燃烧室内的燃烧化学反应建立的羽流流场计算的参数模型，能够客观描述羽流流场的主要特征，快速仿真获得特定条件下的流场数据，为火箭羽流红外特性的计算提供输入参数，提高羽流的计算精度、节省羽流的仿真时间。

Description

一种火箭羽流仿真方法

技术领域

本发明属于目标特性与识别技术领域，涉及一种火箭羽流仿真模型。

背景技术

火箭羽流流场的计算是为羽流的红外辐射场的计算提供基础数据。羽流流场的计算主要包括速度场、温度场、压力场的计算。尾流问题的研究方法多种多样，大致可分为如下三种：半经验和分析法、试验测量法以及数值模拟法。半经验和分析法主要是应用于比较简单的几何模型，在物理过程不太复杂、流动可视为线性的情况下，利用简化之后得到的半经验半理论公式或分析解，研究流动的问题，为飞行器的设计提供理论依据，但是随着人们对飞行器的性能要求越来越高，所需要研究的尾流流场问题就会越来越复杂，同时飞行器尾部的几何形状也更加复杂，流动模型也大多数是湍流的，如此半经验和分析法不能满足设计当今飞行器的要求。

试验测量法是利用实物模型飞行试验的方法模拟出飞行器飞行过程中的尾流，通过在试验中记录所需物理量，从而得到人们关心的数据，为计算设计各种飞行器，提供直接的理论数据，试验测量法在飞行器的设计与研制过程中起着非常重要的作用，可以说，试验是航空航天科技发展的基础，但是由于现代航空航天技术的飞快发展，对飞行器研制的经费、周期等的要求也越来越苛刻，如此说来，试验测量法已不再是设计飞行器的唯一方法，当然，试验测量法在设计飞行器中仍然具有十分重要的作用和重大的意义。

数值模拟法是建立在计算机和计算流体力学(CFD)的发展基础之上，为了满足飞行器设计快、好、省、多的要求而逐渐被人们所接收并采纳的一种设计和研究的方法，数值模拟的方法不仅能提供多种设计研究方案供设计者选择，还能快速提供设计者想要得到的数据结果，在羽流问题的研究上，由于羽流问题的复杂多样，用实物试验无法模拟所有情形，而数值模拟法则能克服并解决这一难题，为设计者提供所需的各种情形下的羽流流场参数，同时可以得到一些用实物试验无法模拟计算得到的数据，这都是实物试验法无法比拟的，除此之外，数值模拟的方法还能节省大量的人力、物力和财力，意义重大。

为使计算简单且很快得到羽流流场参数，可采用通过计算机进行的数值模拟的方法，即计算流体力学CFD方法，该方法分析流体流动和传热等物理现象。通过CFD技术，可利用计算机分析并显示流场中的现象，从而在较短时间内预测流场。CFD模拟能帮助理解流体力学问题，为实验提供指导，为设计提供参考，从而节省人力、物力和时间。然而，传统的假定火箭喷管出口参数已知，直接进行羽流流场参量计算的建模仿真，计算效率低、可靠性差。

发明内容

本发明的目的是提供一种火箭羽流仿真方法，该方法基于火箭羽流流场的形成机理研究，根据燃烧理论、流体力学、气体动力学，考虑燃烧室内的燃烧化学反应建立的羽流流场计算的参数模型，能够客观描述羽流流场的主要特征，快速仿真获得特定条件下的流场数据，为火箭羽流红外特性的计算提供输入参数，提高羽流的计算精度、节省羽流的仿真时间。

本发明的目的是通过以下技术方案实现的：

一种火箭羽流仿真方法，从火箭推力室开始研究，分析火箭羽流流场的形成过程，考虑燃烧室内的燃烧化学反应，根据燃烧理论、流体力学、气体动力学，得到了羽流流场的参数计算公式。具体包括如下内容：

一、利用Gambit软件生成分离比例为0.1-0.5的非结构化网格，并设置如下四种边界条件：入口边设为速度入口边界；出口设为压力输出边界；下边设为对称边界条件；上边设为壁面边界条件。

二、将网格导入到Fluent中，利用Fluent软件求解羽流流场物理模型，主要包括湍流模型、离散相模型和燃烧模型。

三、在Fluent中根据模型需要或实际情况设置边界条件和迭代初始值，模拟仿真获得流场数据。

本发明中在Fluent中进行了2000次以上的迭代仿真运算得到了迭代残差小于1e-5的计算结果，计算得到火箭羽流流场的速度场、温度场、压力场。计算获得的速度场、压力场均符合实际情况，以温度场为例进行计算精度分析，得到温度场仿真结果误差为2％，并且得到了连续的易于观察的场数据分布。与传统的假定火箭喷管出口参数已知，直接进行羽流流场参量计算的建模仿真方法比较，在一定程度上提高了计算效率、精度和可靠性。总之，该方法为火箭羽流流场计算研究提供了一种新思路，可作为工程计算的一种有效模型和算法。

本发明有如下优点：

1、能够快速仿真获得特定条件下的流场数据，提高了计算效率。

2、能够获得连续的易于观察的场数据分布，比起以往成果中的离散数据，本发明得到更加直观、可靠的流场数据。

3、从火箭推力室开始研究，考虑燃烧室内的燃烧化学反应，在一定程度上提高了精度和可靠性。本发明中计算的温度场仿真结果误差为2％。

附图说明

图1为火箭羽流流场计算总体方案；

图2为Gambit生成的网格；

图3为速度场分布；

图4为温度场分布；

图5为温度场计算结果与实验值沿轴线温度对比图；

图6为压力场分布。

具体实施方式

下面结合附图对本发明的技术方案作进一步的说明，但并不局限于此，凡是对本发明技术方案进行修改或者等同替换，而不脱离本发明技术方案的精神和范围，均应涵盖在本发明的保护范围中。

本发明提供了一种火箭羽流仿真方法，具体包括如下内容：

一、火箭羽流流场仿真方法的总体说明

研究思路是首先利用Gambit软件生成计算模拟所用的非结构化网格，并且进行边界条件的定义，以实现Fluent的模拟仿真。利用Fluent计算流体软件对羽流流场进行模拟，得到火箭羽流流场的速度场、温度场、压力场。为后续火箭羽流辐射场的计算分析提供基础数据。具体研究方案如图1所示。

二、火箭羽流流场仿真模型的分布说明

1、Gambit中的网格划分

计算流场数据的第一步是划分网格，用学术语言表达计算区域的离散化，即将空间上连续的计算区域划分成许多个子区域，并确定每个区域中的节点。数学上，生成网格后(即离散化后)，就将连续的控制方程进行了离散，即将描写流动与传热的偏微分方程转化为各个节点上的代数方程组。划分网格的本质是用有限个离散的点来代替原来的连续空间。

本发明采用了非结构化划分网格，并将分离比例分别设为0.5、0.3、0.1的划分网格进行迭代计算对比，采用了0.5分离比例的划分网格的模拟结果，0.5分离比例的划分网格如图2所示。本发明中的网格定义了四种边界条件：入口边设为速度入口边界；出口设为压力输出边界；为减少计算时间，只生成了一半网格，所以下边设为对称边界条件；而上边的不规则边设为壁面边界条件。

2、羽流仿真计算的物理模型

羽流流场物理模型主要包括湍流模型，使用离散相模型来表达液滴喷雾过程，并与连续相流场耦合模拟出了液滴的挥发和燃烧过程。挥发燃烧过程由混合分数平衡化学反应模型模拟。下面将分别详细介绍每个物理模型。

2.1羽流流场的湍流模型

本发明采用湍流模型中的k-ε模型，由于其微分方程形式与N-S方程相同，因而可以和N-S方程进行耦合求解。本发明采用高雷诺数k-ε模型来求解。

设ρ为密度，u、v为x、y向速度，p为压强，e为单位质量的总能量，T为温度，k，ε分别为湍动能和湍流耗散速度。

式中，Re-雷诺数；S_Pφ(φ＝u，v，e)-颗粒相源项，对于液体燃料的火箭羽流，其值为0。

其中各应力项为：

x、y向热流项分别为：

其中，γ为粘性系数μ_e与其密度ρ的比值，Pr、Pr_t分别为层流和湍流的普朗特数，物理意义是流体的动量扩散能力与流体的热量扩散能力的比，μ_e＝μ_l+μ_t为有效粘性系数，包含层流μ_l与湍流μ_t两部分之和。湍流粘性计算如下：

ε_c是考虑可压缩性的湍流耗散速率，根据Sarkar和Lakshmanan模型有：

此处是湍流马赫数。湍动能生成项P_k计算如下：

是考虑压力扩张的湍流相关项。

湍流模型常数如下所示：

C_μ＝0.09，C₁＝1.60，C₂＝1.92，

σ_k＝1.0，σ_ε＝1.3，

α₁＝1.09，α₂＝0.4，α₃＝0.2

计算中还要用到的其它关系式有，完全气体状态方程：

p＝ρRT。

单位质量的总能量：

2.2羽流流场的离散相模型

Fluent中，离散相可用于模拟液滴流动。离散相模型能抓住独立液滴的轨迹。液体燃料和空气流的热量、动量和质量传递可以计算出来。本发明中允许离散相和连续相耦合，这样模拟的离散相能考虑气相运动的影响，也会考虑两相的热量和质量传递。

Fluent采用有限体积法将非线性偏微分方程转变成网格单元上的线性代数方程，然后通过求解线性方程组得到流场的解。从物理角度看，方程的守恒形式反映的是流场变量在网格单元上的守恒关系，即网格单元内某个流场的增量等于各边界上变量的通量的总和。有限体积法的求解策略是用边界面或线上的通量计算出控制点上的变量。

2.3羽流流场的燃烧模型

Fluent中提供了5种模拟反应的方法，本发明中采用非预混合燃烧模型，在这种方法中并不是解每一个组分输运方程，而是解一个或两个守恒标量(质量分数)的输运方程，然后从预测的的混合分数分布推导出每一个组分的浓度。该模型主要用于模拟湍流扩散火焰。对于有限速度公式来说，这种方法有很多优点。在守恒标量方法中，通过概率密度函数或者PDF来考虑湍流的影响。

在非预混燃烧中，燃料和氧化剂以相异流进入反应区。在一定假设条件下，热化学可被简化成一个单一的参数，即混合分数。换句话说，混合分数是在所有组分里，燃烧和未燃烧燃料元素的局部质量分数。因为化学反应中元素是守恒的，所以这种方法极好。反过来，质量分数是一个守恒的数量，因此其控制输运方程不含源项。燃烧被简化成一个混合问题，并且与近非线性平均反应率相关的问题可以避免。

Fluent中关于非预混合燃烧模型解两个关于平均混合分数的方程：平均(时间平均)混合分数方程和平均混合分数均方值的守恒方程。平均(时间平均)混合分数方程为：

式中f为混合分数源项，S_m仅指质量由液体燃料滴或反应颗粒传入气相中，S_user为任何用户定义源项。

Fluent也解一个关于平均混合分数均方值的守恒方程：

式中常数σ_t、C_g和C_d-分别取0.85，2.86和2.0；S_user-用户定义源项。

3、基于Fluent的火箭羽流流场模拟

Fluent软件可以在较短时间内完成流场计算，为实验提供指导，为设计提供参考，从而节省人力、物力和时间。下面是对某型号火箭羽流流场模拟仿真，具体计算结果如图3～图6所示，本发明以温度场为例进行计算精度分析。

速度场计算结果如图3所示。从仿真结果可以看出，火箭羽流速度分布在0-1030m/s区间内，在外界空间中，速度大小呈扩散式变化，查阅资料可以知道火箭羽流的最大速度在2马赫到4马赫之间，本发明仿真计算得到最大速度为3马赫，符合实际情况，可为流场数据提供理论依据，模拟仿真结果比较理想。

温度场计算结果如图4所示。从计算结果可以看出，在燃烧室内富燃区域温度最高，达到2350k，在喷管出口处温度呈扩散式变化，能反映出羽流流场的主要特点。查阅资料可以知道在火箭尾部的燃烧区域能够产生高温2000-3000k，本发明仿真结果的最高温度在这个区间内，并且将所查资料中的热像图实验值相比，其相对误差值为2％左右。

其中T_max1为仿真计算的温度最大值，T_max2为热像图中实验最大值，温度场的数值模拟的结果与实验结果比较，误差在可接受范围内，如图5所示，仿真结果可靠，可为流场数据提供理论依据。

压力场计算结果如图6所示。压力分布场中可以看出，燃烧室内由于气体的燃烧消耗，压强最高，但是总体压力均低于标准大气压，在喷管出口处压力呈扩散式分布，仿真误差在可接受的范围内，仿真模拟结果符合实际情况，能反映出羽流流场的主要特点。

Claims (1)

1.一种火箭羽流仿真方法，其特征在于所述仿真方法为：

一、利用Gambit软件生成非结构化网格，并设置如下四种边界条件：入口边设为速度入口边界；出口设为压力输出边界；下边设为对称边界条件；上边设为壁面边界条件；

二、将网格导入到Fluent中，利用Fluent软件求解羽流流场物理模型：湍流模型、离散相模型和燃烧模型，其中：

所述湍流模型采用高雷诺数k-ε模型来求解，则有：

$\frac{\∂U}{\∂t}+\frac{\∂F}{\∂x}+\frac{\∂G}{\∂y}=\frac{\∂F_v}{\∂x}+\frac{\∂G_v}{\∂y}+S,$

$\begin{array}{ccc}U=\left[\begin{array}{c}\mathrm{\ρ}\\ \mathrm{\ρ}u\\ \mathrm{\ρ}v\\ \mathrm{\ρ}e\\ \mathrm{\ρ}k\\ \mathrm{\ρ}\mathrm{\ϵ}\end{array}\right]& F=\left[\begin{array}{c}\mathrm{\ρ}u\\ {\mathrm{\ρu}}^2+p\\ \mathrm{\ρ}uv\\ \left(\mathrm{\ρ}e+p\right)u\\ \mathrm{\ρ}uk\\ \mathrm{\ρ}u\mathrm{\ϵ}\end{array}\right]& G=\left[\begin{array}{c}\mathrm{\ρ}v\\ \mathrm{\ρ}uv\\ {\mathrm{\ρv}}^2+p\\ \left(\mathrm{\ρ}e+p\right)v\\ \mathrm{\ρ}vk\\ \mathrm{\ρ}v\mathrm{\ϵ}\end{array}\right]\end{array},$

$\begin{array}{ccc}{F}_v=\frac{1}{R_e}\left[\begin{array}{c}0\\ {\mathrm{\τ}}_{xy}\\ {\mathrm{\τ}}_{yy}\\ {\mathrm{\τ}}_{xy}u+{\mathrm{\τ}}_{yy}v-q_y\\ {\mathrm{\τ}}_{yk}\\ {\mathrm{\τ}}_{y\mathrm{\ϵ}}\end{array}\right]& G_v=\frac{1}{R_e}\left[\begin{array}{c}0\\ {\mathrm{\τ}}_{xy}\\ {\mathrm{\τ}}_{yy}\\ {\mathrm{\τ}}_{xy}u+{\mathrm{\τ}}_{yy}v-q_y\\ {\mathrm{\τ}}_{yk}\\ {\mathrm{\τ}}_{y\mathrm{\ϵ}}\end{array}\right]& S=\left[\begin{array}{c}0\\ S_{Pu}\\ S_{Pv}\\ S_{Pe}\\ P_k-\mathrm{\ρ}\left(\mathrm{\ϵ}+{\mathrm{\ϵ}}_c\right)R_e+\stackrel{\‾}{p^{\′\′}{d}^{\′\′}}\\ \frac{\mathrm{\ϵ}}{k}\left(C_1{P}_k-C_2{\mathrm{\ρ\ϵR}}_e\right)\end{array}\right]\end{array},$

式中，Re-雷诺数；S_Pφ(φ＝u,v,e)-颗粒相源项，对于液体燃料的火箭羽流，其值为0；ρ为密度；u、v为x、y向速度；p为压强；e为单位质量的总能量；T为温度；k,ε分别为湍动能和湍流耗散速度；

其中各应力项为：

$\begin{array}{c}{\mathrm{\τ}}_{xx}={\mathrm{\μ}}_e\[-\frac{2}{3}\left(\frac{\∂u}{\∂x}+\frac{\∂v}{\∂y}\right)+2\frac{\∂u}{\∂x}\]\\ {\mathrm{\τ}}_{yy}={\mathrm{\μ}}_e\[-\frac{2}{3}\left(\frac{\∂u}{\∂x}+\frac{\∂v}{\∂y}\right)+2\frac{\∂v}{\∂y}\]\\ {\mathrm{\τ}}_{xy}={\mathrm{\τ}}_{yx}={\mathrm{\μ}}_e\left(\frac{\∂u}{\∂y}+\frac{\∂v}{\∂x}\right)\\ {\mathrm{\τ}}_{xk}=\left({\mathrm{\μ}}_l+\frac{{\mathrm{\μ}}_t}{{\mathrm{\σ}}_k}\right)\frac{\∂k}{\∂x},{\mathrm{\τ}}_{yk}=\left({\mathrm{\μ}}_l+\frac{{\mathrm{\μ}}_t}{{\mathrm{\σ}}_k}\right)\frac{\∂k}{\∂y}\\ {\mathrm{\τ}}_{x\mathrm{\ϵ}}=\left({\mathrm{\μ}}_l+\frac{{\mathrm{\μ}}_t}{{\mathrm{\σ}}_{\mathrm{\ϵ}}}\right)\frac{\∂\mathrm{\ϵ}}{\∂x},{\mathrm{\τ}}_{y\mathrm{\ϵ}}=\left({\mathrm{\μ}}_l+\frac{{\mathrm{\μ}}_t}{{\mathrm{\σ}}_{\mathrm{\ϵ}}}\right)\frac{\∂\mathrm{\ϵ}}{\∂x}\end{array},$

x、y向热流项分别为：

$q_x=-\frac{1}{\mathrm{\γ}-1}(\frac{{\mathrm{\μ}}_l}{\mathrm{Pr}}+\frac{{\mathrm{\μ}}_t}{{\mathrm{Pr}}_t})\frac{\∂T}{\∂x},q_y=-\frac{1}{\mathrm{\γ}-1}(\frac{{\mathrm{\μ}}_l}{\mathrm{Pr}}+\frac{{\mathrm{\μ}}_t}{{\mathrm{Pr}}_t})\frac{\∂T}{\∂y},$

其中，γ为粘性系数μ_e与其密度ρ的比值，Pr、Pr_t分别为层流和湍流的普朗特数，物理意义是流体的动量扩散能力与流体的热量扩散能力的比，μ_e＝μ_l+μ_t为有效粘性系数，包含层流μ_l与湍流μ_t两部分之和；

湍流粘性计算如下：

${\mathrm{\μ}}_t=C_{\mathrm{\μ}}\frac{{\mathrm{\ρk}}^2}{\mathrm{\ϵ}+{\mathrm{\ϵ}}_c};$

ε_c是考虑可压缩性的湍流耗散速率，根据Sarkar和Lakshmanan模型有：

${\mathrm{\ϵ}}_c={\mathrm{\α}}_1{M}_t^2\mathrm{\ϵ},$

此处是湍流马赫数，湍动能生成项P_k计算如下：

$P_k=\frac{{\mathrm{\μ}}_t}{R_e}\{2\[{\left(\frac{\∂u}{\∂x}\right)}^2+{\left(\frac{\∂v}{\∂y}\right)}^2\]+{(\frac{\∂u}{\∂y}+\frac{\∂v}{\∂x})}^2\},$

是考虑压力扩张的湍流相关项；

湍流模型常数如下所示：

$\begin{array}{ccc}{C}_{\mathrm{\μ}}=0.09,& C_1=1.60,& C_2=1.92,\\ {\mathrm{\σ}}_k=1.0,& {\mathrm{\σ}}_{\mathrm{\ϵ}}=1.3,& \\ {\mathrm{\α}}_1=1.09,& {\mathrm{\α}}_2=0.4,& {\mathrm{\α}}_3=0.2\end{array};$

完全气体状态方程：

p＝ρRT；

单位质量的总能量：

$e=\frac{p}{\mathrm{\ρ}(\mathrm{\γ}-1)}+\frac{u^2+v^2}{2};$

所述离散相模型允许离散相和连续相耦合；

所述燃烧模型采用非预混合燃烧模型，非预混合燃烧模型中平均混合分数方程为：

$\frac{\∂}{\∂t}\left(\mathrm{\ρ}\stackrel{\‾}{f}\right)+\▿\·\left(\mathrm{\ρ}\stackrel{\‾}{v}\stackrel{\‾}{f}\right)=\▿\·(\frac{{\mathrm{\μ}}_t}{{\mathrm{\σ}}_t}\▿\stackrel{\‾}{f})+S_m+S_{user},$

式中f为混合分数源项，S_m仅指质量由液体燃料滴或反应颗粒传入气相中，S_user为任何用户定义源项；

平均混合分数均方值的守恒方程为：

$\frac{\∂}{\∂t}\left(\mathrm{\ρ}\stackrel{\‾}{f^{\′2}}\right)+\▿\·\left(\mathrm{\ρ}\stackrel{\‾}{v}\stackrel{\‾}{f^{\′2}}\right)=\▿\·(\frac{{\mathrm{\μ}}_t}{{\mathrm{\σ}}_t}\▿\stackrel{\‾}{f^{\′2}})+C_g{\mathrm{\μ}}_t\left({\▿}^2\stackrel{\‾}{f}\right)-C_d\mathrm{\ρ}\frac{\mathrm{\ϵ}}{k}\stackrel{\‾}{f^{\′2}}+S_{user},$

式中常数σ_t、C_g和C_d-分别取0.85，2.86和2.0；S_user-用户定义源项；

三、在Fluent中根据模型需要或实际情况设置边界条件和迭代初始值，模拟仿真获得流场数据,

所述燃烧模型为非预混合燃烧模型,所述非预混合燃烧模型中，时间平均混合分数方程为：

$\frac{\∂}{\∂t}\left(\mathrm{\ρ}\stackrel{\‾}{f}\right)+\▿\·\left(\mathrm{\ρ}\stackrel{\‾}{v}\stackrel{\‾}{f}\right)=\▿\·(\frac{{\mathrm{\μ}}_t}{{\mathrm{\σ}}_t}\▿\stackrel{\‾}{f})+S_m+S_{user},$

式中，f为混合分数源项，S_m仅指质量由液体燃料滴或反应颗粒传入气相中，S_user为用户定义源项。