郁安吉,薛曉春,余永剛
(南京理工大學(xué) 能源與動力工程學(xué)院,江蘇 南京 210094)
裝藥溫度是影響箭炮武器性能的重要因素,藥溫測量的精度將直接影響火箭彈道解算結(jié)果的正確性,進(jìn)而影響武器性能[1-2]。因此,準(zhǔn)確獲得火箭發(fā)動機內(nèi)部熱環(huán)境參數(shù)是對其熱防護(hù)系統(tǒng)進(jìn)行精細(xì)化設(shè)計的前提[3]。
通??刹捎脗鳠嵯嗨品ê徒馑惴▋煞N方法進(jìn)行裝藥溫度的測量。解算法是利用傳熱學(xué)基本原理,通過測量環(huán)境溫度,對裝藥溫度場進(jìn)行解算[4];傳熱相似法是采用縮比模型,通過相似原理和無量綱分析,保證其傳熱過程的微分方程和單值條件相同[5],可滿足各類火箭實時彈道解算的要求。劉慶才等[6]利用解算法對固體推進(jìn)劑建立了導(dǎo)熱模型,使用非接觸式試驗裝置對推進(jìn)劑的溫度場進(jìn)行了測定,得到的計算結(jié)果與實驗測點溫度值吻合較好。周彥煌等[7]利用全自動間接測量法,通過微機控制的非接觸式實時測量裝置,設(shè)計出的裝藥溫度自動測量系統(tǒng)實現(xiàn)了藥溫自動實時測量,滿足當(dāng)代作戰(zhàn)需要。
目前大口徑火箭發(fā)動機達(dá)到平衡溫度所需時間過長,增加了試驗的能耗,影響其研制效率,而裝藥溫度場的變化特性是確定該時間的前提[8-9]。因此,科研人員在大口徑武器彈藥保溫方面進(jìn)行了深入研究[10]。李媛等[11]為研究環(huán)境溫度對固體火箭發(fā)動機藥柱溫度場的影響,利用有限元傳熱進(jìn)行分析。薛青等[12]根據(jù)發(fā)射藥保溫的原理,通過仿真與試驗相結(jié)合的方法,驗證發(fā)射裝藥保溫過程彈藥溫度接近環(huán)境溫度需要的時間遠(yuǎn)低于現(xiàn)有標(biāo)準(zhǔn)所規(guī)定的時間。
針對火箭發(fā)動機裝藥溫度場平衡溫度所需時間的計算問題,本文采取數(shù)值傳熱學(xué)理論與實驗相結(jié)合的方法,建立了火箭發(fā)動機的傳熱數(shù)學(xué)模型,針對不同環(huán)境溫度條件下的某火箭發(fā)動機的溫度場進(jìn)行試驗測試;通過對比數(shù)值計算與實驗結(jié)果,驗證了所建立的多層材料間傳熱模型的合理性;最后深入分析了火箭發(fā)動機各層材料的溫度變化特性。
由于火箭發(fā)動機是由金屬殼體、絕緣層和推進(jìn)劑裝藥等多層材料組成,試驗中在裝藥上選取多個具有代表性的測點,安裝完傳感器后將實驗?zāi)P头湃氡叵?進(jìn)行保溫測試。因此,為了模擬這種多層材料在高溫或者低溫環(huán)境中的溫度場特性,建立如圖1所示的物理模型,并取a、b、c三點作為特征點記錄其溫度變化特性,其中a點位于裝藥1/3位置,b點位于裝藥1/2位置,c點位于裝藥2/3位置。針對該物理模型,提出以下假設(shè):

