氣氧/氣甲烷針?biāo)▏娮⑵骺仔驮O(shè)計(jì)仿真研究

常一冰，鄒建軍，李清廉

（1.空軍工程大學(xué)航空機(jī)務(wù)士官學(xué)校，信陽(yáng)，464000；2.國(guó)防科技大學(xué)高超聲速?zèng)_壓發(fā)動(dòng)機(jī)技術(shù)重點(diǎn)實(shí)驗(yàn)室，長(zhǎng)沙，410073）

0 引 言

近年來(lái)，隨著商業(yè)航天的興起，國(guó)內(nèi)外航天發(fā)射公司為降低發(fā)射成本，大力拓展發(fā)射裝備的可重復(fù)使用性[1～3]，同時(shí)積極應(yīng)用綠色無(wú)毒、地外天體探索可原位資源利用的液氧/甲烷推進(jìn)劑[4]。在此背景下，具備深度節(jié)流能力、內(nèi)在燃燒穩(wěn)定性且結(jié)構(gòu)簡(jiǎn)單的針?biāo)ㄊ桨l(fā)動(dòng)機(jī)在星球著陸、火箭回收以及軌道控制等場(chǎng)合中顯示出巨大潛力[5]，尤以美國(guó)太空探索技術(shù)公司開(kāi)發(fā)的Merlin 和Raptor 發(fā)動(dòng)機(jī)為典型代表[6]。

然而針?biāo)ㄊ桨l(fā)動(dòng)機(jī)所采用的針?biāo)▏娮⑵髟谌紵鲌?chǎng)中面臨嚴(yán)峻的熱環(huán)境，工作過(guò)程中容易過(guò)熱燒毀[7,8]。目前針對(duì)這一問(wèn)題，提出了3 條解決路徑：a）主動(dòng)冷卻，即在針?biāo)^部開(kāi)孔引出少量推進(jìn)劑并撞擊形成液膜保護(hù)[9]，這種方法簡(jiǎn)單有效，但引出部分推進(jìn)劑用于冷卻會(huì)造成一定的性能損失[10]，同時(shí)對(duì)氧中心式針?biāo)ㄒ鲅趸瘎├鋮s會(huì)造成噴注器材料腐蝕失效[11]；b）采用耐高溫、高導(dǎo)熱性的材料制造噴注器頭部[12,13]，這種方法只是暫時(shí)滿足了工程需要，無(wú)法解釋噴注器頭部受到燒蝕的原因；c）對(duì)噴注器進(jìn)行結(jié)構(gòu)優(yōu)化設(shè)計(jì)[14]，例如在噴孔上沿設(shè)置斜面使軸/徑向推進(jìn)劑撞擊區(qū)域遠(yuǎn)離針?biāo)^以避免其遭受燒蝕，但斜面的存在弱化了推進(jìn)劑的霧化混合[11]。

綜上，本文將以某型500 N 氣氧/氣甲烷針?biāo)òl(fā)動(dòng)機(jī)為研究對(duì)象，在不犧牲燃燒效率的前提下，從噴孔結(jié)構(gòu)優(yōu)化的角度通過(guò)三維數(shù)值仿真考察不同尺寸的圓形噴孔對(duì)針?biāo)▏娮⑵黝^部傳熱環(huán)境和燃燒特性的影響，為針?biāo)▏娮⑵鞣罒嵩O(shè)計(jì)提供有益指導(dǎo)。

1 計(jì)算模型與方法

1.1 500 N 針?biāo)òl(fā)動(dòng)機(jī)模型與孔型變量

500 N 氣氧/氣甲烷針?biāo)òl(fā)動(dòng)機(jī)的剖面結(jié)構(gòu)如圖1所示，其燃燒室設(shè)計(jì)室壓為1.05 MPa，總流量為152 g/s，混合比為3.2，主要幾何參數(shù)列于表1。針?biāo)▏娮⑵髟O(shè)計(jì)參考文獻(xiàn)[15]，采用雙排交錯(cuò)排布的圓形噴孔，每排均勻分布12 個(gè)，沿軸向第1、第2 排噴孔直徑分別記為Df、Ds。

