小型氣氧/煤油火箭發(fā)動(dòng)機(jī)試驗(yàn)及內(nèi)流場(chǎng)數(shù)值模擬

郭康康，聶萬(wàn)勝，劉 瑜，蘇凌宇，石天一

（航天工程大學(xué)，北京，101416）

0 引 言

隨著航天技術(shù)的不斷推進(jìn)，人類(lèi)進(jìn)入太空的探索活動(dòng)越加頻繁。航天不斷成為各國(guó)激烈競(jìng)爭(zhēng)的領(lǐng)域，航天發(fā)射、動(dòng)力先行、航天活動(dòng)對(duì)動(dòng)力系統(tǒng)提出了更高的要求。當(dāng)下航天動(dòng)力系統(tǒng)主要是液體火箭發(fā)動(dòng)機(jī)。液氧/煤油火箭發(fā)動(dòng)機(jī)因推力大、無(wú)污染、經(jīng)濟(jì)和可靠性高成為各國(guó)競(jìng)爭(zhēng)的焦點(diǎn)[1]。鑒于此，本文進(jìn)行了小型氣氧/煤油發(fā)動(dòng)機(jī)的熱試車(chē)，并針對(duì)熱試車(chē)工況分別進(jìn)行基于單步化學(xué)反應(yīng)和多步化學(xué)反應(yīng)的數(shù)值模擬。

中國(guó)新一代大推力液氧/煤油火箭發(fā)動(dòng)機(jī)采用了先進(jìn)的高壓補(bǔ)燃循環(huán)系統(tǒng)，煤油和氧氣形成富氧可燃?xì)鈭F(tuán)在補(bǔ)燃室內(nèi)燃燒，隨后，經(jīng)過(guò)燃燒的混合氣進(jìn)入主燃燒室，由于可燃?xì)鈭F(tuán)高度富氧，可以近似為高溫氧氣進(jìn)入主燃燒室[2]。因此本文進(jìn)行的氣氧/煤油火箭發(fā)動(dòng)機(jī)試驗(yàn)及數(shù)值模擬對(duì)大型發(fā)動(dòng)機(jī)的內(nèi)流場(chǎng)特性研究有一定的參考價(jià)值。

1 試驗(yàn)系統(tǒng)及工況

本文采用的試驗(yàn)發(fā)動(dòng)機(jī)如圖1所示。

圖1 試驗(yàn)發(fā)動(dòng)機(jī)及其三維模型示意Fig.1 Experiment Rocket Engine and It's Three-dimension Model

試驗(yàn)發(fā)動(dòng)機(jī)采用氣氧/煤油作為推進(jìn)劑，點(diǎn)火過(guò)程在點(diǎn)火凹腔中完成，利用火花塞點(diǎn)燃點(diǎn)火劑（氫氣和氧氣）來(lái)啟動(dòng)發(fā)動(dòng)機(jī)。燃燒室收斂段與喉部之間有一套氣膜冷卻系統(tǒng)，發(fā)動(dòng)機(jī)啟動(dòng)后，通過(guò)氣膜冷卻系統(tǒng)向燃燒室側(cè)壁噴入N2，收斂段與喉部之間的內(nèi)壁上形成一層薄氣膜來(lái)阻隔熱量。點(diǎn)火盤(pán)上安裝有6個(gè)噴嘴，噴嘴等距布置在直徑為80 mm的圓周上，噴嘴類(lèi)型為同軸離心式，氧氣和煤油壓降為0.8 MPa和0.6 MPa。壓力傳感器安裝在距離點(diǎn)火盤(pán)面板10 cm的側(cè)壁上。發(fā)動(dòng)機(jī)參數(shù)及試驗(yàn)工況如表1所示。

表1 試驗(yàn)發(fā)動(dòng)機(jī)構(gòu)型參數(shù)及試驗(yàn)工況Tab.1 Parameters of Experiment Rocket Engine and Initial Experiment Conditions

2 試驗(yàn)結(jié)果

