楊 樂，余貞勇，何景軒

(中國航天科技集團(tuán)公司第四研究院第四十一研究所，固體火箭發(fā)動機(jī)燃燒、熱結(jié)構(gòu)與內(nèi)流場國防科技重點實驗室，西安 710025)

0 引言

固體火箭發(fā)動機(jī)點火瞬態(tài)過程是一個機(jī)理復(fù)雜的非定常過程，包括藥柱加熱、局部點燃、火焰沿藥柱表面?zhèn)鞑ズ腿紵以鰤旱冗^程。點火瞬態(tài)過程發(fā)生異常，如出現(xiàn)過大點火延遲、過高壓強(qiáng)峰值等，會導(dǎo)致發(fā)動機(jī)工作失效，甚至可能造成嚴(yán)重后果。所有這些現(xiàn)象與發(fā)動機(jī)幾何形狀、推進(jìn)劑物性參數(shù)、點火器及其位置等有密切關(guān)聯(lián)。定量預(yù)示發(fā)動機(jī)點火瞬態(tài)過程，揭示點火瞬態(tài)的內(nèi)流場特性，有利于對發(fā)動機(jī)點火效果的影響因素進(jìn)行分析，便于調(diào)整參數(shù)和改進(jìn)設(shè)計，從而提高發(fā)動機(jī)質(zhì)量比[1]和安全性。

國內(nèi)外許多學(xué)者[2－4]都對發(fā)動機(jī)的內(nèi)流場進(jìn)行了大量的數(shù)值分析，但對于在FLUENT中影響仿真結(jié)果的具體因素，則較少進(jìn)行深入研究。因此，本文選取了影響點火瞬態(tài)內(nèi)流場仿真的主要因素進(jìn)行分析比較。首先，由于實際點火器結(jié)構(gòu)較為復(fù)雜，需要對其進(jìn)行簡化處理，通常采用等效面積法將側(cè)向噴孔簡化成環(huán)形縫，環(huán)形縫位置不同將會對仿真結(jié)果有不同的影響;其次，在FLUENT軟件設(shè)置中，對于相關(guān)參數(shù)的輸入可采用常值簡化處理，也可利用UDF接口給出更接近實際情形的取值，不同的處理方式也會影響仿真結(jié)果;最后，對于點燃方式的選取，通常認(rèn)為推進(jìn)表面溫度達(dá)到點火溫度即為推進(jìn)劑點燃，這使得仿真結(jié)果所得到的火焰?zhèn)鞑ニ俣韧笥谠囼炛?，從而?dǎo)致最終模擬結(jié)果的精確程度出現(xiàn)較大偏差。鑒于此，采用更為合理的點燃方式則是解決此類問題的有效途徑。

本文采用三維物理模型，模擬大長徑比固體火箭發(fā)動機(jī)不同情形的點火升壓過程。同時，該類固體火箭發(fā)動機(jī)的點火瞬態(tài)過程壓強(qiáng)、推力等變化劇烈，火焰?zhèn)鞑ミ^程較長，對其仿真能更好反映出上述點燃方式對模擬結(jié)果的影響。

1 計算模型

1.1 發(fā)動機(jī)內(nèi)流場結(jié)構(gòu)

計算模型如圖1所示，考慮翼槽的對稱性，取1/8進(jìn)行模擬仿真。

圖1 流場結(jié)構(gòu)Fig.1 Flow field structrue

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

為便于計算，對模型作如下假設(shè):

(1)由于點火過程非常短暫，在計算中不考慮發(fā)動機(jī)內(nèi)流場與結(jié)構(gòu)之間的耦合關(guān)系;

(2)燃燒所生成的混合氣體為理想氣體;

(3)點火瞬態(tài)不計侵蝕燃燒，推進(jìn)劑燃速只與當(dāng)?shù)貕簭?qiáng)有關(guān);

