馬航遠，向紅軍

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

0 引言

帶潛入噴管的固體火箭發(fā)動機在其工作的拖尾階段，由于推進劑即將耗盡，燃燒室內(nèi)的壓力迅速降低。在此過程中，潛入噴管背壁區(qū)沉積的液態(tài)高溫熔渣有可能通過噴管排出[1]，從而導(dǎo)致噴管出口鄰近區(qū)域內(nèi)熱通量的異常升高，甚至損壞附近裝置。本文對熔渣在固體火箭發(fā)動機拖尾階段的排出過程進行了數(shù)值模擬，力圖揭示熔渣排出的相關(guān)機理。

關(guān)于潛入式固體火箭發(fā)動機中熔渣的研究，研究內(nèi)容主要集中在熔渣沉積過程和影響因素兩個方面。Borass等[2]通過模擬航天飛機固體火箭發(fā)動機內(nèi)熔渣沉積過程，提出了一種預(yù)測與沉積和堆積過程有關(guān)的某些參數(shù)分析方法，并通過試驗，定性地證實了其分析結(jié)果。Hopson等[3]通過有限元法結(jié)合空氣動力學(xué)，研究了引起發(fā)動機壓強脈動的熔渣排出機理。Paul Liaw等[4]研究了化學(xué)反應(yīng)和顆粒相對二維軸對稱固體火箭發(fā)動機的影響。此外，他還研究了考慮化學(xué)反應(yīng)、顆粒蒸發(fā)、燃燒和粒子破碎模型時固體火箭發(fā)動機內(nèi)多相流的問題[5]。在他的另一個研究[6]中，Liaw采用歐拉-拉格朗日顆粒軌跡法，結(jié)合VOF(Volume of Fluid)方法研究了固體火箭發(fā)動機背壁區(qū)熔渣在增壓過程中的變化和排出問題。胡春波等[7]利用X射線診斷系統(tǒng)，觀察了背壁區(qū)熔渣粒子的溢流過程，總結(jié)了噴管喉部面積、潛入深度等對熔渣積累和溢出過程的影響。向紅軍等[8]采用歐拉-拉格朗日兩相方法，建立了燃燒室噴管統(tǒng)一流場，結(jié)合顆粒確定軌道模型、Hermsen燃燒模型、粒子破碎模型，模擬了噴管潛入式固體火箭發(fā)動機內(nèi)的兩相流場。魏超等[9]通過試驗和數(shù)值仿真，確定了燃燒室中粒子直徑分布以及熔渣捕獲判據(jù)。李強等[10]等通過對熔渣的沉積過程進行了數(shù)值仿真，得到一些結(jié)論。彭科科等[11]系統(tǒng)介紹了美國的NASA熔渣模型、歐空局的MASTER2005熔渣模型和麻省理工林肯實驗室的MIT/LL模型，分析了熔渣尺寸-數(shù)量和質(zhì)量特性，為進一步分析熔渣對于空間碎片環(huán)境的影響奠定一定基礎(chǔ)，對國內(nèi)自主開展熔渣模型研究也有一定的參考價值。

研究方法分為試驗研究和數(shù)值仿真兩種，試驗研究側(cè)重現(xiàn)象描述，而數(shù)值仿真?zhèn)戎貦C理分析。

綜上所述，關(guān)于熔渣沉積過程的研究已經(jīng)取得了一定的成果，但熔渣排出方面的研究，目前十分欠缺。本文通過數(shù)值仿真方法，對熔渣的排出機理及其影響因素做了初步的探索性研究。

1 模型與算法

相關(guān)研究表明，只有當潛入噴管背壁區(qū)液態(tài)熔渣沉積量達到一定深度后，才出現(xiàn)顯著的熔渣排出現(xiàn)象。熔渣排出過程中主要涉及到氣液兩相間的相互作用以及液相與固壁間的相互作用。氣液界面上存在壓差力、粘性力和表面張力，液相與固壁間存在粘性力與附著力。此外，由于壁面溫度較低還有可能出現(xiàn)液相的凝固現(xiàn)象。震動、外部過載等對熔渣的排出也會有影響，要準確全面地模擬整個熔渣沉積與排出的動態(tài)過程是非常困難的。因此，研究中作了如下假設(shè)：

(1)考慮到拖尾段裝藥基本燃燒完畢，忽略裝藥的三維特征，將發(fā)動機幾何模型簡化為二維軸對稱結(jié)構(gòu)，且只取燃燒室后半段和噴管區(qū)域進行計算；