圖1 500N 氣氧/氣甲烷針?biāo)òl(fā)動(dòng)機(jī)結(jié)構(gòu)Fig.1 Schematic of the 500N GOX/GCH4 Pintle Engine

表1 針?biāo)òl(fā)動(dòng)機(jī)主要結(jié)構(gòu)參數(shù)Tab.1 Major Geometrical Parameters of the Pintle Engine

在針?biāo)▏娮⑵鲊娍自O(shè)計(jì)中，通常以阻塞因子BF和總動(dòng)量比TMR作為控制參數(shù)[16]：

式中N為單排噴孔個(gè)數(shù)；m˙為流量；u為流速；r和a分別代表徑向和軸向。可見(jiàn)BF和TMR分別限制了噴孔的周向尺寸和噴注面積。由500 N 發(fā)動(dòng)機(jī)的設(shè)計(jì)工況和表1 幾何參數(shù)可得BF=0.59，TMR=0.66。

為使所得結(jié)論具有一定通用性，這里取無(wú)量綱孔型參數(shù)孔徑比Df/Ds作為研究變量。同時(shí)為單純考察噴孔尺寸對(duì)噴注器傳熱環(huán)境和燃燒特性的影響，本文在變量取值時(shí)嚴(yán)格保證TMR不變而保持BF近似不變，在此前提下研究變量Df/Ds取4 個(gè)值進(jìn)行對(duì)比研究，取值從約0.5 增至4，如表2 所示，各算例均按設(shè)計(jì)工況計(jì)算。表中Case1 設(shè)為參考算例。

表2 孔型變量取值Tab.2 Dimension of Injection Orifices Variants

1.2 控制方程

為較準(zhǔn)確求解含有燃燒反應(yīng)的氣流流動(dòng)和傳熱問(wèn)題，本文采用雷諾平均的N-S 方程（RANS）來(lái)描述燃?xì)獾娜S可壓穩(wěn)態(tài)流動(dòng)，應(yīng)用Standardk-ε模型描述湍流?？紤]到計(jì)算成本和精度，本文采用Jones-Lindstedt 6 步反應(yīng)機(jī)理[17]來(lái)描述甲烷和氧氣的燃燒反應(yīng)。對(duì)湍流與燃燒反應(yīng)的相互作用則采用渦耗散概念模型（Eddy-Dissipation Concept，EDC）計(jì)算，該模型可在湍流中考慮詳細(xì)的反應(yīng)機(jī)理，從而獲得相對(duì)精確的計(jì)算結(jié)果[18]。計(jì)算中假設(shè)燃?xì)鉃槔硐霘怏w，服從理想氣體狀態(tài)方程。

對(duì)燃?xì)馀c針?biāo)ü瘫诘膶?duì)流換熱采用耦合方法[19]求解，即流體域和固體域的接觸面設(shè)為耦合面，其上滿足溫度連續(xù)和熱流密度連續(xù)：

式中T為溫度；n為法向；下標(biāo)fluid 和solid 分別表示流體域和固體域。針?biāo)ü瘫趦?nèi)的熱傳導(dǎo)采用Fourier方程描述。詳細(xì)控制方程形式見(jiàn)文獻(xiàn)[20]。

1.3 網(wǎng)格生成與邊界條件

以參考算例為例，考慮對(duì)稱(chēng)性，計(jì)算域取模型發(fā)動(dòng)機(jī)1/12 的一半。采用結(jié)構(gòu)化網(wǎng)格，其中流體域和固體域網(wǎng)格分別建立，然后組裝為一體，如圖2 所示。

圖2 計(jì)算域結(jié)構(gòu)化網(wǎng)格Fig.2 Structural Mesh of the Computational Domain

流固耦合面上網(wǎng)格非一致，數(shù)據(jù)交換通過(guò)插值實(shí)現(xiàn)。近壁面處采用標(biāo)準(zhǔn)壁面函數(shù)，耦合面兩側(cè)適當(dāng)加密，以獲得較高的求解精度。固體域材料為304 不銹鋼（牌號(hào)：06Cr19Ni10），其物理性質(zhì)設(shè)為常量，其中熱導(dǎo)率取18.3 W/(m·K)。表3 給出了參考算例Case1的邊界條件，其他算例計(jì)算域選取和推進(jìn)劑流量依其噴孔個(gè)數(shù)N而定。