試驗(yàn)過(guò)程中，燃燒室工作時(shí)間段從點(diǎn)火開(kāi)始到煤油供給結(jié)束，持續(xù)4.1 s。通過(guò)壓力傳感器測(cè)量的壓強(qiáng)如圖2所示。由圖2可知，6.8～8.3 s是燃燒室壓強(qiáng)穩(wěn)定階段，其壓強(qiáng)為2.43 kPa。

圖2 試驗(yàn)壓強(qiáng)值曲線Fig.2 Pressure of Experiment

3 針對(duì)試驗(yàn)工況的數(shù)值模擬

3.1 數(shù)值計(jì)算模型與控制方法

試驗(yàn)發(fā)動(dòng)機(jī)采用同軸離心式噴嘴，氧氣沿軸向進(jìn)入燃燒室，煤油以一定的旋流速度緊貼噴嘴縮進(jìn)段內(nèi)壁噴入燃燒室[3,4]。針對(duì)氣液兩相流，氣相可以用直角坐標(biāo)系下的N-S方程來(lái)描述[5～7]，即：

式（1）分別表示質(zhì)量方程、動(dòng)量方程（包含x，y和z方向）和能量方程。

利用有限體積法對(duì)上述控制方程進(jìn)行離散。針對(duì)N-S方程的強(qiáng)非線性相互耦合性，利用壓力隱式算子分裂算法（Pressure Implicit with Splitting of Operators，PISO）[8,9]實(shí)現(xiàn)壓力和速度的解耦，在離散對(duì)流項(xiàng)和粘性項(xiàng)時(shí)，界面插值分別采用具有二階精度的 Van Leer[10,11]格式和中心差分格式，同時(shí)利用具有二階精度的向后差分格式對(duì)時(shí)間項(xiàng)進(jìn)行離散。

液相通過(guò)顆粒軌道模型進(jìn)行描述；氣液兩相間的質(zhì)量、動(dòng)量和能量交換通過(guò)控制方程源項(xiàng)實(shí)現(xiàn)；通過(guò)給定霧化錐角和液滴直徑分布來(lái)描述液相的霧化過(guò)程；蒸發(fā)過(guò)程由對(duì)流擴(kuò)散來(lái)控制[12～14]。

考慮到計(jì)算成本及所研究問(wèn)題特點(diǎn)，采用的湍流模型為可壓縮模型（Reynolds-Averaged Navier- Stokes，RANS），采用標(biāo)準(zhǔn)為k-ε雙方程[7]。在單步化學(xué)反應(yīng)數(shù)值模擬中，湍流燃燒模型為渦團(tuán)耗散模型（Eddy Dissipation Model，ED）[15,16]，該模型在渦團(tuán)破碎模型（Eddy Break-Up，EBU）的基礎(chǔ)上改進(jìn)得到，認(rèn)為化學(xué)反應(yīng)速率不但受組分輸運(yùn)與湍流脈動(dòng)的影響，而且還與氧化劑、燃料和反應(yīng)產(chǎn)物中濃度值最小的組分有關(guān)；而在多步化學(xué)反應(yīng)模型中，湍流燃燒模型為（Eddy Dissipation Model Concept Model，EDC）渦團(tuán)耗散概念模型[15,16]，EDC模型是在ED模型基礎(chǔ)上考慮湍流流動(dòng)中詳細(xì)的化學(xué)反應(yīng)機(jī)理，假定化學(xué)反應(yīng)進(jìn)行于微小的湍流結(jié)構(gòu)上而建立的。

上述控制方程的離散和求解在開(kāi)源程序平臺(tái)[17]（Open Source Field Operation and Manipulation，OpenFOAM）實(shí)現(xiàn)，OpenFOAM具有良好的開(kāi)源性，利用C++語(yǔ)言把CFD方程組編寫(xiě)為面向?qū)ο蟮膹?qiáng)大底層類(lèi)庫(kù)，通過(guò)調(diào)用特定功能的類(lèi)庫(kù)形成滿足具體應(yīng)用問(wèn)題的求解器，實(shí)現(xiàn)偏微分方程組的數(shù)值求解。