(2)拖尾段開始時熔渣沉積已經(jīng)停止，不考慮殘余燃氣中的凝相成分，按純氣相處理。熔渣池的初始形狀和邊界根據(jù)經(jīng)驗進行簡單設(shè)定；

(3)重點關(guān)注熔渣池主體在燃氣沖刷下的演化過程，忽略熔渣相變現(xiàn)象；

(4)忽略任何外部過載，包括重力作用。

1.1 VOF模型

VOF模型適用于兩種及兩種以上互不穿透的流體間界面的跟蹤，通過求解引入模型的每一相流體在控制單元內(nèi)的體積分數(shù)，確定相間界面。對于文中氣液兩相流動，設(shè)每個控制單元中氣相和液相的體積分數(shù)分別為αa(x，y，t)和αl(x，y，t)，則有

αa+αl=1

(1)

在一個單元格中，每一相體積分數(shù)α的值有三種情況：α=0，該單元格沒有相；α=1，該單元格充滿相；0<α<1，該單元格同時包含兩相，為相間界面。

對于液相，其體積分數(shù)方程為

▽(αlρlvl)=0

(2)

式中αl、ρl和vl分別表示液相的單元格體積分數(shù)、液相密度和速度。

對于每一個控制單元格中的流體，有動量方程：

(3)

其中

ρ=(1-αl)ρa+αlρl

μ=(1-αl)μa+αlμl

(4)

(5)

1.2 湍流模型

采用Launder和Spalding提出的標準k-ε雙方程湍流模型，對流場的湍流進行數(shù)值模擬。

1.3 幾何模型與網(wǎng)格劃分

本研究的主要思路：首先建立熔渣池模型，然后模擬出該發(fā)動機正常工作狀態(tài)的穩(wěn)態(tài)流場，最后以該穩(wěn)態(tài)流場為拖尾階段的起始場進行熔渣排出的模擬。

以國外某潛地導(dǎo)彈的第三級固體火箭發(fā)動機為原型(主要參數(shù)有噴管喉徑125 mm，噴管擴張比34.2，燃燒室平均壓強6.4 MPa)。為加快運算速度，現(xiàn)截取燃燒室后半部分，并按1/4的縮比建立幾何模型。其中，壁面01為穩(wěn)態(tài)場計算中的相間壁面，也是非穩(wěn)態(tài)流動中氣液相間界面的初始位置，該邊界左側(cè)為燃氣，右側(cè)為熔渣池，如圖1所示。

圖1 幾何模型示意圖Fig.1 Geometry configuration of the model

建立幾何模型后，劃分結(jié)構(gòu)化網(wǎng)格。其主要參數(shù)有單元格11 154個，面22 657個，節(jié)點11 521個，沿壁面部分用O型網(wǎng)格局部加密。網(wǎng)格模型如圖2所示。

圖2 網(wǎng)格模型Fig.2 Mesh of the model

1.4 邊界條件

參考該型號的實際工作參數(shù)，給定適當?shù)倪吔鐥l件，模擬出作為拖尾階段初場的穩(wěn)態(tài)場。主要參數(shù)設(shè)置見表1，與該入口邊界條件相匹配的入口靜壓為6.4 MPa。

通過上述條件計算得到的穩(wěn)態(tài)場為t=0時刻的初始場，通過UDF控制入口流率，修改兩相間壁面邊界條件為內(nèi)部邊界，讓熔渣池與穩(wěn)態(tài)流場接觸，模擬熔渣排出的非穩(wěn)態(tài)流動過程。主要參數(shù)設(shè)置見表2，即假定拖尾段入口流率經(jīng)過0.5 s線性下降到0。

表2 非穩(wěn)態(tài)場主要設(shè)置參數(shù)

1.5 物性參數(shù)

本仿真中固體推進劑為TP-1148推進劑，在平均壓強6.4 MPa下，經(jīng)發(fā)動機熱力學(xué)計算得到燃燒室溫度為3434 K，燃氣其他各物性參數(shù)見表3。

熔渣池有液態(tài)Al2O3組成。Glorieux等[12-13]通過實驗測得了液態(tài)Al2O3的密度及表面張力系數(shù)隨溫度變化的關(guān)系式：

ρ=3084.345-0.117 738T(kg/m3)

2000 K≤T≤3100 K

(6)

σ=0.747 5-0.000 039T(N/m)

2300 K≤T≤3200 K

(7)

Paradis等[14-15]通過實驗得到其粘性和定壓比熱容隨溫度T變化的關(guān)系式：

μ=0.003 2 exp(5196.05/T) (Pa·s)

2190 K≤T≤2500 K

(8)

