采用不同气相燃烧模型模拟含硼燃气扩散燃烧过程①

冯喜平，李海波，唐金兰，李贵珠，廖自繁

(1.西北工业大学燃烧、热结构与内流场重点实验室，西安 710072;2.西北工业集团有限公司，西安 710061)

0 引言

含硼富燃燃气与空气的燃烧是固冲发动机技术基础和关键，也是其研究重点［1］。由于硼是固冲发动机理想的能源，硼的点火和燃烧是固冲发动机研究基础。研究者普遍认为，硼粒子的点火和燃烧由点火过程和燃烧过程［2］2个连续的阶段组成，并提出不同的理论。硼粒子点火模型以 King 模型［3－4］和 Li，Williams 和Kuo建立的模型［5］为代表;硼粒子的燃烧模型则以LW模型和Y-K模型为代表，其中Makino和Law建立了无氧化层包覆的硼颗粒表面反应燃烧模型［6］，L-W燃烧模型考虑了扩散和动力学机理［7］。二次燃烧室本质是一能量转换装置，其中的燃烧是研究重点［8－9］，直接关系到发动机性能。

在进行发动机燃烧过程模拟中，需要选择燃烧模型，而不同燃烧模型反映不同燃烧机理，关联燃烧模拟的结果，燃烧模型的确定成为成功进行模拟的关键。针对目前研究普遍存在模型选择的困惑，依靠实验确定正确燃烧模型成为进行二次燃烧研究的有效方法。

本文以实验装置为物理模型，在分析硼粒子King点火模型和L-W燃烧模型反应机理的基础上，湍流燃烧分别采用有限速率/涡耗散模型、涡耗散模型、有限速率模型，建立含硼富燃燃气湍流扩散燃烧模型并开展数值模拟，通过对模拟结果的分析和对比，研究不同模型下的燃烧，为数学模型建立提供支持。

1 燃烧机理分析

燃烧模型能否反映含硼富燃燃气燃烧机理，特别是能否正确反映硼粒子点火和燃烧过程，是正确进行数值模拟的关键。本章分析气相组分的燃烧机理和硼粒子点火与燃烧模型的反应机理，为数值模拟研究做基础［10］。

1.1 气相燃烧机理

由于硼蒸气、B2O2在空气中反应复杂，气相计算中忽略硼与H2O的反应，忽略B2O2与H2O的反应。因此，气相组分燃烧采用以下反应:

关于反应速率的确定，本文分别采用有限速率/涡耗散模型(Finite-Rate/Eddy-Dissipation模型)、涡耗散模型(Eddy-Dissipation模型)及有限速率模型(Laminar Finite-Rate模型)［11］，开展不同模型燃烧特性的数值模拟对比。有限速率模型简化了湍流脉动的影响，反应速率根据Arrhenius公式确定。由于湍流火焰中Arrhenius化学动力学的高度非线性，这一模型对于化学反应剧烈，湍流脉动较大的湍流扩散火焰的计算一般不精确。涡耗散模型称之为湍流-化学反应相互作用模型，反应速率由湍流混合时间尺度k/ε控制，避免Arrhenius化学动力学计算。在非预混火焰反应区发生快速燃烧时，只要湍流出现，反应即可开始不受限制，反应速度往往较快。有限速率/涡耗散模型被广泛应用于湍流扩散燃烧的数值模拟。净反应速率由Arrhenius化学动力学和涡耗散反应速率混合控制，Arrhenius反应速率作为动力学开关，阻止反应在火焰稳定之前发生，延迟了计算中点火的开始，较为符合实际［11］。

1.2 硼颗粒点火模型

King模型是广大研究者普遍接受的硼粒子点火模型。该模型中，假设硼粒子是均质球体，氧化层为厚度均匀的B2O3，环境中的氧经扩散作用通过氧化层扩散至硼与B2O3交界面，与硼反应生成B2O3氧化层。氧化层的消耗由B2O3氧化层的蒸发过程和B2O3在外表面与水蒸气的反应过程组成。当氧化层消耗完毕，点火过程完成。

硼的消耗反应为

反应速率为RB，粒径变化方程为［8］

硼表面氧化层的蒸发过程如下:

B2O3(l)→B2O3(g)+ ΔHvap，反应速率为 RE。B2O3(l)+H2O(g)→HOBO(g)+ ΔHH，反应速率为RH。氧化层的减小速率是由3个过程中B2O3的生成和消耗决定的。