3.2 化學(xué)反應(yīng)模型

本文分別采用單步總包反應(yīng)和多步反應(yīng)對(duì)試驗(yàn)工況進(jìn)行數(shù)值模擬。煤油是一種混合物，含有多種組分，在單步反應(yīng)和9組分14步化學(xué)反應(yīng)中，分別用C12H23和C12H24作為替代燃料。單步反應(yīng)中，化學(xué)反應(yīng)速率計(jì)算如下：

式中 A為指前因子；上標(biāo) Ea為活化能，T為溫度；分別為煤油和氧氣的摩爾濃度，mol/cm3，單步化學(xué)反應(yīng)機(jī)理如表2[18]所示。

表2 單步化學(xué)反應(yīng)機(jī)理Tab.2 Single Step Global Chemical Reaction

多步反應(yīng)參見(jiàn)文獻(xiàn)[19]和文獻(xiàn)[20]，把煤油與氧氣的反應(yīng)簡(jiǎn)化為9組分14步，化學(xué)反應(yīng)速率如下式：

多步化學(xué)反應(yīng)機(jī)理如表3所示[19,20]。

為拓寬基層來(lái)稿渠道，廣泛掌握各地園藝產(chǎn)業(yè)新技術(shù)、新經(jīng)驗(yàn)和新情況，加強(qiáng)與讀者、作者的溝通聯(lián)系，現(xiàn)面向全國(guó)各地縣鄉(xiāng)農(nóng)業(yè)技術(shù)干部征聘本刊通訊員。

表3 多步化學(xué)反應(yīng)機(jī)理Tab.3 nMulti-step Chemical Reactio

續(xù)表3

3.3 網(wǎng)格與邊界條件

計(jì)算網(wǎng)格為非結(jié)構(gòu)網(wǎng)格，因?yàn)榭疾斓氖窃囼?yàn)發(fā)動(dòng)機(jī)的內(nèi)流場(chǎng)特征，為減少計(jì)算量，從噴管段截去一部分。試驗(yàn)網(wǎng)格如圖 3所示，試驗(yàn)發(fā)動(dòng)機(jī)內(nèi)流場(chǎng)邊界條件如表4所示。

圖3 試驗(yàn)發(fā)動(dòng)機(jī)網(wǎng)格示意Fig.3 Grid of Experimental Rocket Engine

表4 發(fā)動(dòng)機(jī)內(nèi)流場(chǎng)邊界條件Tab.4 Boundary Conditions of Internal Flow Field

煤油霧化參數(shù)通過(guò)給定霧化角和霧化粒徑分布的方式給出，其中霧化角為 45°，粒子分布選用Rosin-Rammler模型[5]，平均粒徑為50 μm。監(jiān)測(cè)點(diǎn)是傳感器的安裝位置，位于距離噴注面板10 cm的壁面上。

3.4 數(shù)值計(jì)算結(jié)果分析

針對(duì)單步化學(xué)反應(yīng)模型，利用3套網(wǎng)格進(jìn)行了網(wǎng)格無(wú)關(guān)性的驗(yàn)證，網(wǎng)格量分別為24.7萬(wàn)個(gè)、35.1萬(wàn)個(gè)和 45.3萬(wàn)個(gè)，記為 grid1～3。數(shù)值模擬結(jié)果如圖4所示。當(dāng)網(wǎng)格量大于35.1萬(wàn)個(gè)，數(shù)值計(jì)算對(duì)網(wǎng)格量已經(jīng)不再敏感，因此為了節(jié)省計(jì)算量，選擇grid2作為多步化學(xué)反應(yīng)模型數(shù)值模擬的計(jì)算網(wǎng)格。

