針栓式液氧/煤油發(fā)動機燃燒數(shù)值仿真

俞南嘉，鮑啟林，張 洋，戴 健

(北京航空航天大學 宇航學院，北京 100191)

0 引言

目前，液體火箭發(fā)動機是航天運輸及空間飛行器的主要動力裝置，液氧/煤油火箭發(fā)動機因其制造成本低、發(fā)動機性能高、環(huán)保無污染、推力可控以及安全系數(shù)高等特點，在航天運輸領域得到廣泛應用[1]；同時，液氧/煤油火箭發(fā)動機也是未來可重復使用發(fā)動機發(fā)展的重要方向，可重復使用火箭發(fā)動機意味著火箭的子級可以回收并且重復使用，其實現(xiàn)的重要途徑就是研發(fā)推力可調的火箭發(fā)動機[2]。

針栓式噴注器作為變推力發(fā)動機應用的典型，被廣泛應用到火箭發(fā)動機中。最著名的應用就是美國TRW公司登月艙下降段發(fā)動機，該發(fā)動機先后將12名宇航員送上月球[3]；1995年美國又研制了世界上最大的針栓發(fā)動機TR-106，推力達到了290噸級[4]；美國SpaceX公司于2015年12月完成了獵鷹9號火箭的首次陸上垂直降落回收，其中，獵鷹9號火箭的梅林1D發(fā)動機即采用針栓式噴注器[5]，推力可實現(xiàn)50 %～100 %的調節(jié)；2016年，Min Son等人提出了針栓式噴注器最新的理論設計方法[6]；我國于上世紀70年代開始針栓式噴注器的研究，1992年研制成功機械定位的變推力發(fā)動機，該發(fā)動機采用自燃推進劑N2O4/UDMH[6]；2013年嫦娥三號探測器采用雙調針栓式噴注器實現(xiàn)了月球軟著陸，該發(fā)動機依然采用常溫自燃推進劑[7]。但是，目前我國對于低溫非自燃推進劑的針栓式發(fā)動機的研究很少，開展針栓式液氧/煤油發(fā)動機的研究很有意義。

1 計算模型

1.1 物理模型

針栓式變推力液體火箭發(fā)動機的原理是通過調節(jié)噴注器中的環(huán)形套筒來調節(jié)兩個同軸環(huán)縫的噴注通道的流通面積，實現(xiàn)對燃料和氧化劑噴注壓降、噴注液膜厚度等參數(shù)的調節(jié)，其中，環(huán)形調節(jié)套筒稱為針栓式噴注器的針閥，兩個環(huán)形噴注通道稱為氧化劑噴嘴和燃料噴嘴[8]。針栓式噴注器原理圖如圖 1所示。

圖 1 針栓式噴注器原理圖

1.2 數(shù)學模型

液體火箭發(fā)動機的噴注燃燒一般是氣液兩相的噴注燃燒，推進劑通過氧化劑噴嘴和燃料噴嘴以液體射流或者液膜的形態(tài)噴入燃燒室，氧化劑和燃料之間強烈的相互作用，使得液體推進劑發(fā)生碰撞、霧化、混合等過程，并在對流、輻射等傳熱過程中蒸發(fā)和燃燒。因此，數(shù)學模型涉及到控制方程、湍流模型、燃燒模型等[9]。

1.2.1 控制方程

液體火箭發(fā)動機噴霧燃燒的控制方程包含氣相控制方程和液相控制方程，氣相控制方程一般用歐拉坐標系中的N-S方程來描述，液相控制方程則是由拉格朗日坐標系中N-S方程對單個液滴進行描述，氣液兩相之間的耦合采用氣液兩相的耦合源項來考慮[8]。

對于有N個組元的混合物來說，控制方程由連續(xù)方程、動量方程、能量方程、組元方程以及狀態(tài)方程組成封閉方程組，其中，連續(xù)方程為

(1)

動量方程為

(2)

能量方程為

(3)

組元方程為

(4)

狀態(tài)方程為

p=ρRT

式中：ρ為混合物密度；Sm為質量守恒方程中的源項，代表單位體積中所有液滴蒸發(fā)導致的氣相質量增長率；Sdi,Sdm為動量方程中相互作用的源項，分別為單位體積的液滴反作用力使氣相的動量的增長率和單位體積液滴蒸發(fā)帶給氣相的動量；Se為能量守恒方程中的相互作用源項；Sg為組分方程中的源項；h為混合物靜焓，τij為應力張量，μ為動力粘性系數(shù)，wi為組分i的化學反應速度，R為通用氣體常數(shù)。

1.2.2 湍流模型

湍流模型采用k-ε雙方程模型來計算液體火箭發(fā)動機的燃燒，該模型廣泛應用于各種工程流動和熱交換的數(shù)值模擬中，并得到了較好工程驗證，k-ε模型方程為：

(5)

(6)