式中 rp、δ、f、T0、Tp、TRAD分别表示硼粒子球体半径、硼粒子包覆氧化层厚度、硼粒子融化时液态硼质量分数、环境温度、硼粒子温度及周边辐射温度。

方程(3)、(6)表示固态硼(硼粒子温度低于硼熔点条件下)和液态B(硼粒子温度高于B熔点条件下)条件下硼粒子温度随时间变化。方程(4)表示硼粒子熔化过程中(硼粒子温度等于B熔点)液态B质量分数随时间的增加过程。QRX1、QRX2分别为固态B和液态B的氧化反应热。各变量的定义与取值方法，以及各方程表示的化学过程见文献［9，12］。

1.3 硼颗粒燃烧模型

硼粒子点火过程结束后，硼粒子直接与氧气接触燃烧，进入燃烧过程。硼作为固体单质发生表面燃烧反应，也会有部分硼受热融化，进而有蒸发现象的存在。硼作为单质，没有挥发析出。硼的蒸发量由梯度扩散确定，其蒸发速率见文献［11］。

Makion和Law［10］认为，单质硼的燃烧包括2个阶段:首先，硼与O2在表面反应生成B2O2;随后，气相B2O2与O2在相当高的温度下反应生成B2O3。L-W模型认为，硼粒子燃烧阶段就是指硼与O2生成B2O2的阶段，主导反应为

Yeh和Kuo引入Damkoler数的概念，区分表面反应主导控制机理［5］。本文硼粒子直径较小，符合动力学/扩散控制反应速率模型，表面反应速率同时受到扩散过程和化学反应动力学的影响［12］。

关于硼粒子与环境的换热，认为:

式中 mp、cp、Ap、T∞、h、hfg、εp、σ 、θR分别表示燃烧过程中瞬态硼粒子质量、比热容、表面积和流体温度、对流换热系数、汽化潜热、颗粒黑度、玻耳兹曼常数［11］。

2 计算模型

2.1 物理模型

由于固冲发动机补燃室内燃烧流动非常复杂，包含极度的掺混作用，无法获得单纯的扩散作用对燃烧流动的影响，进一步得出扩散燃烧模型。为简化补燃室内燃烧过程，突出其湍流扩散燃烧本质，通过营造相对简单的燃烧环境，使燃气和空气在补燃室内均匀稳定地发生扩散燃烧，建立均匀分布在燃气与空气的扩散接触面上的三维结构火焰，这种嵌入湍流而掺混不剧烈流场内的火焰是理想的湍流扩散燃烧火焰［13］。为实现这一目的，设计燃气中心进气，空气环周平行进气的冲压发动机实验装置，见图1。

图1 实验装置简图Fig.1 Experimental device sketch

实验装置主要由燃气发生器、补燃室、冲压喷管等组成。文献［14］提出了以实现湍流扩散燃烧为目的的实验方案:雍塞式燃气发生器产生的含硼富燃燃气通过激波减速后，经中心方形管进入补燃室;空气通过集气舱后，经泡沫金属板整流，沿燃气喷口环周进入补燃室;而后，空气与燃气平行均匀流动，并实现扩散燃烧。补燃室采用“方形”三维结构，燃烧室横截面为正方形，边长100 mm，长为627 mm;冲压喷管采用轴对称形式，喉径50 mm;燃气喷嘴为方形结构，喷口尺寸14.2 mm×14.2 mm，伸入补燃室40 mm。燃烧室与喷管简化物理模型见图2。

图2 物理模型示意图Fig.2 Physical model diagram

2.2 湍流输运方程

湍流输运方程包括欧拉坐标系下的三维气相控制方程和k-ε RNG湍流模型方程，具体参数及意义详见文献［8］。

2.3 燃烧模型

为分析不同燃烧模型的影响，气相燃烧分别采用有限速率/涡耗散模型、涡耗散模型、有限速率模型。而颗粒相燃烧采用King点火模型和L-W燃烧模型。

2.4 颗粒轨道模型

两相流燃烧模拟采用颗粒轨道模型［8］。

2.5 燃烧效率表征

燃烧效率反映了燃烧装置实现能量转换的能力，以燃烧生成热占燃料热值比率表征燃烧效率符合燃烧本质。

任意截面硼粒子燃烧效率ηp:

任意截面总燃烧效率ηt:

2.6 边界条件

进口边界条件:

(1)空气入口采用质量入口条件:空气流量0.5 kg/s，空气总温 573 K，空气压强 0.25 MPa;