表3 參考算例邊界條件Tab.3 Boundary Conditions for the Case1

1.4 數(shù)值方法

對(duì)流項(xiàng)空間離散采用二階迎風(fēng)格式，壓力-速度耦合采用SIMPLE 算法。由于同時(shí)進(jìn)行燃燒和傳熱問(wèn)題的穩(wěn)態(tài)求解難以收斂，本文燃燒傳熱計(jì)算中采用瞬態(tài)求解器，最終計(jì)算結(jié)果為時(shí)間收斂解。

耦合傳熱計(jì)算時(shí)，式（4）可離散為

式中λ為熱傳導(dǎo)系數(shù)；s為距交界面的法向距離；Tint為交界面溫度；Δs為網(wǎng)格尺寸。由式（5）可得流固交界面溫度：

1.5 算法驗(yàn)證

首先考察網(wǎng)格無(wú)關(guān)性。本文主要求解燃?xì)馀c針?biāo)▏娮⑵魍獗诿娴膶?duì)流換熱，因此網(wǎng)格無(wú)關(guān)性驗(yàn)證僅改變流體域網(wǎng)格數(shù)量。低、中、高密度網(wǎng)格流體域的六面體單元數(shù)依次為108 264、227 662、485 140。

圖3 為采樣線及待分析對(duì)稱(chēng)面的位置，通過(guò)對(duì)圖3所示的采樣線上溫度沿軸向的變化情況進(jìn)行對(duì)比（以氧氣環(huán)縫所在位置為原點(diǎn)），如圖4 所示，可見(jiàn)中、高密度網(wǎng)格的結(jié)果吻合較好，而低密度網(wǎng)格與前兩者相差較大。

圖3 采樣線及待分析對(duì)稱(chēng)面的位置Fig.3 Position of the Sample Line and the Symmetries Inspected

圖4 不同網(wǎng)格密度下的采樣線溫度分布Fig.4 Temperature on the Sample Line for Different Grid Density

其次，對(duì)本文數(shù)值模型進(jìn)行熱試試驗(yàn)驗(yàn)證。低室壓條件下（pc=0.25 MPa）針?biāo)^溫度Ttip試驗(yàn)值與仿真值對(duì)比結(jié)果如表4 所示。

表4 p c=0.25MPa 時(shí)不同混合比下針?biāo)^溫度試驗(yàn)仿真對(duì)比Tab.4 Experimental and Simulational Value of the Pintle Tip Temperature Tt ip with O/F Varying underp c=0.25MPa

定性來(lái)看，隨混合比O/F 增大，仿真值與試驗(yàn)值同步先減小后增大，變化趨勢(shì)一致性較好；定量來(lái)看，試驗(yàn)值仿真值之差與試驗(yàn)值之比最大為22.97%，最小為1.86%。而對(duì)于設(shè)計(jì)工況（O/F=3.2，理論燃?xì)鉁囟?459 K），限于量程和測(cè)點(diǎn)數(shù)量，測(cè)溫?zé)犭娕家褵o(wú)法給出有效數(shù)據(jù)，僅能從試驗(yàn)后針?biāo)^的燒蝕形貌（如圖5所示）進(jìn)行估計(jì)。試驗(yàn)所用的針?biāo)▏娮⑵鞑牧蠟?04不銹鋼，其熔點(diǎn)為1671～1727 K，結(jié)合圖5 認(rèn)為仿真給出的設(shè)計(jì)工況下針?biāo)^溫度值1573 K 符合實(shí)際。

圖5 設(shè)計(jì)工況熱試5s 后針?biāo)^燒蝕形貌Fig.5 Photo of the Pintle Head Ablated after 5s Hot Test

綜上，本文所采用的數(shù)值模型用于初步定性分析是可接受的。

2 計(jì)算結(jié)果與分析

不同孔徑比下圓孔針?biāo)ǖ牟蓸泳€溫度對(duì)比結(jié)果（以氧氣環(huán)縫位置為原點(diǎn)）如圖6 所示，4 個(gè)算例的變化趨勢(shì)都是溫度沿軸向逐漸升高，針?biāo)^尖部位置達(dá)到最高。