圖1 某火箭發(fā)動機物理模型Fig.1 Physical model of a rocket engine
①由于發(fā)動機溫度沿徑向變化較大,因此采用一維柱坐標(biāo)描述這種多層材料間的非穩(wěn)態(tài)傳熱現(xiàn)象;
②各層材料間無熱阻;
③空氣的導(dǎo)熱系數(shù)和比熱只是溫度的函數(shù);
④發(fā)動機各層材料的物性參數(shù)為常數(shù)。
針對上述假設(shè),建立這種多層材料間的傳熱模型:
(1)
式中:T為發(fā)射藥溫度;r為發(fā)動機截面徑向坐標(biāo);λi、ρi、ci分別是材料(i=1時為推進(jìn)劑、i=2時為絕緣層、i=3時為發(fā)動機殼體)的導(dǎo)熱系數(shù)、密度和比熱容。
初始時刻求解區(qū)域內(nèi)各點的溫度相同,即初始條件為
T(t,r)|t=0=f(0,r)
(2)
當(dāng)火箭發(fā)動機處在高溫或低溫環(huán)境中時,受外界熱源或者冷源的對流傳熱作用,其內(nèi)各層材料的溫度均會發(fā)生變化,并不斷接近環(huán)境溫度。因此,采用第三類邊界條件作為發(fā)動機內(nèi)孔內(nèi)邊界條件和發(fā)動機殼體外邊界條件:
(3)
(4)
式中:Rs、Rn分別為發(fā)動機內(nèi)孔內(nèi)徑和發(fā)動機殼體半徑;Tws為發(fā)射藥內(nèi)邊界溫度;Twn為發(fā)動機殼體外邊界溫度;Tfs和Tfn分別為內(nèi)環(huán)境和外環(huán)境溫度。hs和hn分別為發(fā)動機內(nèi)外邊界的對流換熱系數(shù),通過努塞爾數(shù)Nu計算。如下所示:
Nu=C(GrPr)b
(5)
(6)
式中:C,b為實驗系數(shù);Pr=0.7;αV為空氣體積膨脹系數(shù);g為重力加速度;ΔT為壁面與環(huán)境的溫差;δ為特征長度;μ為空氣動力黏度。
采用控制容積法離散傳熱方程,求解區(qū)域離散采用內(nèi)節(jié)點法,將推進(jìn)劑、絕緣層、殼體和空氣層進(jìn)行網(wǎng)格劃分,并在時間間隔[t,t+Δt]內(nèi)將非穩(wěn)態(tài)導(dǎo)熱方程對控制體作積分:
(7)
假定方程非穩(wěn)態(tài)項中溫度隨時間和空間都是階梯性變化,擴散項隨空間分段線性變化,隨時間階梯性變化,則有:
(8)
式中:下標(biāo)P,S,N表示求解區(qū)域相應(yīng)網(wǎng)格節(jié)點;半徑方向相鄰兩節(jié)點間的距離以δr表示,下標(biāo)s,n表示相應(yīng)控制容積界面;上標(biāo)0表示上一時刻的值。離散后方程可寫為
(9)

為使方程(9)封閉,采用附加源項法進(jìn)一步處理邊界條件,即考慮所有節(jié)點上離散方程系數(shù)形式的一致性,對于內(nèi)外邊界而言,將方程的系數(shù)aP和b改寫如下:
(10)
式中:SC,SP為一般傳熱控制方程形式中源項作局部線性化處理后的系數(shù)。由于裝藥與環(huán)境溫度之間傳熱是無源的,即:SC=0,SP=0,其余參數(shù)保持不變,因此改寫上述系數(shù)不影響離散后傳熱控制方程的解。
對于求解區(qū)域邊界上的節(jié)點而言,由于采用了附加源項法處理邊界條件,其源項部分即SC、SP不為0,因此上述的SC、SP為0僅滿足求解區(qū)域上的內(nèi)部節(jié)點。但由于添加了SC,ad、SP,ad部分,使得發(fā)射裝藥與環(huán)境溫度之間的傳熱關(guān)系式離散后的系數(shù)在求解區(qū)域上對所有節(jié)點具有了統(tǒng)一的形式,有利于程序的編制。
根據(jù)圖2對邊界條件進(jìn)行處理,先對位于邊界上的N節(jié)點表示的邊界條件進(jìn)行研究。