(2)燃气入口亦采用质量入口条件:一次燃气流量0.05 kg/s，一次燃气总温 2 200 K，一次燃气压强0.5 MPa。空燃比为 10。

出口边界条件采用压力出口边界:出口压力为环境压强。壁面条件采用无滑移标准壁面条件。

燃气组分由燃气发生器热力计算获得，燃气组分:气相组分68%(其中，CO气相质量分数47%，H2气相质量分数10%，CO2气相质量分数10%，H2O气相质量分数1%，其余为N2);硼凝相组分32%，颗粒直径为10 μm，表面液态氧化层厚度为1 μm。

3 数值模拟结果及分析

使用上述模型，采用Fluent流场计算软件，对物理模型进行数值模拟。其中，硼粒子点火模型使用用户自定义函数(UDF)。考虑补燃室结构的对称性，选取一半构筑计算区域，计算网格数250 000;为提高模拟结果的精度，燃气出口与补燃室轴线处，采用局部加密处理。

3.1 不同燃烧模型模拟结果对比

图3是采用不同燃烧模型计算的不同截面处温度平均值和组分浓度平均值沿轴向分布。1号线、2号线、3号线分别对应于湍流燃烧为有限速率/涡耗散模型(模型1)、涡耗散模型(模型2)、有限速率模型(模型3)。图4为不同截面处硼燃烧完成效率和总燃烧效率。由图3(a)和3(b)可见，3种模型计算的温度沿轴向均呈上升规律，O2浓度则呈下降规律，表明了补燃室内燃烧的发展过程。在补燃室前段，3种模型计算的温度相差不大，3种模型模拟的补燃室前段燃烧剧烈程度相似，说明3种模型描述补燃室前段的燃烧过程时区别不大。在补燃室后段，温度上升由快到慢，依次为模型3、模型2、模型1。说明3种模型模拟的补燃室后段燃烧剧烈程度由大到小依次为模型3、模型2、模型1。图4表征的总燃烧效率分布规律与温度分布规律一致，补燃室后段硼燃烧效率也表现出相同的规律。这是由于硼在补燃室后段燃烧放热多，其燃烧起主导作用，不同模型对硼燃烧过程的描述决定了其燃烧的剧烈程度。

图3 不同燃烧模型模拟结果对比Fig.3 Simulation results comparison of different combustion models

图3(c)是CO浓度分布图。在3种模型中，CO浓度在补燃室前段沿轴向呈下降趋势，至中部基本消耗完毕。说明气相组分在补燃室前段进行燃烧，到补燃室后段时，气相组分基本燃烧完全。补燃室前段，模型2的CO浓度下降梯度最大，模型1次之，模型3最小。说明在模型2中，气相组分的燃烧更早开始，且反应更为剧烈，而模型1、模型3气相组分燃烧相对开始较晚。形成上述规律的原因是模型2化学反应速率由大涡混合时间尺度k/ε决定，反应起始不可控，气相反应速度快。

图3(d)是B2O2浓度分布图。3种模型计算的B2O2浓度在补燃室前段缓慢上升，而在补燃室后段迅速上升。硼粒子燃烧主要发生在距离头部超过300 mm的补燃室后段。这是由硼粒子King点火模型的机理决定的，硼粒子燃烧存在点火过程，需要一定的点火条件。在补燃室前段，模型1表征的B2O2浓度上升最快，模型3次之，模型2最慢。补燃室后段，模型2、模型3 B2O2浓度上升差别不大，但均比模型1快。图4表征的硼燃烧效率变化也呈现相同规律。说明采用模型1时，硼在补燃室中段即发生较强的点火过程，其点火反应迅速，但在后段的燃烧发展较慢;采用模型2时，其中段点火较晚，但后段燃烧较强;采用模型3时，其点火反应速度适中，燃烧发展也比较强烈。不同模型描述硼燃烧时，硼燃烧发展的剧烈程度与其点火过程的快慢程度存在着较大差异，且均与气相燃烧规律存在不同。

图4 燃烧效率对比Fig.4 Combustion efficiency comparison