cp= 153.5+0.003 1(T-2327) (J/(mol·K))

2175 K≤T≤2435 K

(9)

由于燃氣溫度為3434 K，根據(jù)上述各式外插計算并轉(zhuǎn)換單位得到液態(tài)Al2O3的物性參數(shù)列于表3。其中，熱導(dǎo)率取7.4 W/(m·K)[16]，表面張力系數(shù)為0.613 6 N/m。

表3 燃氣和液態(tài)Al2O3的物性參數(shù)

1.6 算法選擇

采用SIMPLE算法，對壓力場和速度場進行耦合計算，密度、動量和能量方程均采用二階迎風(fēng)格式進行離散，k和ε方程均采用一階迎風(fēng)格式進行離散，體積分數(shù)方程采用幾何重構(gòu)格式進行離散。

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

2.1 穩(wěn)態(tài)場

按表1中參數(shù)設(shè)置算例，得到用于非穩(wěn)態(tài)場計算的穩(wěn)態(tài)初始場。圖3為穩(wěn)態(tài)場兩相的密度分布云圖，熔渣池和氣相穩(wěn)態(tài)場由相間界面隔開分布。圖4為穩(wěn)態(tài)壓強分布云圖，燃燒室壓強達到實際工作狀態(tài)的平均壓強。此時，燃燒室與熔渣池間的氣液邊界為壁面邊界，壓強不受燃燒室影響。圖5為穩(wěn)態(tài)溫度分布云圖，燃燒室溫度達到3434 K。經(jīng)過氣液邊界的熱傳遞，達到穩(wěn)態(tài)時，熔渣池的溫度與燃燒室溫度達到平衡，亦為3434 K。圖6為穩(wěn)態(tài)時的流線圖。

圖3 穩(wěn)態(tài)場的密度分布云圖Fig.3 Density distribution of the steady flow-field

圖4 穩(wěn)態(tài)場的壓強分布云圖Fig.4 Pressure distribution of the steady flow-field

圖5 穩(wěn)態(tài)場的溫度分布云圖Fig.5 Temperature distribution of the steady flow-field

圖6 穩(wěn)態(tài)場的流線圖Fig.6 Streamline of the steady flow-field

從圖6可看到，氣體從入口流入，在熔渣池的氣液相界面左側(cè)沿徑向向外位置形成渦流，在徑向向內(nèi)位置沿界面貼壁流動，并經(jīng)噴管排出。

2.2 拖尾階段

以穩(wěn)態(tài)場為進入非穩(wěn)態(tài)計算的初始場，按照表2數(shù)據(jù)修改邊界條件，計算拖尾階段的熔渣排出過程。計算完成后，選取不同時刻液相體積分數(shù)分布云圖，可顯示出熔渣的流動、排出狀況，見圖7。

從圖7可看到，進入拖尾段后，熔渣池在氣流作用下產(chǎn)生變形和流動，氣液界面由中心向熔渣池內(nèi)側(cè)產(chǎn)生凹陷?？拷鼭撊雵姽苋肟谔幍娜墼N著壁面流入噴管，最后通過喉部并排出噴管。熔渣的整個流動、排出過程，為小量、貼壁、間歇性排出的過程?？拷紵彝獗谝粋?cè)，熔渣在氣流的作用下，形成如圖7所示形態(tài)，并保持相對穩(wěn)定。

經(jīng)分析，從t=0時刻開始，液相熔渣同時受到壓差力、粘性力和氣液表面張力的作用。沿徑向取t=0時刻兩相界面左側(cè)的壓強分布如圖8所示?？煽闯?，作用于氣液相間界面的壓強沿徑向向外先增大、后減小，即氣液間界面中心位置的壓強大于兩側(cè)位置的壓強。又由于液相不可壓，熔渣池在中心位置受到相對較大的壓強作用時，產(chǎn)生內(nèi)凹，將熔渣擠向上下兩側(cè)。在上側(cè)位置形成熔渣池沿壁面爬升狀態(tài)，在下側(cè)靠近噴管口的位置，熔渣沿壁面流入噴管，進而隨氣流排出。

(a)t=0.000 1 s

(b)t=0.035 6 s

(c)t=0.205 8 s

(d)t=0.309 6 s

(e)t=0.5 s圖7 拖尾階段的熔渣排出過程Fig.7 Expulsion of the slag during tail-off

另一方面，在氣液兩相界面之間存在的粘性力和表面張力的共同作用下，液相熔渣池在相間界面受到沿界面向上下兩側(cè)的力，表面液體被氣流帶向兩側(cè)。被帶向兩側(cè)的熔渣，與由壓強擠壓產(chǎn)生流動的熔渣一起，分別向上沿壁面爬升和沿噴管排出。