圖2 附加源項圖示說明Fig.2 Illustration of additional source items
系數(shù)aN表示為
(11)
式中:rp為邊界點N到P點的距離,λb為邊界節(jié)點的導(dǎo)熱系數(shù)。
邊界條件可改寫為
(aP-aN)TP=aN(TN-TP)+aSTS+b
(12)
利用已知的邊界條件把TN消去未知的邊界溫度,由Fourier定律得:
(13)
得:
(14)
于是,有
(15)
式中:h為表面對流換熱系數(shù);SC,ad、SP,ad為附加源項。
火箭發(fā)動機是由多層材料組成,多層材料間物性均采用調(diào)和平均方法進(jìn)行計算,各層材料(發(fā)動機殼體、絕緣層及推進(jìn)劑)的導(dǎo)熱系數(shù)、密度和比熱容,其物性參數(shù)見表1。

表1 材料物性參數(shù)Table 1 Physical property parameters of materials
為了驗證上述數(shù)值模型的合理性,本文對20 ℃環(huán)境條件下的初溫-40 ℃的某火箭發(fā)動機溫度變化特性進(jìn)行了實驗研究,溫度測點位置按圖1所示布置。發(fā)動機內(nèi)部溫度的數(shù)值計算結(jié)果與實驗測量變化曲線如圖3所示??芍?a,b,c測點處的數(shù)值計算與實驗最大誤差分別為6.72%、7.84%和7.96%,且在整個升溫過程中,通過式(16)可得數(shù)值計算與實驗的溫度誤差最大為1.423 ℃;

圖3 實驗與數(shù)值計算結(jié)果對比圖Fig.3 Comparison of experimental data and calculation results
(16)
式中:TE,k為實驗測量值,TC,k為理論計算值,n為統(tǒng)計次數(shù)。
2.2.1 升溫過程分析
為了探究升溫過程中,發(fā)動機初始溫度對其溫度場的影響,針對環(huán)境溫度為20 ℃,裝藥初始溫度分別為-30 ℃、-40 ℃、-50 ℃和-60 ℃時的發(fā)動機溫度場變化特性進(jìn)行多工況分析,結(jié)果如圖4所示。

圖4 初始溫度對發(fā)動機內(nèi)部溫度的影響Fig.4 Temperature changes at different initial temperatures
由圖4可知,當(dāng)環(huán)境溫度相同時,升溫過程中裝藥初始溫度越高,各點溫度越高,升溫速率越小。在起始的72 h內(nèi),裝藥溫度上升速度較快,各工況間溫度差較大。隨著時間推移,發(fā)動機溫差越來越小,當(dāng)t=96 h時,不同工況下a測點和外徑位置處的溫度均接近環(huán)境溫度,溫差最小分別為0.037 9 ℃和0.021 8 ℃。當(dāng)裝藥初溫增加10 ℃時,裝藥特征點達(dá)到19 ℃所需時間最多減少了2.44 h。表2和表3分別為發(fā)動機不同初溫條件下,a測點和外徑位置處各個時刻的溫升速率,可進(jìn)一步分析發(fā)動機內(nèi)溫度場的變化規(guī)律。

表2 不同初溫下a測點溫升速率Table 2 Temperature rise rate at the a measurement point with different initial temperatures