呈现以上规律的原因是不同模型定义化学反应速率的机理不仅对气相反应影响不同，对硼粒子King点火模型和L-W燃烧模型也存在不同影响。不同模型可影响气相组分的反应，同时控制硼粒子的反应。对于硼的点火过程，模型2单独使用湍流混合时间尺度k/ε决定反应速率，使氧化层外表面氧组分迅速反应，但却使氧浓度迅速降低，影响了氧通过氧化层的扩散过程，使得氧化层内侧硼反应较慢，影响了整个点火过程的进行;模型3单独使用Arrhenius反应速率，难以描述湍流脉动，影响氧化层外表面的蒸发反应，减缓点火过程;模型1使用Arrhenius反应速率和涡耗散反应速率对反应进行混合控制，使得氧化层外侧反应速度适中，保证了氧在氧化层中扩散作用的进行，从而在整体上促进了点火反应的进行。对于硼的燃烧过程，LW模型采用化学动力学过程和扩散过程共同控制反应速率，模型2使用涡耗散速率计算，表面反应剧烈，各组分分布梯度大，提高了硼粒子表面扩散程度，加快了燃烧反应的进行。模型3使用Arrhenius反应速率，对于固相硼与气相燃烧的描述与单纯的气相反应不同，单一控制方式加强了硼的表面燃烧。模型1使用Arrhenius反应速率和涡耗散反应速率对反应进行混合控制，反而降低了L-W燃烧模型中的燃烧反应速率，使硼燃烧发展较慢。

3.2 不同燃烧模型流场对比

图5～图8分别是采用不同燃烧模型计算的流线、温度等值线和O2、B2O2浓度等值线的分布。自上至下对应于湍流燃烧分别使用有限速率/涡耗散模型(模型1)、涡耗散模型(模型2)、有限速率模型(模型3)。图5中的3种流线分布趋势一致，燃气和空气流动平行程度高，说明设计的物理模型实现了平行进气的目标。由图6可知，3种模型在补燃室前段，均形成了由内而外温度先升高后降低的纺锤型气相燃烧火焰，在补燃室后段轴线附近，形成了由硼粒子燃烧而产生的核心高温区。

图5 不同燃烧模型流线图Fig.5 Flow lines of different combustion models

图6 不同燃烧模型温度场分布Fig.6 Temperature field comparison of different combustion models

图7 不同燃烧模型O2浓度场分布Fig.7 O2concentration distribution comparison of different combustion models

图8 不同燃烧模型B2O2浓度场分布Fig.8 B2O2concentration distribution comparison of different combustion models

图7反映了O2组分的变化规律。由于O2从燃气喷口外围喷入，该处O2浓度最高，随着流动与扩散的进行，O2得到消耗，其浓度由前向后，由四周向轴线逐渐降低，在补燃室轴线附近达到最低值。由于结构对称性，O2浓度分布沿轴对称。由于扩散作用，O2在火焰面发生化学反应，其浓度下降梯度大，消耗剧烈，对应着O2浓度等值线与温度等值线密集的环向剪切层。补燃室前段对应温度为2 800 K的6号温度等值线反映了气相火焰的长度。可知，模型2气相火焰最长，模型1次之，模型3气相火焰最短，验证了3种模型气相组分燃烧难易程度的规律。

图8反映了硼的燃烧规律。B2O2的浓度变化反映了硼的含量变化，B2O2浓度等值线的疏密反映了硼的消耗强度。由图8可知，3种模型所得B2O2浓度等值线呈环形结构，与补燃室后段，O2浓度等值线对应，反映了硼燃烧的环形火焰面。此处O2浓度等值线分布相对稀疏，说明硼的燃烧不如气相组分燃烧剧烈，其燃烧效率相对偏低。其规则的火焰形状，便于实验验证3种模型对硼粒子点火和燃烧的影响规律。

表征B2O2最低浓度的1号环形浓度等值线的出现位置，与对应温度为3 200 K的7号温度等值线的出现位置对应，反映了硼开始点火的位置。模型1最靠前，距头部0.22 m;模型3其次，距头部0.26 m;模型2最靠后，距头部0.32 m。模型1中硼粒子点火反应提前于模型2、模型3，其补燃室前段硼的燃烧效率高。模型1使用的反应机理与硼粒子King点火模型发生了耦合，使硼粒子点火反应更快开始。而模型1 B2O2环形浓度等值线分布较模型2、模型3稀疏，说明模型1硼粒子点火后燃烧，相对模型2、模型3发展较慢，其补燃室后段硼的燃烧效率低。其原因是模型1采用混合方式控制化学反应速度，为反应提供了动力学开关，使L-W模型的表面反应在火焰稳定后开始，延缓了表面反应时间，从而控制了点火后燃烧的发展。在模型2和模型3的计算结果中，冲压喷管处出现B2O2浓度等值线密集区和温度等值线密集区，说明使用这2种模型时，喷管内硼燃烧较剧烈。冲压喷管内速度高，温度下降快，硼粒子滞留时间短，不利于硼的燃烧。模型2、模型3使用单一方法控制反应的机理较难约束L-W模型硼燃烧的发展。

