張　濤，李國岫，陳　君,，虞育松，魏延明，劉旭輝，周廣文

(1.北京交通大學(xué) 機(jī)械與電子控制工程學(xué)院，北京　100044；2.北京控制工程研究所，北京　100190)

?

AP/HTPB復(fù)合固態(tài)微推力器燃燒過程仿真

針對高氯酸銨/端羥基聚丁二烯(AP/HTPB)混合固體推力器燃燒過程，采用燃速公式和氣相化學(xué)反應(yīng)動力學(xué)模型相結(jié)合的方法建立AP/HTPB固體推進(jìn)劑的燃燒模型；同時利用用戶自定義函數(shù)(UDF)來模擬固體推進(jìn)劑的燃面退移過程，研究得到固體微推力器內(nèi)復(fù)雜的燃燒反應(yīng)過程，以及AP/HTPB推進(jìn)劑配比對燃燒過程的影響規(guī)律。結(jié)果表明，AP/HTPB固態(tài)推力器進(jìn)行著復(fù)雜的物理化學(xué)過程，所得到的計算結(jié)果與實驗具有較好的一致性；同時，在推進(jìn)劑配比逐漸變化過程中，平均推力呈現(xiàn)出先增后減的趨勢，得到當(dāng)AP/HTPB=85/15時，推力器獲得最佳性能。

固態(tài)微推力器；高氯酸銨/端羥基聚丁二烯；推進(jìn)劑配比；燃燒

0　引言

微納衛(wèi)星由于其具有較傳統(tǒng)大尺寸衛(wèi)星設(shè)計周期短、成本低、可靠性高等優(yōu)點(diǎn)，受到各國專家和學(xué)者的關(guān)注[1-3]。這其中基于AP/HTPB推進(jìn)劑的固態(tài)微推力器中，可為微納型衛(wèi)星的姿態(tài)進(jìn)行連續(xù)精確調(diào)節(jié)[4]。因此，對其進(jìn)行數(shù)值模擬研究具有重要的意義。Alexeenko A A等[5]采用數(shù)值方法研究了高雷諾數(shù)、熱邊界條件、推力器三維幾何高度等對推力器性能的影響。Orieux S等[6]提出了一種預(yù)測固體微推力器性能的集總參數(shù)簡化模型。Hu等[7-8]獲得了不同固體推進(jìn)劑在不同配比情況下的微推力器內(nèi)燃燒過程仿真計算結(jié)果，之后他們探討了推力器內(nèi)燃燒過程中熱損失對推力器性能的影響。Sebasti?o等[9]仿真計算得到了固體推力器微噴管內(nèi)熱量傳遞以及壓力分布等結(jié)果。在國內(nèi)，胡松啟等[10]針對一種基于MEMS技術(shù)的固體微推力器，分析了瞬態(tài)燃燒效應(yīng)和推力器喉部尺寸對推進(jìn)性能的影響。林博穎等[11]針對微推力器采用零維瞬態(tài)模型，結(jié)合詳細(xì)化學(xué)反應(yīng)動力學(xué)機(jī)理，研究了微推力器的微尺度特性和啟動特性，同時對比研究了多孔介質(zhì)對微推力器工作性能的影響?？梢钥闯?，目前對于固體推進(jìn)劑的燃燒方面，多采用瞬態(tài)燃燒模型計算。

本文利用數(shù)值模擬的方法對固體微推力器的工作過程進(jìn)行了仿真研究。采用燃速公式和氣相化學(xué)反應(yīng)動力學(xué)機(jī)理相結(jié)合的方法建立了固體推進(jìn)劑AP/HTPB的燃燒模型。利用FLUENT軟件中的動網(wǎng)格技術(shù)，結(jié)合UDF進(jìn)行二次開發(fā)，模擬固體推進(jìn)劑的燃面退移過程。主要討論了AP/HTPB復(fù)合固態(tài)推力器內(nèi)進(jìn)行的復(fù)雜的物理化學(xué)過程，以及推進(jìn)劑配比對燃燒過程影響規(guī)律。

1　計算模型

1.1幾何模型

圖1為微推力器幾何結(jié)構(gòu)圖。采用文獻(xiàn)[12]中的計算模型，其中燃燒室長度為1 100 μm，噴管收縮段和擴(kuò)張段長度均為190.95 μm，喉部直徑為40 μm，微推進(jìn)器底部壁面厚度為100 μm，燃燒室直徑為350 μm，燃燒室段壁面厚度為350 μm，噴管收縮及擴(kuò)張半角均為35.26°。固態(tài)燃料點(diǎn)火后迅速完全燃燒，將燃料的化學(xué)能轉(zhuǎn)化為動能，經(jīng)噴管產(chǎn)生推力。

