固體火箭發動機內彈道多參數辨識及精準預示方法①

張 凌,王德友,李世鵬,王寧飛

(北京理工大學 宇航學院,北京 100081)

0 引言

固體火箭發動機內彈道學的基本任務是在發動機各種工作條件下計算燃燒室壓強隨時間和空間的變化規律。零維內彈道計算方法假定燃燒室內參數均勻一致,與一維及多維內彈道計算方法相比更加簡便、直觀,是工程上應用較多的求解方法[1-2]。發動機的內彈道特性直接決定著整機性能,進而影響著火箭導彈的飛行特性。因此,進一步提高內彈道預示的精度,對保障飛行器工作可靠性有著重要意義。

工程上通常根據單臺或多臺發動機的實測數據和燃面幾何退移規律,對燃速特性等參數進行辨識,以獲得更精確的內彈道性能預示結果。目前,業內學者針對發動機內彈道性能精確預示問題,從理論算法修正、參數辨識方法改進等角度開展了研究。如劉楊等[3]基于對發動機實測數據處理誤差的理論分析,提出了修正的燃速系數處理方法,該方法可提高具有較高壓強比的發動機內彈道預示精度。李春艷等[4]針對短時間工作固體發動機特點,從理論角度對現有內彈道計算方法進行了適應性修改,所得預示結果與實測數據相比具有較高精度。

在辨識參數選擇方面,CAVALLINI等[5-6]通過研究Vega固體運載火箭發動機實際飛行性能影響因素,發現噴管喉部燒蝕率對推力和燃燒室壓強的影響較大。BIANCHI等[7]采取數值仿真方法研究噴管喉部燒蝕率和燃速變化對燃燒室壓強的影響,發現這兩個參數的改變對壓強-時間曲線影響十分明顯。李曉斌等[8]基于發動機工作過程的內彈道數學模型及實測燃燒室壓強時變數據,通過參數辨識方法可對燃速、喉襯燒蝕率等參數進行計算。

在參數辨識方法改進方面,由于最小二乘法、最大似然法等傳統辨識方法常存在初始點敏感等問題,因而如遺傳算法等具有初值穩健性和全局搜索性的方法在內彈道參數辨識領域逐漸得到深入應用。樊超等[9]將遺傳算法用于固體火箭發動機參數辨識,獲得了推進劑燃速模型和喉徑變化模型的全局最優辨識值,所得內彈道性能預示結果與試驗測量值吻合良好。莊建華等[10]對發動機貯存老化內彈道性能的影響因素進行參數建模,基于拉丁超立方設計方法進行了瞬態參數辨識,結合非定常內彈道模型對發動機貯存老化內彈道進行了預示。蒲曉航等[11]基于部分實測壓強數據,通過遺傳算法逆向反算獲得預示參數,再結合正向計算實現了內彈道的預示自修正,對復雜環境剖面下的發動機內彈道性能進行了準確預示。孫瑞陽等[12]建立了基于經典AlexNet卷積神經網絡的內彈道參數辨識方法,對長期貯存發動機推進劑裝藥的燃速特性參數進行了辨識。張小平等[13]建立了高能固體推進劑的配方表征方法,利用遺傳神經網絡(GA-BP)對推進劑高壓燃燒性能及其變化規律進行了預示。

傳統發動機的實測數據處理常通過燃燒時間內的平均壓強和藥柱平均燃速計算獲得燃速系數等參數并反算燃面,而實際上平均壓強與平均燃速間并不符合燃速模型對應關系[3],因而計算存在一定誤差。而目前的研究工作大多采用遺傳算法和神經網絡等方法對關鍵內彈道參數進行辨識,并沒有根據實測數據對藥柱燃面以及內彈道參數進行修正與二次辨識,內彈道參數辨識以及性能預示的精度還可以進一步提高。因此,有必要對關鍵內彈道參數開展二次辨識,得到更加精確的辨識結果,進而提高發動機內彈道預示精度。