表3 不同初溫下裝藥外徑溫升速率Table 3 Temperature rise rate at charge outer diameter with different initial temperatures
由表2可知,t=12 h時,各工況的溫升速率均大于1.5 ℃/h,初始溫度越低,溫升速率越大。隨時間推移,各工況的溫升速率開始減小,t=72 h時溫升速率將至0.072 ℃/h以下,此時初始溫度為-60 ℃的溫升速率仍為最大。t=96 h時,不同初溫工況下的溫升速率均為0.018 ℃/h,此時各工況測點溫度均接近環(huán)境溫度,因此溫升速率也接近。
由表3可知,裝藥外徑的溫度變化規(guī)律與a測點相似。由于裝藥外徑處殼體材料的導(dǎo)熱系數(shù)較大,傳熱量越大,溫升速率更快。結(jié)合圖4(b)可知,在t=12 h時,裝藥外徑處4種工況下與環(huán)境最大溫差達(dá)到14.31 ℃,a測點4種工況下與環(huán)境最大溫差達(dá)到39.42 ℃,裝藥外徑處溫度在變化初期迅速上升,縮小與環(huán)境的溫差,使得溫升速率在溫度變化后期較小。
為了探究升溫過程中,環(huán)境溫度對其溫度場的影響規(guī)律。針對裝藥溫度為-40 ℃,環(huán)境溫度分別為30 ℃、20 ℃、10 ℃和5 ℃時的發(fā)動機溫度場變化特性進(jìn)行分析,結(jié)果如圖5所示。

圖5 不同環(huán)境溫度下不同測點的溫度變化Fig.5 Temperature variation at different ambient temperatures
由圖可知,當(dāng)裝藥溫度相同時,環(huán)境溫度越高,升溫過程中各點溫度越高,溫升速率越大??梢钥闯?在起始72 h內(nèi)裝藥溫度上升速度較快,各工況間溫差較小。隨著時間的推移,不同工況下的發(fā)動機內(nèi)溫度均已經(jīng)接近環(huán)境溫度。通過圖5(c)可知,當(dāng)t=2 h、環(huán)境溫度分別為10 ℃和30 ℃時,a測點的溫升速率分別為1.906 ℃/h和2.805 ℃/h,此時裝藥外徑處對應(yīng)的溫升速率分別為2.769 ℃/h和3.83 ℃/h,明顯大于a測點的溫升速率。這是由于溫度變化初期,裝藥外徑處殼體材料的導(dǎo)熱系數(shù)大于推進(jìn)劑的導(dǎo)熱系數(shù),傳熱量大,溫度上升速率快。隨著升溫過程的繼續(xù),裝藥外徑與環(huán)境間的溫度差減小,溫升速率明顯小于a測點。當(dāng)環(huán)境溫度為30 ℃時,兩處裝藥特征點達(dá)到29 ℃所需時間相差14.42 h,當(dāng)環(huán)境溫度為10 ℃時,兩處裝藥特征點達(dá)到9 ℃所需時間相差3.82 h。表4和表5分別為不同環(huán)境溫度條件下發(fā)動機內(nèi)溫度場的變化規(guī)律。

表4 不同環(huán)境溫度下a測點溫升速率Table 4 Temperature rise rate at the measurement point of a with different ambient temperatures

表5 不同環(huán)境溫度下裝藥外徑溫升速率Table 5 Temperature rise rate of charge outer diameter with different ambient temperatures
由表4可以看出,t=12 h時,溫升速率均大于1.5 ℃/h,環(huán)境溫度越高,溫升速率越快。隨試驗時間推移,各工況下溫升速率開始減小;t=72 h時,各工況下的溫升速率都小于0.054 ℃/h,此時環(huán)境溫度為30 ℃的發(fā)動機的溫升速率仍為最大;t=96 h時,不同環(huán)境溫度工況下的溫升速率均較小。
結(jié)合圖5可以看出,溫度變化初期溫差越大,溫升速率越快,經(jīng)過較長時間,各工況的裝藥溫度與環(huán)境溫度的溫差減小,溫升速率減小直至為0 ℃/h,最終與環(huán)境溫度一致。
在t=12 h時,裝藥外徑處4種工況下與環(huán)境的最大溫差為13.136 ℃,a測點4種工況下與環(huán)境的最大溫差為35.446 ℃;t=96 h時,裝藥外徑處的裝藥溫度與環(huán)境溫度最大溫差僅為0.23 ℃,此時a測點的裝藥溫度與環(huán)境溫度最大溫差為0.47 ℃。
2.2.2 溫降過程分析
為了探究降溫過程中,發(fā)動機內(nèi)裝藥初始溫度和環(huán)境溫度對其溫度場的影響規(guī)律。針對不同裝藥初始溫度不同環(huán)境溫度進(jìn)行分析,結(jié)果如圖6所示。