圖1　微推力器幾何模型示意圖

1.2數(shù)學(xué)模型

分別求解動網(wǎng)格區(qū)域內(nèi)連續(xù)性方程、動量方程以及能量方程，推力計算公式：

(1)

式中m為噴管出口流量；c為噴管出口流速；pe為噴管出口壓力；pa為大氣壓力；Ae為噴管出口截面積。

根據(jù)文獻(xiàn)[12]，AP/HTPB的燃速公式采用指數(shù)公式形式：

r=apcn

(2)

式中r為推進(jìn)劑燃面的退移速度；a為燃速系數(shù)，a=6.73×10-5kg-nmn+1s2n-1；n為壓強(qiáng)指數(shù)，n=0.4；pc為燃燒室壓力。

對于推進(jìn)劑的燃燒過程分為兩個部分。固體推進(jìn)劑被點(diǎn)燃后，推進(jìn)劑在燃面附近經(jīng)凝聚相反應(yīng)區(qū)后進(jìn)入氣相反應(yīng)區(qū)。凝聚相反應(yīng)采用燃速公式和動網(wǎng)格UDF結(jié)合的方法，模擬推進(jìn)劑的燃面退移，其首先產(chǎn)生高溫氣態(tài)的HClO4、NH3、C2H4和C4H6，反應(yīng)過程如公式(3)所示。

AP/HTPB=0.68HClO4+0.12NH3+0.05C4H6+

0.15C2H4

(3)

上述反應(yīng)產(chǎn)生的4種氣相組分作為氣相反應(yīng)區(qū)的反應(yīng)物繼續(xù)分解燃燒。采用Jeppson M B等[13]提出的72步反應(yīng)和36種組分的化學(xué)反應(yīng)機(jī)理。

1.3初始條件及求解

入口邊界采用質(zhì)量流量入口，流量大小和燃面推移規(guī)律由推進(jìn)劑燃速公式和UDF的動網(wǎng)格同時確定。出口采用壓力邊界，且壁面采用無滑移邊界條件。

2　計算結(jié)果與討論

2.1模型驗證

根據(jù)文獻(xiàn)[12]中的固體微推力器結(jié)構(gòu)，利用上述方法對固體微推進(jìn)器的工作過程進(jìn)行了仿真研究，并將計算結(jié)果與文獻(xiàn)中的實驗結(jié)果、瞬態(tài)燃燒模擬結(jié)果進(jìn)行了比較，如表1所示。

表1　3種研究方法的結(jié)果對比

從表1可見，由于燃速公式結(jié)合化學(xué)反應(yīng)動力學(xué)機(jī)理的方法相對更全面地考慮了推進(jìn)劑的燃燒過程，因此本文計算結(jié)果相對于瞬態(tài)燃燒模擬結(jié)果更接近實驗值。其原因是，當(dāng)對微推力器中的模型依然采用瞬態(tài)燃燒模擬方法時，誤差會很大，而燃速公式結(jié)合化學(xué)反應(yīng)動力學(xué)機(jī)理的方法是針對固體推進(jìn)劑的整個燃燒過程而言的，準(zhǔn)確性進(jìn)一步提高。

2.2固體推進(jìn)劑燃面推移過程

圖2是不同時刻(從上至下分別為5、20、40、60、90 ms)燃燒室壓力、速度及HClO4組分分布云圖。從圖2可看出，隨反應(yīng)的進(jìn)行，推進(jìn)劑燃面逐漸向左退移，燃燒壓力逐漸升高，且較均勻，噴管擴(kuò)張段軸線上出現(xiàn)了多個低壓區(qū)域。燃?xì)庠趪姽軆?nèi)逐漸加速，因為粘性邊界層的影響，沿中軸線方向上的流速更高，所以壓力就更低，甚至小于邊界層的壓力，在擴(kuò)張段軸線上出現(xiàn)了低壓區(qū)域。

(a) 5 ms,壓力　　　(b) 5 ms,速度　　　(c) 5 ms,HClO4

(d) 20 ms,壓力　　　(e) 20 ms,速度　　　(f) 20 ms,HClO4

(g) 40 ms,壓力　　　(h) 40 ms,速度　　　(i) 40 ms,HClO4