分別記單步總包反應(yīng)與9組分14步化學(xué)反應(yīng)模型下的算例為case1和case2，監(jiān)測(cè)點(diǎn)壓強(qiáng)如圖5所示。

圖4 網(wǎng)格無(wú)關(guān)性驗(yàn)證結(jié)果曲線Fig.4 Numerical Simulation of Grid Independence

圖5 單步反應(yīng)與多步反應(yīng)壓強(qiáng)結(jié)果對(duì)比曲線Fig.5 Pressure Results of Single-step Reaction and Multi-step Reactions

單步總包反應(yīng)和多步反應(yīng)壓強(qiáng)分別為2.33 kPa和2.27 kPa，與試驗(yàn)穩(wěn)定室壓相比，誤差分別為 4.1%和6.6%，均在合理誤差范圍內(nèi)，驗(yàn)證了數(shù)值計(jì)算模型的正確性。數(shù)值模擬計(jì)算值小于試驗(yàn)值的原因是：在數(shù)值模擬中，氧氣入口質(zhì)量流量按主氧路來(lái)給定而忽略點(diǎn)火路氧氣流量的影響。試驗(yàn)中，點(diǎn)火路氧氣供給存在響應(yīng)遲滯，因此氧氣質(zhì)量流量偏多是富燃試驗(yàn)工況誤差的主要來(lái)源。case1壓強(qiáng)稍微高于 case2，這是由于單步化學(xué)反應(yīng)模型忽略了中間的吸熱反應(yīng)。由于兩者通過(guò)數(shù)值模擬得到的壓強(qiáng)差別很小，因此在計(jì)算室壓時(shí)可以互相替代。若單純考察燃燒室室壓，選擇單步化學(xué)反應(yīng)模型更加高效。

圖6為case1和case2的溫度分布云圖，兩者的溫度分布存在較大差異。其中，case1的最高溫度約為4550 K，case2的最高溫度約為3530 K。如前所述，由于單步總包化學(xué)反應(yīng)模型忽略了中間吸熱反應(yīng)，化學(xué)反應(yīng)釋放的熱量偏多，導(dǎo)致溫度的數(shù)值計(jì)算結(jié)果較多步化學(xué)反應(yīng)模型偏高，因此多步化學(xué)反應(yīng)模型計(jì)算的溫度值更符合實(shí)際。

圖6 基于單步與多步化學(xué)反應(yīng)對(duì)稱(chēng)面溫度云圖Fig.6 Temperature Contours of Chamber Based on Single-step Reaction and Multi-step Reactions

在試驗(yàn)中，很難準(zhǔn)確地測(cè)量溫度值。因此通過(guò)多步化學(xué)反應(yīng)模型的數(shù)值計(jì)算結(jié)果，可以在一定程度上反應(yīng)內(nèi)流場(chǎng)的溫度分布。由圖6可知，燃燒室圓柱端靠近壁面的區(qū)域溫度較中心高溫區(qū)域低，這是因?yàn)榭拷紵冶诿?C12H24質(zhì)量分?jǐn)?shù)較高，如圖 7所示。C12H24蒸發(fā)吸熱，使得近壁區(qū)溫度降低，這對(duì)燃燒室內(nèi)壁面冷卻有利。由圖6還可看出，高溫區(qū)主要集中在燃燒室收斂段和喉部之間，該區(qū)域熱負(fù)荷較大。為了降低此區(qū)域的熱負(fù)荷，如圖1所示，在燃燒室收斂段和喉部之間安裝有N2膜冷卻系統(tǒng)，通過(guò)氣膜冷卻系統(tǒng)以一定速度向燃燒室噴入N2，使得高熱負(fù)荷區(qū)域形成N2薄膜，起熱防護(hù)作用。在確定氣膜冷卻系統(tǒng)的冷卻區(qū)域時(shí)，可以參考內(nèi)流場(chǎng)溫度云圖，使得發(fā)動(dòng)機(jī)在工作狀態(tài)下，整個(gè)收斂段內(nèi)壁面上的高溫區(qū)域均能形成N2薄模，以此達(dá)到較好的冷卻效果。