本文提出了一種發動機內彈道多參數辨識及性能精準預示方法,根據發動機實測壓強數據,通過遺傳算法對燃速系數、壓強指數、推進劑密度和噴管喉襯燒蝕速率系數開展初步辨識,基于辨識結果修正燃面-肉厚關系,在此基礎上開展二次參數辨識,最終確定關鍵內彈道參數。根據第二次參數辨識結果和修正后的燃面-肉厚關系,通過正向預示可獲得發動機的內彈道性能。以美國C1xb型發動機[14-15]和國內某型長通道發動機[16]為案例進行計算分析,驗證了該方法內彈道性能預示的精度和有效性,相關結果可為固體火箭發動機的優化設計提供技術參考。

1 內彈道預示模型

本文以零維內彈道模型[17-18]為基礎進行發動機性能正向預示,該模型有如下假設:(1)燃氣為理想氣體,氣體流動為零維定常流動;(2)裝藥點燃過程瞬時完成;(3)噴管喉襯燒蝕速率與燃燒室壓強正相關。在考慮點火藥和噴管喉襯燒蝕作用對內彈道計算影響的情況下,結合燃燒室內質量守恒定律、能量守恒定律、燃燒室內體積守恒方程、燃速方程以及氣體狀態方程,推導出的內彈道微分方程組如下:

(1)

其中,

(2)

(3)

(4)

(5)

2 多參數辨識及精準預示方法

考慮到發動機試驗成本問題,基于一次或幾次試車數據開展發動機內彈道關鍵參數的精確辨識研究,對降低固體發動機設計成本和提高設計效率有著重要意義。基于工程中的內彈道參數辨識方法,本文提出一種發動機內彈道多參數辨識及預示方法,可提高內彈道預示精度。

2.1 工程中內彈道參數辨識方法

(6)

(7)

(8)

式中W為裝藥肉厚;t為燃燒時間;Ab為燃面面積;At為噴管喉部面積;C*為特征速度。

2.2 多參數辨識及精準預示方法流程

上述平均燃速和平均壓強并不具有實際對應關系,計算存在一定誤差,在內彈道預示精度上仍有較大優化空間。本文提出了一種內彈道多參數辨識及精準預示方法,該方法是對工程辦法的改進。以單臺發動機的實測壓強數據為基礎,在發動機穩定工作段選取部分時間段進行計算,基于所建立辨識修正方法獲得燃速模型參數和燃面-肉厚數據等內彈道預示所需全部關鍵參數,進行內彈道性能預示,該方法能提高內彈道預示精度,具體流程如圖1所示。

圖1 多參數辨識及精準預示方法流程Fig.1 Multi-parameter identification and prediction method process

從發動機實測壓強數據出發,基于GA遺傳算法對燃速系數a、壓強指數n、推進劑密度ρp及噴管喉襯燒蝕速率系數as進行第一次參數辨識。通過參數辨識結果和燃面-壓強公式對燃面-肉厚曲線進行反算,進而更新多參數辨識的計算模型,再通過所得燃面數據和遺傳算法開展第二次參數辨識。根據二次辨識結果對初次辨識參數進行修正,在第一章介紹的內彈道預示模型的基礎上,根據所得參數及燃面數據對發動機內彈道性能開展正向預示。

2.3 基于遺傳算法的內彈道多參數辨識

內彈道多參數辨識主要根據發動機工作時間內的實測壓強數據進行逆向反算,該工作基于GA遺傳算法完成,對燃速系數a、壓強指數n、推進劑密度ρp及噴管喉襯燒蝕速率系數as四個參數進行辨識,具體步驟如下:

(1)根據實測數據信息及設計經驗,選定時間間隔Δt,作為計算兩條曲線擬合程度的最小時間微元;

(2)給定參數初值范圍,其中燃速系數a和壓強指數n還需滿足燃速模型關系,故需施加額外限制條件;

(3)在給定初值范圍內,隨機給出N組數據(a,n,ρp,as)作為模型初始群體,每組數據為初始群體的個體;

(4)根據零維內彈道理論建立正向預示數學模型,基于任意一組數據可進行一次內彈道性能預示,即N組數據可獲得N條壓強曲線;

(5)由于發動機穩定工作段數據更具代表性,故在穩定工作段內以Δt為間隔總計選取M個節點,建立內彈道預示模型的適應度函數,其為M個節點處預示壓強與實測壓強的差值組成的行向量,該行向量的范數作為衡量該組數據適應度高低的判據;

