謝愛元，武曉松，王向東

（1.南京理工大學(xué) 機(jī)械工程學(xué)院，南京210094；2.晉西工業(yè)集團(tuán)有限責(zé)任公司，太原030027）

固體燃料沖壓發(fā)動(dòng)機(jī)（Solid Fuel Ramjet，SFRJ）的燃面退移速率（簡稱燃速）不僅與固體燃料、來流空氣質(zhì)量流率、總溫有關(guān)，還和燃燒室結(jié)構(gòu)尺寸等息息相關(guān)[1－4]。目前國內(nèi)外的科研工作者通過大量的實(shí)驗(yàn)研究，擬合出單一因素下的燃速公式[1－4]，成本高昂。同時(shí)，試驗(yàn)表明，在其它條件均相似的情況下，尺寸越小的固體燃料沖壓發(fā)動(dòng)機(jī)具有更高的燃面退移速率，且差別較大[2－3]；即針對某一尺寸發(fā)動(dòng)機(jī)所得的實(shí)驗(yàn)結(jié)果，不能用于與其幾何相似的另一尺寸發(fā)動(dòng)機(jī)的工程設(shè)計(jì)中。由此可見，建立能夠?qū)Σ煌叽纭砹鳁l件下固體燃料沖壓發(fā)動(dòng)機(jī)燃速進(jìn)行準(zhǔn)確預(yù)示的模型顯得尤為重要。

假定工作過程中固體燃料通道內(nèi)腔的表面溫度（以下簡稱為內(nèi)表面溫度）為裝藥的液化溫度，且不變；流場中不存在回流區(qū)。在上述簡化的基礎(chǔ)上，文獻(xiàn)[5]建立了固體燃料沖壓發(fā)動(dòng)機(jī)的燃速模型。由于其簡化得太多，導(dǎo)致其預(yù)示的結(jié)果誤差很大，在來流空氣總溫Ta＜450K的工況中，誤差超過30%。

本文在文獻(xiàn)[5]的基礎(chǔ)上，采用回流區(qū)內(nèi)平均燃速公式，結(jié)合控制固體燃料燃速的阿雷尼烏斯公式來修正燃料內(nèi)表面溫度。建立了固體燃料沖壓發(fā)動(dòng)機(jī)燃速模型，編寫燃面預(yù)示程序，進(jìn)行了燃速計(jì)算并與已有實(shí)驗(yàn)進(jìn)行對比，驗(yàn)證了該方法的準(zhǔn)確性；首次通過數(shù)值計(jì)算發(fā)現(xiàn)工作過程中固體燃料內(nèi)腔體積隨時(shí)間呈線性關(guān)系變化并做了理論分析。

1 燃速模型

固體燃料沖壓發(fā)動(dòng)機(jī)的物理模型如圖1所示。

圖1 固體燃料沖壓發(fā)動(dòng)機(jī)燃燒室流場示意圖

實(shí)驗(yàn)研究表明，在低壓下，固體燃料沖壓發(fā)動(dòng)機(jī)中流場與固體裝藥表面間的輻射換熱量占總換熱量不足5%。忽略輻射，在燃料表面上，存在能量平衡：

式中：h為對流換熱系數(shù)，Tf為火焰溫度，Ts為燃料表面溫度，ρf為固體燃料密度，為局部燃速，hv為單位質(zhì)量固體燃料的有效汽化熱。

基于湍流邊界層傳熱理論對對流換熱系數(shù)進(jìn)行求解，文獻(xiàn)[5]建立了固體燃料沖壓發(fā)動(dòng)機(jī)平均燃速計(jì)算公式：

如圖1所示，由于突擴(kuò)臺階的作用，固體燃料沖壓發(fā)動(dòng)機(jī)燃燒室內(nèi)流場存在流動(dòng)的分離與再附著，以再附著點(diǎn)為界可將流場分為2個(gè)區(qū)域：燃燒室入口至再附著點(diǎn)為回流區(qū)，該區(qū)域內(nèi)渦旋流動(dòng)為主要流動(dòng)；再附著點(diǎn)以后為湍流再發(fā)展區(qū)，該區(qū)域內(nèi)湍流附面層不斷增厚；回流區(qū)與湍流再發(fā)展區(qū)的換熱量不能按同一種方法計(jì)算。式（2）是基于湍流邊界層傳熱理論而得的，故只能用于湍流再發(fā)展區(qū)的燃速計(jì)算，而不適用于回流區(qū)。

式中：Lr為回流區(qū)長度。

由此可見，當(dāng)其他參數(shù)不變時(shí)，再發(fā)展區(qū)中，平均退移速率正比于（Ld＝L－Lr，為再發(fā)展區(qū)長度）。