(j) 60 ms,壓力　　　(k) 60 ms,速度　　　(l) 60 ms,HClO4

(m) 90 ms,壓力　　　(n) 90 ms,速度　　　(o) 90 ms,HClO4

同時從速度分布中可看出，燃燒室內(nèi)的流速很低，在40 m/s以內(nèi)。通過拉法爾噴管后，燃?xì)饬魉偌彼偕?，可觀察到在噴管出口處的中軸線上出現(xiàn)了馬赫盤。這是因為粘性作用，壁面附近存在速度邊界層。邊界層內(nèi)的燃?xì)饬魉佥^低，而沿固體微推力器軸線上的燃?xì)饬魉佥^高，軸線附近的高速燃?xì)庖蚴苓吔鐚訁^(qū)域的流體粘性力的作用而減速，但是其上游的燃?xì)馊砸愿咚傧蛳聸_擊，而下游的燃?xì)庖蛩俣容^低就對上游燃?xì)庑纬烧系K，導(dǎo)致在噴管擴(kuò)張段出現(xiàn)了質(zhì)量壅塞，上游燃?xì)獗粔嚎s，從而形成激波。產(chǎn)生激波以后，激波面后面的高壓通過邊界層內(nèi)的亞聲速區(qū)向上游傳播到激波前區(qū)，形成了逆向壓力梯度。由于逆向壓力梯度的存在，導(dǎo)致了邊界層內(nèi)流場狀態(tài)的改變，激波的形態(tài)因此而發(fā)生變化，在微噴管的軸線附近出現(xiàn)了連續(xù)的激波，激波間也就形成了二次馬赫盤。而從不同時刻固體微推進(jìn)器內(nèi)HClO4(來源于AP)的濃度分布，可得到HClO4進(jìn)入燃燒室后，快速分解，在距燃面約80 μm處，HClO4已經(jīng)完全分解。隨著HClO4的分解，產(chǎn)生大量的氧化性物質(zhì)，驅(qū)動了整個燃燒反應(yīng)的進(jìn)行。

2.3推進(jìn)劑配比的影響

圖3是不同AP/HTPB推進(jìn)劑質(zhì)量配比下，固體微推進(jìn)器工作過程中燃燒室內(nèi)壓力變化曲線對比。從圖3可看出，燃燒室壓力曲線整體呈先快速升高而后緩慢降低的趨勢。其原因是，反應(yīng)初期燃?xì)馀c壁面接觸面積小，隨著固體推進(jìn)劑燃面的退移，燃?xì)馀c壁面的接觸面積加大，導(dǎo)致熱損失的增加，燃燒溫度降低、壓力降低。AP/HTPB配比為85/15時，平均壓力最高，但與配比為80/20工況的壓力計算結(jié)果差別不大。當(dāng)推進(jìn)劑配比為70/30時，燃燒室內(nèi)的壓力最低。比最高壓力(配比為85/15時)低0.05 MPa。結(jié)果表明，配比對推進(jìn)劑的燃燒過程有較大影響。當(dāng)配比為85/15時，固體微推進(jìn)器性能較好。其原因是，隨著AP含量的增加，與HPTB反應(yīng)的氧化劑增加，燃燒壓力升高，但是當(dāng)氧化劑含量過多，而燃料HPTB相對不足時，燃燒反應(yīng)中CO、H2等小分子降低，壓力也隨之降低。因此，存在一個最佳的氧化劑和燃料的配比。

圖4是推進(jìn)劑燃燒終了時刻不同AP/HTPB推進(jìn)劑質(zhì)量配比下的推進(jìn)劑的燃燒溫度曲線。從圖4可見，AP/HTPB配比中，隨著AP的增多，燃燒溫度逐漸升高。由1 231、1 320、1 420、1 500 K(AP/HTPB=70/30)增加到1 585 K(AP/HTPB=90/10)，基本呈線性增加的趨勢。這是因為，AP屬強(qiáng)氧化劑，反應(yīng)速率快；而HTPB是碳?xì)漕惛叻肿泳酆衔?，其分解燃燒需要AP分解提供氧化劑。當(dāng)AP含量升高時，氧化劑含量高，反應(yīng)速率加快，燃燒相對更充分且放熱率高，所以燃燒溫度也更高。

圖5是不同AP/HTPB推進(jìn)劑質(zhì)量配比下，推進(jìn)劑的燃燒時間對比。

圖3　不同配比下燃燒室內(nèi)壓力變化對比曲線