圖9為拖尾段t=0.035 6 s時刻的溫度分布。對比同一時刻的熔渣排出圖可看到，在噴管擴張段熔渣經(jīng)過的位置，溫度比周圍氣相高。這種“高溫包”排出噴管時，會加大對周邊環(huán)境的熱輻射。

圖8 初始時刻兩相界面左側(cè)的徑向壓強分布Fig.8 Pressure distribution along radial direction on leftside of the initial interface between the two phases

圖9 t=0.035 6 s時刻的溫度分布云圖Fig.9 Temperature distribution at t=0.035 6 s

2.3 影響因素

分別改變氣相粘性、氣液間表面張力系數(shù)和液相熔渣的密度，考察它們對液相熔渣排出的影響。

圖10為不同氣相粘性(圖例中的數(shù)字表示粘性大小，如0.000 184 kg/(m·s))下熔渣排出的質(zhì)量累積量隨時間的變化情況。從圖10可看出，熔渣排出量隨時間的走勢一致，而且隨著氣相粘性的增大，液相熔渣的最終排出量依次減小。增大氣相粘性后，氣相流場會在熔渣池左側(cè)上部產(chǎn)生更大的渦旋，峰值壓強的位置下移，使產(chǎn)生壁面爬升的熔渣量增多，進入噴管的熔渣量減少，因而導(dǎo)致熔渣的最終排出量減少。不同粘性下的熔渣最終排出比例如表4所示。

圖11為不同表面張力系數(shù)下熔渣排出的質(zhì)量累積量隨時間的變化情況。從圖11可看出，隨著相間表面張力系數(shù)的增大，熔渣排出量隨時間走勢一致，最終排出量依次增大。經(jīng)分析，大的表面張力系數(shù)產(chǎn)生大的表面作用力，使得氣流帶走液相熔渣的能力增強，因而導(dǎo)致熔渣的最終排出量增大。

圖10 不同氣相粘性下熔渣的排出累積量隨時間變化曲線Fig.10 Mass accumulation of slag dischargedwith gas viscosity over time

氣相粘性/[kg/(m·s)]4．6×10-59．2×10-51．84×10-4百分比/%62．9743．3436．44

不同表面張力系數(shù)下的熔渣最終排出比例如表5所示。圖12為不同液相密度下熔渣排出的體積累積量隨時間的變化情況。從圖12可看出，隨著相間液相密度的增大，熔渣排出量隨時間走勢和排出累積量基本一致。液相密度對熔渣的排出過程影響不大，對其最終排出量的影響也沒有明顯規(guī)律。

圖11 不同表面張力系數(shù)下熔渣排出累積量隨時間變化曲線Fig.11 Mass accumulation of slag dischargedwith tension coefficient over time

圖12 不同液相密度下熔渣的排出累積量隨時間變化曲線Fig.12 Volume accumulation of slag dischargedwith liquid density over time

表面力系數(shù)/(N/m)0．31360．61360．9136百分比/%39．3343．3453．35

3 結(jié)論

(1)液相熔渣在壓差力、粘性力和表面張力的共同作用下，氣液相間界面產(chǎn)生變形，進而隨氣流貼壁面排出噴管。其中，沿徑向變化的壓差力主要使熔渣池產(chǎn)生上下不均勻的變形，熔渣在相界面受壓較小的下側(cè)位置開始溢出。粘性力和表面張力主要使相界面產(chǎn)生向下的切向流動，并與受壓溢出的液體一起脫離熔渣池，隨氣流緊貼壁面排出。

(2)氣相粘性對熔渣的排出過程及最終排出量有明顯影響：氣相粘性越大，熔渣在背壁區(qū)產(chǎn)生的渦旋越明顯，進而導(dǎo)致熔渣的最終排出量越小。

(3)表面張力系數(shù)對熔渣的排出過程及最終排出量亦有明顯影響：表面張力系數(shù)越大，產(chǎn)生的表面作用力越大，使得氣流帶走液相的能力越強，進而導(dǎo)致熔渣的最終排出量越大。

(4)液相密度對熔渣的最終排出量影響不明顯，對排出過程的影響沒有明顯規(guī)律。

本研究作為初步的探索性研究，尚有一些不足。對于熔渣池大小、液相粘性等影響因素的討論尚不夠充分，這些將在后續(xù)探究中不斷進行完善。

參考文獻：

[1] Mehta M,Dufrene A T,Seaford M,et al.Space launch system base heating test:environments and Base Flow Physics[C]//Aiaa Aerospace Sciences Meeting,2015.