回流區(qū)內(nèi)燃速的平均值可表示為[6]

式中：a為與燃?xì)庹承浴⑷剂掀療嵊嘘P(guān)的比例系數(shù)；Ta為來流空氣總溫；Dp，Din分別為裝藥通道直徑及燃燒室入口直徑。

試驗(yàn)研究發(fā)現(xiàn)，Lr與燃燒室入口臺階高度H（H＝（Dp－Din）/2）有關(guān)，無量綱化后的回流區(qū)長度為

至此，燃燒室內(nèi)固體燃料瞬時(shí)平均燃速可表示為

對于給定的氧化劑/燃料混合物，其燃燒狀態(tài)的影響因素有：環(huán)境壓力，氧化劑/燃料的質(zhì)量比，氧化劑/燃料的初溫等。而最大的影響因素為氧化劑/燃料的質(zhì)量比：

式中：MR為空燃比，Ab為燃料的燃面面積。且有

式中：A，E分別為燃料退移的前指因子、活化能；Ru為通用氣體常數(shù)。

至此，假定固體燃料壁面溫度Ts，聯(lián)立方程（3）～方程（8）便可求出燃速，將其代入方程（9）中，反解出對應(yīng)的壁面溫度，判斷收斂性。未收斂則修正壁面溫度，再次求解；收斂則更新相關(guān)參數(shù)，進(jìn)入下一個(gè)時(shí)間步求解。

同時(shí)，對于固體燃料燃速，存在關(guān)系式：

2 燃速計(jì)算

本文假設(shè) Polymethy/methacrylate（PMMA）、Polyethylene（PE）的熱解產(chǎn)物分別為單體 MMA（C5H8O2）[1]和 C2H4。

2.1 影響燃速的參數(shù)選取及計(jì)算流程

國內(nèi)外大量文獻(xiàn)表明，固體藥柱的燃速受藥柱通道直徑，燃燒室入口直徑，燃燒室壓力，來流空氣質(zhì)量流率、總溫等因素影響。國內(nèi)外研究結(jié)果表明，在低壓下（燃燒室壓力pc＜0.8MPa），燃燒室壓力變化對燃 速 的 影 響 很 小[1，4]；本 文 所 計(jì) 算 的 工 況 均在0.8MPa以下，因此僅考慮藥柱通道直徑，燃燒室入口直徑，來流空氣質(zhì)量流率、總溫對燃速的影響。

計(jì)算中取PE的有效汽化熱、密度、前指因子及活化能分別為5.5 MJ/kg[2]，965kg/m3[8]，640 m/s和133.66kJ/mol[9]；PMMA 的有效汽化熱、密度、前指因子及活化能分別為2.0 MJ/kg，1 180kg/m3[8]，0.072 1m/s和53.086kJ/mol[9]。固體燃料沖壓發(fā)動(dòng)機(jī)中常用由溫升定義的燃燒效率ηΔT來表征固體燃料的燃燒效率。成功點(diǎn)燃的固態(tài)碳?xì)淙剂希铅平均值約為0.9[1，4]，本文取該值作為 PMMA、PE的燃燒效率。計(jì)算步長為0.005～0.05s。收斂準(zhǔn)則為：固體燃料質(zhì)量消耗率的初判值與二次判斷值的相對誤差小于0.5%，且壁面溫度相對誤差小于1%。

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

本文對PE、PMMA為燃料的固體燃料沖壓發(fā)動(dòng)機(jī)燃速共進(jìn)行了34次計(jì)算，結(jié)果如表1所示。表中試驗(yàn)數(shù)據(jù)的來源為燃速公式一欄中所對應(yīng)的參考文獻(xiàn)，燃速公式為國外進(jìn)行單一因素試驗(yàn)所得結(jié)果為實(shí)驗(yàn)所得的平均燃速為本文計(jì)算所得的平均燃速，Dp0為裝藥初始內(nèi)徑，Δ為誤差。其中，文獻(xiàn)[3]擬合的公式，燃速與壓力pc有關(guān)，計(jì)算其退移速率的試驗(yàn)值時(shí)，取pc＝0.7MPa（本文計(jì)算產(chǎn)物物性時(shí)選取的壓力）。

表1 計(jì)算結(jié)果

由表1可知計(jì)算值與實(shí)驗(yàn)值的誤差在－14.4%～13.5%之間；當(dāng)燃料為PE時(shí)，誤差絕對值的平均值為4.94%；燃料為PMMA時(shí)，誤差絕對值的平均值為8.89%。