(4)點火器燃?xì)馀c推進(jìn)劑燃?xì)饩哂邢嗤再|(zhì)，忽略各組分間化學(xué)反應(yīng)，忽略比定壓熱容隨溫度的變化，取為常數(shù);

(5)采用動態(tài)溫度點火方式，即相對于通常取燃面附近流體單元溫度達(dá)到恒定點火溫度作為點燃判據(jù)而引入的定義。其以固體推進(jìn)劑薄層內(nèi)部表面達(dá)到點火溫度為依據(jù)，在考慮燃面與流場的對流換熱特性及推進(jìn)劑燃燒時自身的化學(xué)反應(yīng)過程中吸、放熱等各個因素共同作用的基礎(chǔ)上，引入隨點燃處推進(jìn)劑燃速變化的轉(zhuǎn)換因子將此恒定的點火溫度轉(zhuǎn)化為在數(shù)值模擬中容易獲取的燃面附近流體單元的動態(tài)點火溫度，進(jìn)而將其作為推進(jìn)劑的點燃判據(jù)。

1.3 控制方程

流場燃?xì)獠捎梅嵌ǔ？蓧嚎sN-S方程，以連續(xù)、動量和能量方程為基礎(chǔ)，考慮氣體粘性和熱擴(kuò)散率隨溫度的變化關(guān)系;湍流模型采用RNGk-ε兩方程模型，該模型可更好地處理高應(yīng)變率及流線彎曲程度較大的流動;近壁區(qū)由于Re數(shù)較低，湍流發(fā)展并不充分，湍流的脈動影響不如分子粘性大，所以采用標(biāo)準(zhǔn)壁面函數(shù)法處理。

為考察脈動影響，目前對湍流廣泛采用的方法是時間平均法。時均化的N-S方程通用表達(dá)式[5]為

式中 φ為通用變量(φ=1，u，v，T分別對應(yīng)連續(xù)、動量和能量方程);Γφ為與φ相對應(yīng)的廣義擴(kuò)散系數(shù);湍流脈動能附加項是與φ相對應(yīng)的湍流擴(kuò)散系數(shù));S為廣義源項，包含有兩項間的相互作用。

1.4 初始及邊界條件

初始條件:該發(fā)動機(jī)為地面靜止試驗發(fā)動機(jī)，且燃燒室初始沖入壓強(qiáng)為0.1 MPa，環(huán)境壓強(qiáng)為0.094 8 MPa，因此整個流場區(qū)域的初始狀態(tài)取為:T=300 K，p=0.194 8 MPa，3個方向的初始速度為零。

邊界條件:

(1)點火器出口采用質(zhì)量流率邊界條件;

(2)發(fā)動機(jī)頭部為絕熱邊界;

(3)藥柱表面為熱耦合邊界，即燃面點燃前按加熱表面處理;點燃之后按側(cè)壁加質(zhì)邊界進(jìn)行處理;

(4)噴管出口在堵蓋打開后設(shè)為壓力出口邊界條件;堵蓋在打開前其為固體壁面邊界條件，打開后為內(nèi)部邊界條件。

1.5 計算方法

文中利用UDF接口編程進(jìn)行二次開發(fā)[6－7]，采用側(cè)壁加質(zhì)的方法設(shè)定固體推進(jìn)劑燃面邊界，點火判據(jù)采用動態(tài)溫度點火方式。由于發(fā)動機(jī)燃燒室在點火瞬態(tài)期間壓強(qiáng)、溫度迅速上升，會造成燃?xì)鉄釋?dǎo)率遠(yuǎn)大于固體推進(jìn)劑的熱導(dǎo)率[8]。本文采用的轉(zhuǎn)換因子的變化規(guī)律為與點燃處推進(jìn)劑燃速呈線性關(guān)系，隨燃速增加而減小，變化范圍約為0.96～0.85。根據(jù)點火器的壓強(qiáng)-時間曲線，在保證總質(zhì)量流量的前提下對點火器燃?xì)饬髁窟M(jìn)行了假定。