圖7 基于多步化學(xué)反應(yīng)對(duì)稱(chēng)面煤油質(zhì)量分?jǐn)?shù)云圖Fig.7 Mass Fractions Contour of C12H24 Based on Multi-step Reactions

圖8 和圖9分別為case1和case2的O2、CO2、H2O質(zhì)量分?jǐn)?shù)云圖分布。對(duì)比圖8和圖9所示，在相同的氧氣質(zhì)量流量和噴射速度下，case2的 O2噴射長(zhǎng)度較case1長(zhǎng)，因?yàn)樵趩尾椒磻?yīng)中，O2一經(jīng)噴入便快速消耗；而在多步反應(yīng)中，氧氣在多個(gè)中間反應(yīng)中逐漸消耗，因此氧氣在燃燒室中的留存時(shí)間較長(zhǎng)，表現(xiàn)為O2的噴射長(zhǎng)度較長(zhǎng)。

圖8 基于單步化學(xué)反應(yīng)模型的組分云圖Fig.8 Mass Fractions Contour of Components Based on Single-step Reaction

續(xù)圖8

圖9 基于9步14組分化學(xué)反應(yīng)模型的組分云圖Fig.9 Mass Fractions Contour of Components Based on Multi-step Reactions

case1和case2的CO2和H2O組分分布差異較大，在case1中，化學(xué)反應(yīng)一步完成，同時(shí)產(chǎn)生CO2和H2O，所以?xún)烧哔|(zhì)量分?jǐn)?shù)云圖分布幾乎相同；而在case2中，CO2和H2O在不同的中間反應(yīng)過(guò)程生成，生成CO2和H2O的中間反應(yīng)過(guò)程在時(shí)間和空間上均不一致。因此，CO2和H2O的質(zhì)量分?jǐn)?shù)云圖分布不一致且和單步總包反應(yīng)計(jì)算結(jié)果存在較大差異。

4 結(jié) 論

a）單步總包反應(yīng)和多步化學(xué)反應(yīng)模型下，數(shù)值計(jì)算壓強(qiáng)值與試驗(yàn)壓強(qiáng)值一致，驗(yàn)證了模型的正確性。若單純考察內(nèi)流場(chǎng)的壓強(qiáng)分布，為提高計(jì)算效率，單步總包反應(yīng)可以替代多步化學(xué)反應(yīng)。

b）單步與多步化學(xué)反應(yīng)模型溫度場(chǎng)的數(shù)值計(jì)算結(jié)果存在較大偏差，兩種化學(xué)反應(yīng)模型下的內(nèi)流場(chǎng)最高溫度分別為4550 K和3530 K，單步化學(xué)反應(yīng)忽略了很多中間吸熱反應(yīng)，因此內(nèi)流場(chǎng)溫度值偏高失真。

c）單步與多步化學(xué)反應(yīng)模型組分場(chǎng)（O2、CO2和H2O）的數(shù)值計(jì)算結(jié)果存在較大偏差，多步化學(xué)反應(yīng)下O2在燃燒室的噴射長(zhǎng)度較單步長(zhǎng)，且CO2和H2O組分的空間分布差異較大。在考察燃燒室的組分分布時(shí)，多步化學(xué)反應(yīng)模型的計(jì)算結(jié)果在一定程度上能更接近真實(shí)情況。

d）燃燒室圓柱段近壁區(qū)域煤油質(zhì)量分?jǐn)?shù)比較高，煤油蒸發(fā)吸熱，使得該區(qū)域溫度較中心區(qū)低；而燃燒室的高溫區(qū)主要集中在燃燒室收斂段和喉部之間，驗(yàn)證了設(shè)置氣膜冷卻系統(tǒng)的必要性，此外，氣膜冷卻系統(tǒng)冷卻區(qū)域選取可以參考內(nèi)流場(chǎng)的溫度數(shù)值計(jì)算結(jié)果。