由1＃～4＃，9＃的結(jié)果可知，燃速隨著來流空氣總溫增大而增大。這是因?yàn)殡S著來流空氣總溫增大，燃燒產(chǎn)物的溫度將變高，流場與固體裝藥間的換熱量增大，使得燃速增大。

由5＃，6＃，10＃的結(jié)果可知，燃速隨著來流空氣質(zhì)量流率的增大而增大。來流空氣質(zhì)量流率的增加，使得發(fā)動(dòng)機(jī)內(nèi)流動(dòng)速度增大，對流換熱系數(shù)變大，導(dǎo)致燃速變大。

由7＃，8＃的結(jié)果可知，燃速隨著入口直徑的減小而增大。來流空氣質(zhì)量流率相同時(shí)，入口直徑越小，燃燒室內(nèi)主流區(qū)的雷諾數(shù)、努賽爾數(shù)、對流換熱量越大[10]，使得燃速變大。

對比9＃和10＃，來流空氣總溫均為300K的條件下，9＃的裝藥內(nèi)徑大于10＃，來流空氣質(zhì)量流率、突擴(kuò)臺階高度的相對量均大于10＃中的3個(gè)工況，燃速小于10＃的3個(gè)工況，說明燃速隨著裝藥內(nèi)徑的增大而減小。內(nèi)徑的增大使得裝藥通道內(nèi)的流動(dòng)速度降低，對流換熱系數(shù)減小。

以上結(jié)論均和已有試驗(yàn)結(jié)果相符，證明了本文所建立的燃速計(jì)算模型是準(zhǔn)確可靠的。

圖2為工作過程中藥柱內(nèi)腔體積V隨時(shí)間變化的對比圖，其中試驗(yàn)數(shù)據(jù)取自文獻(xiàn)[7]。

圖2 燃料內(nèi)腔的瞬時(shí)體積

通過線性擬合發(fā)現(xiàn)，34次計(jì)算所得的藥柱內(nèi)腔瞬時(shí)體積V均為

式中：b，c為系數(shù)；t為時(shí)間。

采用上述標(biāo)準(zhǔn)，所得結(jié)果的相關(guān)性最小值為0.998 4。這表明，工作過程中燃料的質(zhì)量消耗率、空燃比不變，由此得出，火焰溫度不變，產(chǎn)物的物性不變。

由式（10）可得：

為解釋這一現(xiàn)象，本文做如下假設(shè)：①工作過程中火焰溫度及產(chǎn)物的物性不變，由文獻(xiàn)[11]的試驗(yàn)可得，工作過程中燃料質(zhì)量消耗率、空燃比不變，同時(shí)來流空氣總溫不變，在此前提下進(jìn)行這一假設(shè)是合理的；②忽略工作過程中燃料壁面溫度的微小變化，以圖3中2個(gè)工況為例，計(jì)算結(jié)果表明，工作過程中壁面溫度變化值/平均值分別為80 K/1 250K，100 K/1 380K；③不考慮再附著點(diǎn)的移動(dòng)，圖3中a＝8.5g/s工況試驗(yàn)過程中，再附著點(diǎn)位置僅變化8mm/70 mm（再附著點(diǎn)位置變化量/藥柱長度）；④不考慮燃料壁面質(zhì)量加入引起同一時(shí)刻不同截面質(zhì)量通量的變化，工作過程中，與來流空氣相比，固體燃料的質(zhì)量流率較低，約在5%～8%之間。

Heaton等人的研究結(jié)果表明，管內(nèi)流動(dòng)入口區(qū)不同軸向位置的Nu僅與Pr及無量綱特征長度x/（DReDPr）相關(guān)[11]，其中，D為圓管內(nèi)徑，ReD為按內(nèi)徑所定義的雷諾數(shù)。如圖3所示，本文涉及的問題對應(yīng)于圖中的混合入口段問題，解依賴于x/（DReDPr）和Pr，圖中僅給出了Pr＝0.7的解。

圖3 圓管內(nèi)流入口段努賽爾數(shù)的解[11]

式中：u為流速。

由圖3可得，裝藥通道再發(fā)展區(qū)內(nèi)努賽爾數(shù)Nu可表示為2個(gè)函數(shù)的乘積：

由假設(shè)④可得：

由方程（1）可得再發(fā)展區(qū)內(nèi)不同點(diǎn)的燃速：

式中：C＝（Tf－Ts）/（ρfhv），k為導(dǎo)熱系數(shù)。

由假設(shè)①②得，工作過程中，式（15）中C為常數(shù)。

結(jié)合假設(shè)①及式（14）、式（15）可得，工作過程中再發(fā)展區(qū)內(nèi)各點(diǎn)的燃速與內(nèi)徑Dp成反比，故有