式中：C1ε和C2ε為經驗常數(shù)，C1ε=1.44，C2ε=1.92；σk和σε分別是與湍流動能k和耗散率ε對應的普讓特數(shù)；模型常數(shù)取值：Cμ=0.09，σk=1.0，σε=1.3；μt為湍流粘性；Gk是由于速度梯度引起的湍流動能k的產生項。

1.2.3 組分輸運模型

燃燒的組分輸運模型采用有限速率-渦耗散模型，反應考慮了組分之間的相互擴散過程和推進劑進口的擴散作用；該模型的化學反應速率由大渦混合時間尺度k/ε控制，并對Arrhenius和渦耗散反應速率都進行計算，并取兩者之間化學反應速率較小的一個[9]；液氧和煤油的化學反應機理采用單步總包反應，其化學反應關系式為：

C12H23+17.75O2=12CO2+11.5H2O

由阿倫尼烏斯計算的燃燒反應速率為：

(7)

由渦耗散控制的化學反應速率為：

(8)

(9)

2 設計工況

所設計發(fā)動機的總質量流量為381 g/s，理論比沖為2 806.59 m/s，特征速度為1 776.6 m/s，噴注方案為燃料中心式，即燃料在內側噴注，噴注角度為90o；氧化劑在外側噴注，噴注角度為52.36o，發(fā)動機的參數(shù)如表 1所示。

2.1 網格劃分

由于發(fā)動機計算域結構的對稱性，本文采用1/8模型進行網格的劃分，使用ICEM軟件劃分計算域網格。由于燃燒室內存在劇烈的化學反應，并伴隨著復雜的流動過程，因此為更準確地捕捉流場信息，對噴嘴、燃燒室頭部和燃燒室壁面進行網格加密，計算域網格數(shù)量為118萬，圖 2為燃燒室頭部的計算域網格。

圖2 燃燒室頭部計算域網格

2.2 邊界條件設置

1)入口邊界條件：對于液體火箭發(fā)動機的數(shù)值仿真，一般采用離散相模型作為入口邊界。本文設置液滴每隔20步化學反應噴注一次，同時考慮顆粒對連續(xù)相的影響，顆粒軌道的迭代步數(shù)設置為1 000，離散物模型選擇壓力梯度模型、隨機碰撞模型以及裂解模型；液氧和煤油顆粒噴射源的噴注類型選擇SURFACE，顆粒形式為DROPLET，同時分別設置液滴的入口溫度、入口速度、質量流量、液滴尺寸等，其中，煤油和液氧的液滴尺寸均設置為0.4 mm。

2)出口邊界條件：發(fā)動機噴管出口一般為超音速的氣流，各參數(shù)計算采用外插公式得到。出口邊界條件設置湍流粘度和水利直徑，出口壓力為環(huán)境壓力。

3)壁面邊界條件：壁面邊界采用絕熱無滑移的壁面條件，即u=w=v=0；并且壓力梯度導數(shù)、組分質量分數(shù)梯度導數(shù)為零。

4)對稱面邊界：流場是對稱分布的，因此對稱面采用對稱面邊界條件，對稱面處參數(shù)的法向導數(shù)為零[10-11]。

3 計算結果

仿真計算采用Fluent軟件完成，計算獲得了發(fā)動機的燃燒流場分布。圖 3為發(fā)動機燃燒室內頭部附近流線分布，根據(jù)圖可知，在燃燒室頭部以及針栓頭部附近存在兩塊回流區(qū)，即外側回流區(qū)和中心回流區(qū)，這是針栓式噴注器噴霧混合模式的典型特征，這種形式是實現(xiàn)高效燃燒和保證燃燒穩(wěn)定性的基礎。外側回流區(qū)主要是推進劑向頭部的回流；中心回流區(qū)主要是推進劑向針栓頭部回流，加劇推進劑液滴的二次破碎和摻混燃燒效果，提高燃燒效率。

圖 3 燃燒室內頭部附近流線分布

圖4為動量比1工況下的溫度分布以及Y=100 mm位置截面的溫度分布，根據(jù)軸向的溫度分布可知，燃燒室壁面上游區(qū)域存在由氧化劑形成的液膜低溫區(qū)，這是因為液氧和煤油從噴注通道流出后發(fā)生碰撞，氧化劑獲得徑向速度分量，向上運動碰撞到燃燒室壁面后，分成兩部分，一部分向低壓區(qū)(燃燒室頭部)流動，形成燃燒室頭部回流區(qū)；一部分沿燃燒室壁面向下游流動，形成燃燒室壁面的液膜，該液膜對發(fā)動機壁面起到冷卻的作用，對發(fā)動機冷卻是有利的。流動過程中煤油與氧的不斷燃燒消耗，使得液膜厚度變得越來越小，到燃燒室中下游附近消失；中心回流區(qū)內的主要組分是未參與燃燒的煤油，在主流燃氣的帶動下，未燃燒的煤油不斷向下游流動，使得中心軸線附近出現(xiàn)低溫區(qū)，這個低溫區(qū)有利于降低針栓頭部的溫度，起到保護針栓的作用。

圖 4 溫度分布