圖6 不同邊界條件下不同測點的溫度變化Fig.6 Temperature variation at different measurement points under different boundary conditions
由圖6可知,在起始的72 h內(nèi),裝藥溫度下降速度較快,隨著裝藥溫度的下降,兩測點處溫度變化均較為平緩。這是由于各工況下發(fā)動機內(nèi)溫度均接近環(huán)境溫度,溫差越小傳熱越慢。通過圖6(c)可知,溫度變化初期兩種工況下a測點溫度變化較為平緩,裝藥外徑溫度變化迅速。這是由于裝藥外徑的殼體導(dǎo)熱系數(shù)大,溫度下降快,溫度變化后期裝藥溫度非常接近環(huán)境溫度,使得發(fā)動機內(nèi)裝藥溫度需要較長時間與環(huán)境溫度達(dá)到完全平衡。當(dāng)環(huán)境溫度為-20 ℃時,兩處裝藥特征點達(dá)到-19 ℃所需時間相差3.38 h,當(dāng)環(huán)境溫度為-45 ℃時,兩處裝藥特征點達(dá)到-44 ℃所需時間相差3.28 h。表6和表7為不同邊界條件下,裝藥1/3和外徑位置處各時刻的溫升速率。可知,同一時刻同一測點處,裝藥初始溫度與環(huán)境溫度溫差越大,溫度變化速率越大。對比10 ℃~-20 ℃和10 ℃~-45 ℃這兩種工況下同一時刻不同測點處溫升速率可知,由于裝藥外徑處殼體材料的導(dǎo)熱系數(shù)大于推進(jìn)劑的導(dǎo)熱系數(shù),t=2 h時裝藥外徑的溫度變化速率大于a測點;隨著降溫過程的繼續(xù),當(dāng)t=96 h時,溫度變化速率接近于0 ℃/h,這是由于后期裝藥溫度非常接近環(huán)境溫度,傳熱量減少,溫升速率較慢,因此發(fā)動機內(nèi)各點溫度與環(huán)境溫度達(dá)到完全平衡需要較長時間。

表6 不同邊界條件下a測點溫度變化速率Table 6 Temperature change rates of measurement point a under different boundary conditions

表7 不同邊界條件下裝藥外徑溫度變化速率Table 7 Temperature change rate of charge outer diameter under different boundary conditions
①通過某口徑火箭發(fā)動機多層材料在不同邊界條件下裝藥溫度場與環(huán)境之間的數(shù)值計算與實驗結(jié)果的對比,驗證了所建立理論模型能夠準(zhǔn)確描述裝藥溫度場的溫度變化特性。
②由于殼體與推進(jìn)劑的導(dǎo)熱系數(shù)的不同,裝藥外徑的殼體導(dǎo)熱系數(shù)大于裝藥內(nèi)側(cè),使得溫度變化初期裝藥外側(cè)溫度變化波動較大,裝藥內(nèi)側(cè)溫度變化較為平緩,隨測試時間的增長,整個裝藥溫度場終將趨于平衡狀態(tài)。
③同一時刻同一測點,固定裝藥初溫或者環(huán)境溫度時,溫差越大,溫度變化速率也就越大;而同一時刻不同測點的溫度變化速率與材料的導(dǎo)熱系數(shù)相關(guān)。
④溫差是熱量傳遞的動力,熱量的傳遞與物體內(nèi)部溫度分布狀況密切相關(guān),多層材料發(fā)動機裝藥和環(huán)境的初始溫差越大,發(fā)動機內(nèi)各點與環(huán)境溫度達(dá)到完全平衡的時間越長。