計算過程以堵蓋打開為分界線，將整個工作過程分為兩部分，堵蓋打開前設(shè)為壁面邊界，打開后設(shè)為內(nèi)部邊界條件，同時將之前的計算結(jié)果作為初始條件進(jìn)行計算，直到發(fā)動機(jī)點火瞬態(tài)過程完成。取1/8整翼進(jìn)行計算，由于發(fā)動機(jī)翼槽呈周向均勻分布，此簡化方式可較準(zhǔn)確地接近實際情形。利用有限體積法建立離散方程，壓強(qiáng)和密度耦合方式在堵蓋打開前選用PISO算法，堵蓋打開后選用SIMPLE算法，守恒方程采用一階迎風(fēng)格式，時間步長設(shè)置為1×10－5s。

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

在以上物理模型的基礎(chǔ)上，選取了4種組合方式進(jìn)行分析，如表1所示。所模擬的燃燒室頭部壓強(qiáng)和試驗所測壓強(qiáng)的結(jié)果如圖2所示。其中由于點火器側(cè)向噴孔呈周向非均勻分布，三維建模較為復(fù)雜，且劃分網(wǎng)格時不容易生成高質(zhì)量的三維網(wǎng)格，則采用等效面積法對其進(jìn)行2種方式的簡化。

2.1 不同點火器簡化方式的影響

對比表1中方式1和方式2，可看到前者升壓速率明顯高于后者，表明其火焰?zhèn)鞑ニ俣容^快，在0.12 s時刻已經(jīng)接近穩(wěn)態(tài)工作壓強(qiáng)，與試驗值不一致。這是因為初始階段采用方式1簡化處理的點火器會使引燃主藥柱面積較大，從而在隨后的升壓過程中，其累積效應(yīng)更加明顯，使所模擬的升壓過程中不同時刻的壓強(qiáng)值均高于方式2。圖3為2種方式下4、50 ms時刻的燃燒室頭部溫度云圖。由圖3(a)、(b)可見，2種方式下引燃的主藥柱初始點火位置明顯不同，方式1的初始點燃位置相對方式2要遠(yuǎn)離頭部，這是因為雖然燃?xì)馔瑯釉陬^部發(fā)生聚集，但在進(jìn)口燃?xì)獾臎_擊點，溫度升高會更快速，從而在正對點火器出口處的單元會先被點燃。另外，圖3(c)、(d)也顯示出方式1火焰?zhèn)鞑ニ俣纫哂诜绞?。

表1 參數(shù)組合方式Table 1 Method of parameter combination

圖2 計算值與試驗值對比Fig.2 Comparison between calculation and experiment value

圖3 不同時刻溫度云圖對比Fig.3 Comparison of temperature contour at different times

2.2 不同點火器出口參數(shù)選取方式的影響

對比表1中方式2和方式3，二者的升壓速率相差不大。整個模擬過程以0.02 s時刻為分界點，之前二者壓強(qiáng)相差較大，之后差距則有所降低。這是因為初始階段引起燃燒室壓強(qiáng)升高的因素中，點火器燃?xì)獾募尤胝紦?jù)主導(dǎo)作用。方式3以點火器出口燃?xì)獾钠骄导尤?，顯然在達(dá)到點火器峰值時間之前所加入的燃?xì)饪傎|(zhì)量相對較多，導(dǎo)致模擬的升壓速率會大于前者，隨著時間的推進(jìn)，方式2在點火器峰值期間引燃主藥柱面積會比后者顯著增加，升壓速率大于后者，從而在點火器峰值附近時刻二者燃燒室壓強(qiáng)相差會很小。此后，燃燒室壓強(qiáng)的增加主要由主藥柱的加質(zhì)來決定，點火器的貢獻(xiàn)相對則很小，且二者的其余邊界條件相同，升壓規(guī)律基本相同。圖4為3、45 ms時刻的溫度云圖。由圖4(a)、(b)可見，方式3的初始引燃面積明顯大于方式2。另外，圖4(c)、(d)顯示出二者的主藥柱點燃面積已經(jīng)相差很小。