(6)淘汰適應度低的初始個體,保留適應度高的初始個體,在保留個體間進行數據的交換和變異,形成初始群體的子代,再次通過適應度函數對新群體進行選擇;

(7)重復上步操作直至結果收斂,將適應度最高的數據代入內彈道預示模型中,即可獲得與實測數據吻合度最高的壓強曲線。

上述為初次參數辨識的步驟。基于初次參數辨識結果進行燃面反算,根據反算燃面數據更新內彈道預示數學模型,再進行二次參數辨識修正。第二次辨識的參數初值范圍可適當縮小,以加快收斂速度并提高效率。

3 案例計算分析

為驗證內彈道性能預示方法的準確性和有效性,本文選取C1xb型發動機和一種長通道發動機作為研究案例。

3.1 C1xb型發動機案例

C1xb型發動機的裝藥組成為79% AP和21% HTPB的復合推進劑,該發動機結構參數和試車數據較為詳盡[15],許多關于發動機燃燒穩定性和內流場特性的研究均選取該發動機作為驗證案例。該型發動機的具體參數如表1所示,結構如圖2所示。

表1 C1xb發動機參數Table 1 Parameters of C1xb SRM

圖2 C1xb發動機結構圖Fig.2 Schematic of C1xb SRM

3.1.1 第一次參數辨識

(1)取最小時間微元Δt=0.01 s進行曲線擬合度對比計算;

(2)在穩定工作段內,以Δt為時間間隔選取300個點作為曲線擬合度計算時間節點;

(3)計算所有節點上預測壓強值與實測壓強值的差向量,將該行向量的Euclid范數作為適應度函數。

C1xb發動機實測壓強曲線存在初始尖峰,根據工程經驗可知其為點火藥與裝藥共同燃燒導致產氣量劇增所致[17-18]。點火過程影響因素較為復雜,但其對內彈道曲線整體趨勢影響并不顯著。辨識計算中取發動機穩定工作部分時間段作為計算區間,避開初始點火過程和壓強下降拖尾段,以提高多參數辨識的精度。

在第一次多參數辨識的過程中,設定4項待辨識參數的初始取值范圍分別為3.5 mm/(s·MPan)≤a≤5 mm/(s·MPan),0.4≤n≤0.5,1550 kg/m3≤ρp≤1650 kg/m3,0.1 mm/(s·MPa)≤as≤0.2 mm/(s·MPa)。基于遺傳算法進行多參數辨識,在經過多次迭代計算后結果收斂,求解過程如圖3所示。初次參數辨識結果為:燃速系數a=4.48 mm/(s·MPan),壓強指數n=0.487 6,推進劑密度ρp=1556 kg/m3,噴管喉襯燒蝕速率系數as=0.138 7 mm/(s·MPa)。

圖3 第一次遺傳算法優化求解過程Fig.3 Calculation progress of first multi- parameter identification

由于原始文獻未給出發動機點火裝置參數,且點火過程對宏觀內彈道特性影響較小,故將點火藥假定為球形顆粒黑火藥,再基于遺傳算法和考慮點火過程的內彈道預示模型對點火藥質量和顆粒粒徑進行辨識。該過程所需其余參數采用初次參數辨識的結果,最終獲得點火藥質量為2.5 g,粒徑為1.5 mm。根據參數辨識結果進行內彈道預示,可獲得初步預示的壓強-時間曲線,如圖4所示。

圖4 第一次辨識計算結果與實測數據對比Fig.4 Comparison between first prediction and measured pressure data

3.1.2 燃面反算及第二次參數辨識

基于初次參數辨識結果,根據燃面-壓強公式進行燃面反算。為避免點火過程對反算燃面數據產生影響,在計算時不考慮點火藥燃燒的影響,可獲得如圖5所示的反算燃面-肉厚曲線。

圖5 反算燃面-肉厚曲線Fig.5 Inverse calculation of burning area-web curve