圖5為燃燒室中心軸線的溫度曲線，根據(jù)曲線可知，在針栓頭部附近溫度大約為1 445 K，在噴管的喉部附近溫度達到最高，約為2 500 K；圖 6為Y=100 mm，Z=0 mm位置處X方向的湍流強度分布，據(jù)圖可知，在徑向X=8 mm左右位置處湍流強度最大，該位置對應溫度流場中的燃燒高溫區(qū)，該區(qū)域存在劇烈的化學反應，湍流脈動速度大，導致推進劑混合劇烈，可以提高燃燒效率。

圖 5 軸向溫度曲線

圖 6 湍流強度分布

3.1 動量比的影響

動量比是指推進劑的徑向射流動量值與軸向射流動量值的比，即

動量比作為針栓式噴注器的重要性能參數(shù)，代表了噴注器的霧化混合程度，推進劑的霧化混合主要是通過徑向射流液膜與軸向射流液膜碰撞實現(xiàn)，針栓式噴注器可以通過改變動量比的大小，提高或者降低射流的穿透能力，增強或者減弱推進劑的霧化混合程度，從而改變燃燒效率高低同時。因此，本文分別設計動量比為0.45，0.8，1.0，1.2和1.5的工況進行數(shù)值計算，得出如下結果。

圖7為不同動量比的溫度分布，根據(jù)圖可知，隨著動量比增加，燃燒室壁面液膜厚度越來越小，壁面溫度越來越高；動量比為1.5時，燃燒高溫區(qū)已經貼近燃燒室壁面，說明隨著徑向射流動量的增加，燃料和氧化劑的燃燒面向壁面移動，使得燃燒高溫區(qū)靠近壁面。根據(jù)圖 8和圖 9的壓力分布可知，隨著動量比增加燃燒室壓力呈現(xiàn)先增加后減小的趨勢，在動量比為1時，燃燒室壓力最高，這是因為動量比小于1時，隨著動量比增加，推進劑徑向相對射流動量增加，一定程度上提高了煤油和氧的混合程度，有利于燃燒，使的燃燒室壓力有增加的趨勢，但是，動量比大于1時，隨著動量比的增加，燃燒室壓力有減小的趨勢，說明動量比大于1時，徑向射流動量比較大，軸向射流動量相對減小，使得軸向推進劑組元穿透徑向射流的難度增加，導致碰撞、混合的效果減弱，發(fā)動機燃燒效果略有減小，發(fā)動機燃燒室壓力有所降低。圖10為不同動量比的煤油質量分數(shù)分布云圖，根據(jù)圖可知，隨著動量比增加，徑向射流動量增強，使得煤油組分在燃燒室內占據(jù)的區(qū)域增加，氧和煤油摻混的區(qū)域向燃燒室壁面靠近，導致燃燒的高溫區(qū)向燃燒室壁面靠近，提高了發(fā)動機壁面的熱載荷，不利于發(fā)動機燃燒室的熱防護。

圖 7 不同動量比溫度分布

圖 8 不同動量比壓力分布

圖 9 不同動量比軸向壓力分布曲線

圖 10 不同動量比煤油質量分數(shù)分布

3.2 燃燒效率分析

為了評估燃燒室內能量轉換過程中的損失大小，引入燃燒效率的概念，火箭發(fā)動機的燃燒效率是指實際特征速度與理論特征速度的比值[11]，即：

根據(jù)數(shù)值仿真結果計算出不同動量比工況下的燃燒效率如表 2所示。根據(jù)表可知，動量比為1時，發(fā)動機的燃燒效率最高為96.65%。為了分析發(fā)動機燃燒效率和動量比之間的變化關系，探究發(fā)動機具有最佳燃燒性能的工況點，根據(jù)表 2數(shù)據(jù)采用多項式擬合的方法擬合出燃燒效率與動量比的關系曲線如圖 11所示，根據(jù)圖可知，隨著動量比增加，燃燒效率呈現(xiàn)先增加后減小的趨勢，并且動量比為1時燃燒效率最高，說明該狀態(tài)下推進劑混合效果最好，具有最佳的燃燒性能。因此，在設計針栓式噴注器時，為了提高發(fā)動機性能，盡量選擇動量比為1的工況點進行設計。

表 2 不同動量比燃燒效率

圖 11 燃燒效率與動量比關系曲線

4 結論

本文針對液氧/煤油推進劑的針栓式發(fā)動機進行數(shù)值仿真，得出如下結果：

1)針栓式發(fā)動機在燃燒室內形成兩個回流區(qū)，即燃燒室頭部回流區(qū)和針栓頭部回流區(qū)，回流區(qū)的存在有利于發(fā)動機噴注器的熱防護；

2)針栓式發(fā)動機的噴注形式能夠在燃燒室壁面形成液膜，有利于燃燒室冷卻；

3)隨著動量比增加，燃燒該高溫區(qū)向壁面靠近，增加了壁面熱載荷，對發(fā)動機的熱防護是不利的；

4)動量比接近1時，霧化混合效果最好，燃燒效率最高。