2.3 不同點燃方式的影響

對比表1中方式2和方式4，方式4的升壓速率相對較高，其達(dá)到穩(wěn)態(tài)工作壓強(qiáng)的時間也必然早于方式2。其原因在于點火瞬態(tài)期間，氣體溫度的升高速率會明顯大于固體推進(jìn)劑表面溫度的升高速率，方式4以燃面附近氣體單元溫度達(dá)到點火溫度為依據(jù)，會使所模擬的火焰?zhèn)鞑ニ俣雀哂趯嶋H值。而方式2考慮了對流傳熱及各種影響因素對火焰?zhèn)鞑サ挠绊懀x擇固體推進(jìn)劑薄層內(nèi)單元溫度達(dá)到點火溫度作為點燃判據(jù)，同時利用轉(zhuǎn)換因子將此恒定的點火溫度轉(zhuǎn)化為燃面附近流體單元的動態(tài)點火溫度，以滿足與實際點火過程的流場狀態(tài)更加吻合，因而模擬的點火瞬態(tài)過程的壓強(qiáng)值與試驗值基本一致。圖5為40 ms時刻的溫度云圖。由圖5可見，方式2的主藥柱點燃面積要小于方式方式4，反映出方式2相對方式4有效地降低了火焰?zhèn)鞑ニ俣?，使得仿真結(jié)果更為準(zhǔn)確。

圖4 不同時刻溫度云圖對比Fig.4 Com parison of tem perature contour at different times

圖5 40 m s時刻溫度云圖對比Fig.5 Com parison of tem perature contour at 40 m s

3 結(jié)論

(1)點火器側(cè)向噴孔沿軸向簡化為多條環(huán)形縫結(jié)構(gòu)時，主藥柱初始點火位置更靠近頭部，仿真值與試驗值基本吻合。

(2)流量入口邊界條件中各數(shù)值采用平均值和采用實測值在點火器壓強(qiáng)達(dá)到峰值附近時刻之前升壓規(guī)律不同，其后則基本一致。

(3)燃面加質(zhì)邊界條件中點火溫度判據(jù)采用恒定值會導(dǎo)致仿真的火焰?zhèn)鞑ニ俣让黠@加快，而考慮了傳熱影響的動態(tài)點火溫度判據(jù)所模擬出的結(jié)果則更好地反映了實際情形。

[1] 方丁酉，等.固體火箭發(fā)動機(jī)內(nèi)彈道學(xué)[M].長沙:國防科技大學(xué)出版社，1997.

[2] 劉君，郭正，等.非結(jié)構(gòu)網(wǎng)絡(luò)技術(shù)應(yīng)用于固體火箭發(fā)動機(jī)數(shù)值模擬[J].固體火箭技術(shù)，2001，24(4):9-11.

[3] Cang S T，Han S，Joh C.Radiation effect ignition on I-D transient anlysis of SRM[R].AIAA 96-3055.

[4] 陳軍濤，蹇澤群，陳林泉.固體火箭發(fā)動機(jī)點火瞬態(tài)內(nèi)流場軸對稱數(shù)值分析[J].固體火箭技術(shù)，2004，27(3):183-176.

[5] 王福軍.計算流體動力學(xué)分析[M].北京:清華大學(xué)出版社，2004.

[6] 周志清，高雙武.含裝藥缺陷固體火箭發(fā)動機(jī)點火瞬態(tài)內(nèi)流場的數(shù)值研究[D].第二炮兵工程學(xué)院，2006.

[7] 孫淑霞，肖陽春，等.C/C++程序設(shè)計教程[M].北京:電子工業(yè)出版社，2009.

[8] 匹茨D R，西遜姆L E.傳熱學(xué)的理論和習(xí)題[M].上海:機(jī)械工業(yè)出版社，1983.