郭康康，聶萬勝，劉 瑜，蘇凌宇，石天一

（航天工程大學，北京，101416）

0 引 言

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

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

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

本文采用的試驗發(fā)動機如圖1所示。

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

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

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

2 試驗結(jié)果

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

圖2 試驗壓強值曲線Fig.2 Pressure of Experiment

3 針對試驗工況的數(shù)值模擬

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

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

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

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

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

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

上述控制方程的離散和求解在開源程序平臺[17]（Open Source Field Operation and Manipulation，OpenFOAM）實現(xiàn)，OpenFOAM具有良好的開源性，利用C++語言把CFD方程組編寫為面向?qū)ο蟮膹姶蟮讓宇悗?，通過調(diào)用特定功能的類庫形成滿足具體應用問題的求解器，實現(xiàn)偏微分方程組的數(shù)值求解。

3.2 化學反應模型

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

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

表2 單步化學反應機理Tab.2 Single Step Global Chemical Reaction

多步反應參見文獻[19]和文獻[20]，把煤油與氧氣的反應簡化為9組分14步，化學反應速率如下式：

多步化學反應機理如表3所示[19,20]。

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

表3 多步化學反應機理Tab.3 nMulti-step Chemical Reactio

續(xù)表3

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

計算網(wǎng)格為非結(jié)構(gòu)網(wǎng)格，因為考察的是試驗發(fā)動機的內(nèi)流場特征，為減少計算量，從噴管段截去一部分。試驗網(wǎng)格如圖 3所示，試驗發(fā)動機內(nèi)流場邊界條件如表4所示。

圖3 試驗發(fā)動機網(wǎng)格示意Fig.3 Grid of Experimental Rocket Engine

表4 發(fā)動機內(nèi)流場邊界條件Tab.4 Boundary Conditions of Internal Flow Field

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

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

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

分別記單步總包反應與9組分14步化學反應模型下的算例為case1和case2，監(jiān)測點壓強如圖5所示。

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

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

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

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

圖6 基于單步與多步化學反應對稱面溫度云圖Fig.6 Temperature Contours of Chamber Based on Single-step Reaction and Multi-step Reactions

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

圖7 基于多步化學反應對稱面煤油質(zhì)量分數(shù)云圖Fig.7 Mass Fractions Contour of C12H24 Based on Multi-step Reactions

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

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

續(xù)圖8

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

case1和case2的CO2和H2O組分分布差異較大，在case1中，化學反應一步完成，同時產(chǎn)生CO2和H2O，所以兩者質(zhì)量分數(shù)云圖分布幾乎相同；而在case2中，CO2和H2O在不同的中間反應過程生成，生成CO2和H2O的中間反應過程在時間和空間上均不一致。因此，CO2和H2O的質(zhì)量分數(shù)云圖分布不一致且和單步總包反應計算結(jié)果存在較大差異。

4 結(jié) 論

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

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

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

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