圖4　不同配比下推進(jìn)劑的燃燒溫度對比曲線

圖5　不同配比下推進(jìn)劑燃燒時間對比曲線

由圖5可看出，當(dāng)AP/HTPB配比為85/15時，推進(jìn)劑燃燒時間最短，為89.5 ms；而當(dāng)配比為70/30時，燃燒時間最長，為95 ms?？煽闯?，當(dāng)AP/HTPB配比從70/30變化到85/15的過程中，推進(jìn)劑燃燒溫度逐漸升高，燃燒時間逐漸變短，這是因為當(dāng)推進(jìn)劑中AP含量高時，反應(yīng)更加充分，放熱量更高，所以燃燒溫度、壓力等更高，燃燒速度也更快。但當(dāng)配比達(dá)到AP/HTPB=90/10時，雖然燃燒溫度更高，但由于推進(jìn)劑中HTPB含量降低，導(dǎo)致了燃燒產(chǎn)物中CO、H2等小分子減少，燃速減慢。因此，配比為AP/HTPB=90/10時，雖燃燒溫度更高，但燃燒時間更長。

圖6是推進(jìn)劑燃燒終了時刻，重要中間氣相產(chǎn)物沿固體微推力器軸線上的質(zhì)量分?jǐn)?shù)分布。

(a) OH

(b) ClO

(c) NO

(d) CH2O

由圖6可看出，反應(yīng)物進(jìn)入推力器后很快分解燃燒。隨著AP/HTPB配比從70/30到90/10的變化過程中，OH存在反應(yīng)區(qū)域變厚，且OH的含量也逐漸加大。這也說明了隨著配比中AP含量的增加，氧化劑含量高。而組分ClO，只存在于燃面附近0.1 mm以內(nèi)的空間區(qū)域，而后消失，這是因為HClO4反應(yīng)活性高，含Cl反應(yīng)相對較快，反應(yīng)較為迅速。NO在不同配比條件下基本被氧化完全，但是當(dāng)氧化劑含量較高時(配比為90/10)，NO并不能充分被氧化，濃度較高。而隨著AP濃度的逐漸降低，CH2O質(zhì)量分?jǐn)?shù)維持在較高水平，隨著AP濃度的升高逐漸降低。

圖7是不同AP/HTPB推進(jìn)劑質(zhì)量配比下，平均推力的變化曲線。從圖7可看出，平均推力在AP/HTPB配比從70/30到90/10的變化過程中，呈現(xiàn)出先升高后降低的趨勢，當(dāng)配比為AP/HTPB=85/15時，平均推力最大，達(dá)到了1.75 mN。原因是HTPB的含量并不是越少越好，因為HTPB參與反應(yīng)后，會生成CO、H2、CO2等小分子組分，能夠使燃燒室壓力升高，進(jìn)而提升推力。

圖7　不同配比下平均推力的變化曲線

3　結(jié)論

(1)在AP/HTPB復(fù)合固態(tài)推力器中，進(jìn)行著復(fù)雜的物理化學(xué)過程。采用將燃速公式和氣相化學(xué)反應(yīng)動力學(xué)機(jī)理相結(jié)合方法獲得的計算結(jié)果，比文獻(xiàn)中提出的瞬態(tài)燃燒模擬的方法得到的燃燒時間和推力計算結(jié)果，更加接近于實驗值。

(2)隨著反應(yīng)的進(jìn)行，推進(jìn)劑燃面逐漸向左退移，燃燒壓力逐漸升高，在噴管擴(kuò)張段軸線上出現(xiàn)了多個低壓區(qū)域，并形成激波。同時HClO4在距燃面約80 μm處，HClO4已經(jīng)完全分解。

(3)推進(jìn)劑配比影響中間產(chǎn)物時間和空間分布。隨AP含量增加，OH存在反應(yīng)區(qū)域變厚且OH峰值濃度也逐漸加大，CH2O被氧化，濃度逐漸降低。

(4)當(dāng)氧化劑AP含量高時，即AP/HTPB配比從70/30到90/10的變化過程中，推進(jìn)劑的燃燒溫度是逐漸升高的，從70/30時的1 231 K增加到90/10時的1 585 K，增幅為28.8%。

(5)隨著AP/HTPB配比的逐漸變化中，燃燒時間呈現(xiàn)先減后增的趨勢。在AP/HTPB=85/15時，燃燒時間最短，為89.5 ms。同時推力也呈現(xiàn)出先升高后降低的趨勢，當(dāng)配比為AP/HTPB=85/15時，平均推力最大，達(dá)到了1.75 mN。