至此，式（10）～式（12）所對應(yīng)的現(xiàn)象得到了解釋。

表2為同一工況下不同裝藥長度的燃速計(jì)算結(jié)果。對其進(jìn)行指數(shù)擬合，結(jié)果為

這與式（12）理論推導(dǎo)的結(jié)論相符，造成差別的原因?yàn)椋寒?dāng)只改變裝藥長度時(shí)，回流區(qū)長度Lr不變，導(dǎo)致再發(fā)展區(qū)長度Ld的相對增加量大于裝藥長度L的相對增加量，平均燃速減小量比L－0.2稍大。由此可見，長度對燃速的影響不可忽略。

3 結(jié)論

本文建立了固體燃料沖壓發(fā)動(dòng)機(jī)的燃速預(yù)示模型。編寫了燃速預(yù)示程序，以PE、PMMA為固體燃料進(jìn)行了34次的計(jì)算驗(yàn)證，主要結(jié)論如下：

①取燃燒效率ηΔT＝0.9，計(jì)算值與試驗(yàn)值吻合很好，誤差在±15%以內(nèi)，誤差絕對值的平均值為5.87%；

②固體燃料的燃面退移速率隨來流空氣總溫、質(zhì)量流率、突擴(kuò)臺階高度增大而增大，隨通道內(nèi)徑增大而減小，這些結(jié)論與國內(nèi)外實(shí)驗(yàn)相同；

③34次計(jì)算均表明，工作過程中固體燃料內(nèi)腔體積隨時(shí)間呈線性關(guān)系，這是繼文獻(xiàn)[7]的實(shí)驗(yàn)之后首次通過數(shù)值計(jì)算發(fā)現(xiàn)這一現(xiàn)象，并對形成此現(xiàn)象的原因做了分析；

④在給定工況下的計(jì)算結(jié)果表明，固體燃料沖壓發(fā)動(dòng)機(jī)中燃速正比于L－0.21，這是首個(gè)將工作過程中突擴(kuò)臺階高度、再發(fā)展區(qū)長度變化考慮在內(nèi)的固體燃料沖壓發(fā)動(dòng)機(jī)燃速模型。

[1]夏強(qiáng).固體燃料沖壓發(fā)動(dòng)機(jī)工作過程研究[D].南京：南京理工大學(xué)，2011：81－108.XIA Qiang.Investigation on the working process of solid fuel ramjet[D].Nanjing：Nanjing University of Science and Technology，2011：81－108.（in Chinese）

[2]ELANDS P J M，KORTING P A O G.Combustion of polyethylene in a solid fuel ramjet：a comparison of computational and experimental result，AIAA88－3043[R].1988.

[3]SCHULTE G.Fuel regression and flame stabilization studies of solid fuel ramjets[J].Journal of Propulsion and Power，1986，2（4）：301－304.

[4]KORTING P A O G，GELD C W M.Combustion of PMMA in a solid fuel ramjet，AIAA86－1401[R].1986.

[5]謝愛元，武曉松，夏強(qiáng).基于焓平衡法的固體燃料沖壓發(fā)動(dòng)機(jī)燃面退移速率模型[J].固體火箭技術(shù)，2012，35（5）：613－618.XIE Ai－yuan，WU Xiao－song，XIA Qiang.Prediction model of propellant regression in SFRJ based on enthalpy－balance[J].Journal of Solid Rocket Technology，2012，35（5）：613－618.（in Chinese）

[6]KRISHNAN S，GEORGE P.Solid fuel ramjet combustor design[J].Aerospace Science，1998，34（4）：219－256.

[7]ZVULONI R，GANY A，LEVY Y.Geometric effects on the combustion in solid fuel ramjets[J].Journal of Propulsion and Power，1989，1（5）：32－37.

[8]TIMNAT Y M，KORTING P A O G.Hybrid rocket motor experiments[R].Delft：Department of Aerospace Engineering，Delft University of Technology，1985.

[9]WILDE J P D.Fuel pyrolysis effects on hybrid rocket and solid fuel ramjet combustor performance[M].Delft：Delft University Press，1991：68－78.

[10]NETZER A，GANY A.Burning and flame holding characteristics of a miniature solid fuel ramjet combustor[J].Journal of Propulsion and Power，1991，7（3）：357－363.

[11]INCROPERA F P，DEWITT D P.傳熱和傳質(zhì)基本原理[M].葛新石，葉宏譯.北京：化學(xué)工業(yè)出版社，2009：317.INCROPERA F P，DEWITT D P，BERGMAN T L，et al.Fundamentals of heat and mass tranfer[M].Beijing：Chemical Industry Press，2009：317.（in Chinese）