[2] Boraas S,徐濤清.航天飛機固體火箭發(fā)動機內(nèi)的熔渣沉積模擬[J].固體火箭技術(shù),1983,6(3):57-74.

Boraas S,XU Taoqing.Modeling slag deposition in the space shuttle solid rokets motor[J].Journal of Solid Rocket Technology,1983,6(3):57-74.

[3] Hopson C.Space shuttle solid rocket motor slag expulsion mechanisms[C]//AIAA Paper 95-2725,31st AIAA/ASME/SAE/ASEE Joint Propulsion Conference and Exhibit,San Diego,CA,Jul.10-12,1995.

[4] Pual Liaw,Chen Yensen,Shang Huanmin.Numerical analysis of particulate multi-phase flowfield for a 3-D advanced solid rocket motor[C]//AIAA Paper 94-0866,32nd Aerospace Sciences Meeting and Exhibit,Reno,NV,Jan.10-13,1994.

[5] Pual Liaw,Chen Yensen,Shang Huanmin.Particulate multi-phase flowfield calculation with combustion/breakup models for solid rocket motor[C]//AIAA Paper 94-2780,30th AIAA/ASME/SAE/ASEE Joint Propulsion Conference,Indianapolis,IN,Jun.27-29,1994.

[6] Pual Liaw,Chen Yensen,Shang Huanmin.Numerical investigation of the slag behavior in the aft-end cavity of the solid rocket motors[C]//AIAA Paper 95-0815,33rd Aerospace Sciences Meeting and Exhibit,Reno,NV,Jan.9-12,1995.

[7] 胡春波,何國強,劉佩進,等.潛入噴管背壁區(qū)熔渣溢流沉積實驗研究[J].推進技術(shù),2001,22(3):204-206.

HU Chunbo,HE Guoqiang,LIU Peijin,et al.Experimental study on slag accumulation in back well region of SRM submerged nozzle[J].Journal of Propulsion Technology,2001,22(3):204-206.

[8] 向紅軍.固體火箭發(fā)動機熔渣沉積數(shù)值模擬與實驗研究[D].北京:北京航空航天大學(xué),2000.

XIANG Hongjun.Numerical simulation and experimental study of slag accumulation for solid rocket motors[D].Beijing:Beijing University of Aeronautics and Astronautics,2000.

[9] 魏超,侯曉.潛入噴管背壁區(qū)熔渣沉積的機理分析與數(shù)值模擬[J].航空動力學(xué)報,2006,21(6):1109-1115.

WEI Chao,HOU Xiao.Mechanism investigation and numerical simulation of slag deposition for submerged nozzle[J].Journal of Aerospace Power,2006,21(6):1109-1115.

[10] 李強,甘曉松,劉佩進,等.大型固體發(fā)動機潛入式噴管背壁區(qū)域熔渣沉積數(shù)值模擬[J].固體火箭技術(shù),2010,33(2):148-151.

LI Qiang,GAN Xiaosong,LIU Peijin,et al.Numerical simulation of slag accumulation for submerged nozzle of large scale SRM[J].Journal of Solid Rocket Technology,2010,33(2):148-151.

[11] 彭科科,龐寶君,肖偉科,等.基于不同模型的固體火箭發(fā)動機噴射物熔渣特性分析[J].航天器環(huán)境工程,2012,29(6):616-620.

PENG Keke,PANG Baojun,XIAO Weike,et al.Characteristic analysis of slags based on different solid rocket motor firing models[J].Spacecraft Environment Engineering,2012,29(6):616-620.

[12] Glorieux B,Millot F,Rifflet J C,et al.Density of superheated and undercooled liquid alumina by a contactless method[J].International Journal of Thermophysics,1999,20(4):1085-1094.

[13] Glorieux B,Millot F,Rifflet J C.Surface tension of liquid alumina from contactless techniques[J].International Journal of Thermophysics,2002,23(5):1249-1257.

[14] Paradis P F,Ishikawa T.Surface tension and viscosity measurements of liquid and undercooled alumina by containerless techniques[J].Japanese Journal of Applied Physics,2005,44(7A):5082-5085.

[15] Paradis P F,Ishikawa T,Saita Y,et al.Non-contact thermophysical property measurements of liquid and undercooled alumina[J].Japanese Journal of Applied Physics,2004,43(4A):1496-1500.

[16] Sinn H,Glorieux B,Hennet L,et al.Microscopic dynamics of liquid aluminum oxide[J].Science,2003,299(5615):2047-2049.