[1]余文革, 鐘先信, 李曉毅，等. 皮衛(wèi)星發(fā)展展望[J]. 壓電與聲光, 2004, 26(4): 289-292.

[2]Granger R, Dalgleish B. Developing nano/micro satellite technology for government and military applications[C]. Milsatcoms 2015, 2015.

[3]Underwood C, Pellegrino S, Lappas V J, et al. Using cubesat/micro-satellite technology to demonstrate the autonomous assembly of a reconfigurable space telescope (AAReST)[J]. Acta Astronautica, 2015, 114: 112-122.

[4]Cai W, Thakre P, Yang V. A model of AP/HTPB composite propellant combustion in rocket-motor environments[J]. Combustion Science and Technology, 2008, 180(12): 2143-2169.

[5]Alexeenko A A, Levin D A, Fedosov D A, et al. Performance analysis of microthrusters based on coupled thermal-fluid modeling and simulation[J]. Journal of Propulsion and Power, 2005, 21(1): 95-101.

[6]Orieux S, Rossi C, Esteve D. Compact model based on a lumped parameter approach for the prediction of solid propellant micro-rocket performance[J]. Sensors and Actuators A: Physical, 2002, 101(3): 383-391.

[7]Hu S, Chen J, Hua Y, et al. Combustion behaviour of solid propellants for micro-thruster[J]. Materials Testing, 2014, 56(5): 399-407.

[8]Hu S, Chen J, Wu G, et al. Influences of heat loss on micro thruster performance[J]. Combustion Science and Technology, 2015 (just-accepted).

[9]Sebasti?o I B, Santos W F N. Numerical simulation of heat transfer and pressure distributions in micronozzles with surface discontinuities on the divergent contour[J]. Computers & Fluids, 2014, 92: 125-137.

[10]胡松啟，高勝靈, 劉凱，等. 固體微推進(jìn)器工作過程數(shù)值分析[J].西北大學(xué)學(xué)報, 2012, 42(1): 21-24.

[11]林博穎，張根烜，劉明侯，等. 微小型化學(xué)推進(jìn)器的性能分析[J].中國科學(xué)技術(shù)大學(xué)學(xué)報, 2008, 38(1): 94-99,104.

[12]Zhang K L, Chou S K, Ang S S. Performance prediction of a novel solid-propellant microthruster[J]. Journal of Propulsion and Power, 2006, 22(1): 56-63.

[13]Jeppson M B, Beckstead M W, Jing Q. A kinetic model for the premixed combustion of a fine AP/HTPB composite propellant[D]. Brigham Young University. Department of Chemical Engineering,1998.

(編輯：呂耀輝)

Numerical simulation of combustion process of AP/HTPB solid micro thruster

ZHANG Tao1, LI Guo-xiu1, CHEN Jun1, 2, YU Yu-song1,WEI Yan-ming2, LIU Xu-hui2, ZHOU Guang-wen1

(1.School of Mechanical, Electronic and Control Engineering, Beijing Jiaotong University, Beijing100044, China；2.Beijing Institute of Control Engineering, Beijing100190, China)

Within the present investigation, the combustion characteristics of ammonium perchlorate/hydroxyl terminated polybutadiene (AP/HTPB)-based micro-scale thruster were numerically studied by using the method of the combination of burning rate formula and chemical reaction kinetics model. Meanwhile, the combustion surface of solid thruster was described by the user-defined function (UDF). The effect of AP/HTPB ratio on the combustion process was investigated. Based on the calculated results, it can be found that, the complex physical and chemical processes were discussed, and the results show a good consistency with the experimental result. Moreover, the thruster performance showed a trend of first increasing then decreasing with the increasing of AP/HTPB ratio. When AP/HTPB ratio was 85/15, thruster could acquire better performance.

solid micro thruster；AP/HTPB；propellant mass ratio；combustion

2015-09-25；

2015-11-05。

張濤(1987—)，男，博士生，研究方向為航天器推進(jìn)系統(tǒng)。E-mail:tzhang903@bjtu.edu.cn

張濤1，李國岫1，陳君1,2，虞育松1，魏延明2，劉旭輝2，周廣文1

(1.北京交通大學(xué) 機(jī)械與電子控制工程學(xué)院，北京100044；2.北京控制工程研究所，北京100190)

V439

A

1006-2793(2016)05-0632-06

10.7673/j.issn.1006-2793.2016.05.006