圖6 不同孔徑比下的采樣線上溫度對(duì)比Fig.6 Temperature Distribution on the Sample Line for Different Df/Ds

圖7 不同孔徑比下對(duì)稱(chēng)面B 上的流場(chǎng)結(jié)構(gòu)Fig.7 Flow Field in the Symmetry B for Different Df/Ds

圖8 為不同孔徑比下對(duì)稱(chēng)面B 上針?biāo)^部區(qū)域氧氣質(zhì)量分?jǐn)?shù)，由圖8 可知，這一小回流區(qū)是氧氣氣膜穿透甲烷射流進(jìn)而貼覆針?biāo)ū诿媪髦玲標(biāo)^下方形成，且氧氣質(zhì)量分?jǐn)?shù)不低于0.7。這樣該區(qū)域不會(huì)充分燃燒，所以其本身溫度相對(duì)較低，同時(shí)小回流區(qū)的存在也避免了針?biāo)^受到燃燒室中心大回流區(qū)高溫燃?xì)饣亓鞯闹苯訜g。而Case2 的氧氣氣膜在與甲烷撞擊后消耗較多，針?biāo)^下方區(qū)域氧氣質(zhì)量分?jǐn)?shù)低于0.7，針?biāo)^也就無(wú)法得到保護(hù)。因此針?biāo)^下方小回流區(qū)的存在是有利于避免針?biāo)^過(guò)熱燒毀的。

圖8 不同孔徑比下對(duì)稱(chēng)面B 上針?biāo)^部區(qū)域氧氣質(zhì)量分?jǐn)?shù)Fig.8 Distribution of Oxygen Mass Fraction in the Region Around the Pintle Head in the Symmetry B for Different Df/Ds

圖7 表明燃燒室室壁附近的溫度高達(dá)3000 K 以上，這是在氣氣燃燒以及絕熱條件下的計(jì)算結(jié)果，實(shí)際應(yīng)用時(shí)通常會(huì)利用液態(tài)甲烷的高比熱進(jìn)行再生冷卻以及采用耐高溫材料制造推力室，這方面已有諸多研究證明是可行的，本文不再贅述。

為評(píng)估孔徑比對(duì)發(fā)動(dòng)機(jī)燃燒特性的影響，引入特征速度燃燒效率：

按照孔徑比增大的順序，Case2、Case3、Case1、Case4 的燃燒效率依次為96.38%、96.33%、96.44%、96.55%，可見(jiàn)孔徑比變化對(duì)燃燒效率的影響非常微小。

3 結(jié) 論

本文基于500 N 氣氧/氣甲烷針?biāo)òl(fā)動(dòng)機(jī)在阻塞因子BF=0.59、總動(dòng)量比TMR=0.66 和其他幾何條件不變的前提下，就圓形噴孔孔徑比對(duì)針?biāo)▏娮⑵黝^部傳熱環(huán)境和燃燒特性的影響進(jìn)行仿真研究，得出以下結(jié)論：

a）當(dāng)0.55≤Df/Ds≤4 時(shí)，隨孔徑比Df/Ds增大，針?biāo)▏娮⑵鞅砻孀罡邷囟戎饾u減小，而發(fā)動(dòng)機(jī)燃燒效率變化微小且均在96%以上。

b）三維仿真給出的流場(chǎng)特征表明除燃燒室中2 個(gè)大型回流區(qū)外，當(dāng)孔徑比Df/Ds≥1 時(shí)針?biāo)^下方還存在一小回流區(qū)，且該區(qū)域軸向推進(jìn)劑質(zhì)量分?jǐn)?shù)較大，對(duì)保護(hù)針?biāo)▏娮⑵黝^部防止過(guò)熱燒毀具有重要作用。

但從試驗(yàn)結(jié)果來(lái)看，大孔徑比下針?biāo)▏娮⑵黝^部仍然存在較嚴(yán)重的燒蝕，因此僅對(duì)孔型結(jié)構(gòu)進(jìn)行優(yōu)化是不夠的，后續(xù)針對(duì)針?biāo)▏娮⑵鞯姆罒g設(shè)計(jì)將進(jìn)行更深入的研究。