4 结论

(1)在设计的平行进气物理模型中，燃气和空气沿各自流道呈平行流动，燃烧发生在各种等值线分布密集的环向剪切层内，燃烧火焰沿中心轴呈对称结构，气相组分的燃烧主要发生在补燃室前段，硼的燃烧发生于补燃室后段。

(2)涡耗散模型化学反应速率由大涡混合时间尺度k/ε决定，反应起始不可控，其气相燃烧速度最为剧烈。有限速率模型反应速率根据Arrhenius公式确定，难以描述湍流脉动对反应的影响，其化学反应速率较低，其气相燃烧最为缓慢。有限速率/涡耗散模型使用混合方法控制反应速度，其气相燃烧剧烈程度居于另二者之间。

(3)采用不同模型描述扩散燃烧时，硼粒子的点火和燃烧呈现不同规律，这种规律与气相燃烧规律不尽不同。使用有限速率/涡耗散模型时硼粒子点火最快，其混合控制反应速度的机理，使得氧化层外侧反应速度适中，保证了氧在氧化层中扩散作用的进行，从而在整体上促进了点火反应的进行;对于硼与氧直接发生表面反应的燃烧过程，混合控制方式降低了硼的表面反应速度，控制了硼燃烧过程的发展，使得其点火后燃烧最缓慢，燃烧效率低。使用有限速率模型硼粒子点火速度居于另二者之间，其原因是Arrhenius反应速率难以描述湍流脉动，影响氧化层外表面的蒸发反应，减缓点火过程;单独使用Arrhenius反应速率描述硼的燃烧，有利于提高L-W模型反应速度，硼燃烧较有限速率/涡耗散模型剧烈。使用涡耗散模型硼粒子点火最慢，其原因是氧化层外侧反应过于剧烈，阻碍了氧在氧化层内的扩散，反而降低了硼的点火速率;描述硼的燃烧过程时，增强了硼粒子表面扩散作用，促进L-W模型的燃烧剧烈程度。

以上规律由数值模拟获得，需要进一步通过实验进行验证研究。

［1］张炜，朱慧.冲压发动机发展现状和关键技术［J］.固体火箭技术，1998，21(3):24-30.

［2］Macek A，Semple J.Combustion of boron particles at atmospheric pressure［J］.Combustion Science and Technology，1969，1(3):181-191.

［3］King M K.Boron particle ignition in hot gas streams［J］.Combustion Science and Technology，1974，8(4):243-255.

［4］King M K.Ignition and Combustion of boron particles and clouds［J］.Journal of Spacecraft and Rockets，1982，19(4):294-306.

［5］Yeh C L，Kuo K K.Ignition and combustion of boron particles［J］.Progress in Energy and Combustion Science，1996，22(6):511-541.

［6］Makion A，Law C K.A simplified model for the combustion of uncoated boron particles［J］.Combustion Science and Technology，1988，61(4-6):155-168.

［7］Li S C，Williams F A.Ignition and combustion of boron particles combustion of boron-based solid propellant and solid fuels［M］.Boca Raton，FL，CRC Press，Inc.，1993:248-271.

［8］胡建新.含硼推进剂固体火箭冲压发动机补燃室工作过程研究［D］.长沙:国防科技大学，2006.

［9］胡建新，夏志勋，罗振兵，等.固体火箭冲压发动机补燃室内硼颗粒点火计算研究［J].含能材料，2004，12(6):342-345.

［10］敖文，杨卫娟，韩志江，等.硼颗粒点火燃烧模型研究进展［J］.固体火箭技术，2012，35(3):361-366.

［11］Fluent Inc.Fluent user＇s guide［M］.http://ceclub.cn .2010.

［12］King M K.A review of studies of boron ignition and combustion phenomena at atlantic research corporation over the past cecade［C］//Kuo K K.Combustion of boron-based solid propellants and solid fuels，Boca:CRC press，1993:1-86.

［13］周力行.多相湍流反应流体力学［M].北京:国防工业出版社，2002.

［14］Yang A S，Kuo K K，Hsieh W H.Supercritical evaporation and combustion of liquid oxygen in an axisymmetric configuration［C］//Kuo K K.Recent advances in spray combustion:spray atomization and drop burning phenomena volumeⅠInc.，1996:439-480.