采用不同氣相燃燒模型模擬含硼燃氣擴散燃燒過程①

馮喜平，李海波，唐金蘭，李貴珠，廖自繁

(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.