在初次辨識基礎上開展二次辨識,基于反算燃面數據對已有數學模型進行更新。為加速收斂,第二次辨識的參數初值范圍相較第一次辨識適當縮小。設定4項待辨識參數的初始取值范圍分別為4.2 mm/(s·MPan)≤a≤4.8 mm/(s·MPan),0.45≤n≤0.5,1550 kg/m3≤ρp≤1600 kg/m3,0.12 mm/(s·MPa)≤as≤0.2 mm/(s·MPa)。如圖6所示,經過多次迭代計算后結果收斂,第二次參數辨識結果為:燃速系數a=4.648 mm/(s·MPan),壓強指數n=0.455 1,推進劑密度ρp=1589 kg/m3,噴管喉襯燒蝕率系數as=0.142 3 mm/(s·MPa)。根據反算燃面數據和第二次參數辨識結果再次進行內彈道預示,可獲得圖7所示的修正壓強-時間曲線。由圖7可知,第二次性能預示精度較第一次更高,且預示結果與實測數據高度吻合。

圖6 第二次遺傳算法優化求解過程Fig.6 Calculation progress of second multi- parameter identification

圖7 第二次辨識計算結果與實測數據對比Fig.7 Comparison between second prediction and measured pressure data

3.2 某型長通道發動機案例

測試案例二為某型長通道固體發動機,該發動機裝藥為端面包覆的內外燃管形裝藥,采用雙基推進劑并以黑火藥為點火藥,其相關參數較為詳盡[16],具體結構參數如表2所示,結構圖如圖8所示。

表2 長通道發動機參數Table 2 Parameters of the long passage SRM

圖8 長通道發動機結構圖Fig.8 Schematic of the long passage SRM

應用上述介紹的內彈道多參數辨識方法對該長通道發動機進行第一次參數辨識,第一次參數辨識結果為:燃速系數a=3.01 mm/(s·MPan),壓強指數n=0.444 1,推進劑密度ρp=1563 kg/m3,噴管喉襯燒蝕率系數as=0.112 2 mm/(s·MPa)。用第一次多參數辨識結果進行燃面反算,根據內彈道預示模型進行性能預示,計算得到如圖9所示的燃燒室壓強預示曲線和如圖10所示反算燃面曲線。由于燃面數據計算采用式(8)所示的工程方法,故上升段和下降段與實際燃面變化的物理過程存在一定差異。

圖9 長通道發動機第一次辨識預示曲線Fig.9 First prediction of the long passage SRM

圖10 長通道發動機反算燃面-肉厚曲線Fig.10 Inverse calculation of burning area-web curve of the long passage SRM

根據反算燃面數據更新模型后,進行第二次多參數辨識,得到最終參數辨識結果為:燃速系數a=3.09 mm/(s·MPan),壓強指數n=0.440 5,推進劑密度ρp=1564 kg/m3,噴管喉襯燒蝕率系數as=0.113 1 mm/(s·MPa)。由于辨識計算過程的主要依據為發動機穩定工作段壓強數據,故所得辨識參數結果具有較好的合理性。把反算的燃面數據以及修正后的參數代入內彈道正向預示模型進行計算,可獲得如圖11所示的燃燒室壓強曲線。由圖11可知,第二次性能預示精度較第一次更高,且預示結果與實測數據高度吻合。

圖11 長通道發動機第二次辨識預示曲線Fig.11 Second prediction of the long passage SRM

4 結論

(1)提出了固體火箭發動機內彈道多參數辨識及精準預示方法,并以C1xb型發動機和某型長通道發動機為研究案例進行計算分析,基于二次參數辨識修正和燃面反算,獲得了與實測壓強數據高度吻合的預示結果,驗證了該方法對提高內彈道性能預示精度的有效性。

(2)除本文案例參數外,還可針對發動機使用環境剖面下具有較強散布特性的多項預示參數開展辨識修正,相關方法也適用于具有復雜裝藥構型和燃面退移規律的固體發動機,其應用具有較強的普遍性。

(3)基于該內彈道參數辨識與預示方法,可充分利用有限次發動機試車實測數據開展內彈道性能精準預示,掌握復雜環境剖面下的內彈道性能散布情況,提高發動機設計效率并降低成本。