冀 鵬, 梁樹強， 肖明杰， 雷凡培

(1.中國航天科技集團有限公司西安航天動力研究所液體火箭發(fā)動機重點實驗室， 西安 710100； 2.航天推進技術(shù)研究院， 西安 710100； 3.中國船舶集團有限公司， 北京 100044)

床載荷是單組元發(fā)動機的重要設(shè)計參數(shù)[1]。床載荷越高，發(fā)動機推力室的結(jié)構(gòu)尺寸越小，結(jié)構(gòu)質(zhì)量越輕。尤其在飛行器用于變軌的大推力單組元發(fā)動機提高床載時，該輕質(zhì)化作用對整個飛行器就更為明顯。然而不同于其他單組元發(fā)動機，單組元肼催化分解發(fā)動機可以實現(xiàn)室溫下的自發(fā)起動，即冷起動，大大拓寬了其任務(wù)適用性，但在冷起動過程中常伴隨著明顯的壓力峰，極大地影響著發(fā)動機的工作可靠性和飛行器的變軌穩(wěn)定性[2]。因此，開展高床載單組元肼發(fā)動機的冷起動過程研究顯得尤為重要。

單組元肼催化分解發(fā)動機冷起動過程的現(xiàn)象描述以及冷起動對推力器性能壽命的影響關(guān)系主要是由試驗研究完成的。Kagawa等[3]采用中子放射線照相技術(shù)實現(xiàn)了催化床內(nèi)液態(tài)肼流動狀態(tài)的可視化，試驗表明噴注壓降并不影響肼的分解區(qū)域，而冷起動影響了在發(fā)動機催化床內(nèi)肼的流動區(qū)域。劉昌國等[4]進行了單組元300 N發(fā)動機低溫試驗，發(fā)現(xiàn)溫度差異、催化劑活性、低溫下電磁閥響應(yīng)對該發(fā)動機工作特性的影響。由于缺乏充分的過程認識和可靠的數(shù)值仿真手段，阿里安5號運載火箭[5]、美國航空航天局(National Aeronautics and Space Administration，NASA)空間飛行器[6-7]等所采用的單組元發(fā)動機也主要采用試驗的手段來鑒定發(fā)動機的冷起動能力并評判其對發(fā)動機性能的影響。但依靠試驗的方法難以對影響發(fā)動機起動過程的諸多因素及其組合開展全面系統(tǒng)的分析。

而在單組元催化分解發(fā)動機起動過程數(shù)值仿真方面，Schmitz等[8]將肼推力器的起動過程假設(shè)為氣體充填和催化床溫升換熱兩個子過程，并分別建立相應(yīng)數(shù)學(xué)模型得到了發(fā)動機的起動過程數(shù)學(xué)模型，但該數(shù)學(xué)模型中采用修正因子來描述催化床入口處推進劑的分解，只能通過試驗來獲取修正因子限制了該模型的適用性。孫得川等[9]在此基礎(chǔ)上進一步考慮相變過程通過計算流體力學(xué)的方法分析了硝酸羥胺基單組元發(fā)動機的起動過程，發(fā)現(xiàn)考慮相變過程后計算值與試驗值吻合更好，但零維模型的計算結(jié)果仍令人滿意。

上述學(xué)者通過試驗和數(shù)值仿真的方法都將重點放在了起動邊界溫度和推力室內(nèi)部復(fù)雜物理化學(xué)過程的研究上，而液體火箭發(fā)動機的起動過程是一個能量高密度釋放和工況復(fù)雜變化的耦合過程[10]，起動控制對發(fā)動機起動過程也極具影響[11]。但目前尚未見從發(fā)動機系統(tǒng)的角度對單組元冷起動過程開展研究。為此，從減小起動瞬時流量的角度出發(fā)，建立反映冷起動延遲與分解過程平均效應(yīng)的發(fā)動機冷起動過程數(shù)學(xué)模型，對某型高床載的單組元肼發(fā)動機冷起動過程開展仿真計算，分析系統(tǒng)結(jié)構(gòu)對冷起動過程的影響規(guī)律，旨在為單組元肼催化分解發(fā)動機抑制冷起動壓力峰提供理論參考與支持。

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

單組元姿控動力系統(tǒng)是由氣路增壓系統(tǒng)、液路推進劑供應(yīng)系統(tǒng)以及推力裝置組成。對于恒壓擠壓的單組元姿控動力系統(tǒng)，貯箱增壓壓力基本保持不變，貯箱內(nèi)液注靜壓和慣性壓力的影響相對較小，忽略流體流出貯箱進入管路的損失[12]，假設(shè)貯箱內(nèi)氣體壓力恒定，因此建立液路推進劑供應(yīng)系統(tǒng)相關(guān)組件以及推力裝置的動力學(xué)模型。

1.1 管路模型

在建立液體管路模型時不考慮管路的變形，認為管道為剛性的；管路內(nèi)液體是一維的，管道橫截面上的速度分布是均勻的。采用分段集中參數(shù)的管路模型。

對于特定的一段長度L、直徑d、截面積A、體積V和粗糙度h的管路，考慮慣性、黏性時管路方程為

(1)

式(1)中：j為管路的慣性系數(shù)；qm為管路內(nèi)流量；p1、p2分別為管路的入口壓力與出口壓力；ρp為管內(nèi)液體密度；t為時間；ξ為總的流阻系數(shù)，由管路的沿程損失和局部損失組成[13]。

同時考慮液體的壓縮性方程：

(2)

式(2)中：p為分段管路等效容積的壓力；qm1、qm2為管路等效容積出入口流量；z=V/a2為管路的流容系數(shù)，表征管路流體的壓縮性，其中a為液體推進劑的聲速。

1.2 電磁閥模型

對于系統(tǒng)仿真中所關(guān)注的電磁閥響應(yīng)時間和閥門開度建立相應(yīng)的數(shù)學(xué)模型為[14]

(3)

式(3)中:qval為通過閥門流量；Aval為閥門流通截面積；pval、pv分別為閥門入口壓力與出口壓力；Γ為電磁閥開度與流量系數(shù)的乘積(即Γ=τCv，其中Cv為閥門的流量系數(shù))，它主要由實驗結(jié)果來確定。但當(dāng)缺乏實驗數(shù)據(jù)時，通常采用如下方法來確定。

當(dāng)閥門打開時，有

(4)

當(dāng)閥門關(guān)閉時，有

(5)

式中：tc,k和tc,g分別為電磁閥開啟與關(guān)閉的響應(yīng)時間；m為試驗曲線擬合指數(shù)，主要由Γ-t曲線形狀而定，一般情況下m為1～3。

1.3 催化床模型

采用反映平均效應(yīng)的零維催化床模型。在建立催化床模型時做如下假設(shè)：分解產(chǎn)物為理想氣體；不考慮氣體的不均勻性和波動過程。

1.3.1 冷起動延遲模型

在單組元發(fā)動機冷起動時，由于初始溫度較低、催化劑上微孔直徑較小(微米量級)，液態(tài)推進劑會滲透到催化劑顆粒微孔內(nèi)部進行分解。催化劑微孔內(nèi)部的氣液分界面會隨著分解逐步外移直至分解氣體溢出微孔，在整個催化床內(nèi)建壓。基于上述冷起動過程分析，對催化劑顆粒內(nèi)部進行建模，計算氣液分界面在顆粒微孔內(nèi)部壓力的作用下逐漸推移到微孔外的動態(tài)過程，得到單組元發(fā)動機冷起動延遲時間為[15]

(6)

式(6)中：M為氣體產(chǎn)物的平均摩爾質(zhì)量；a0為催化劑顆粒直徑；Pm=2σcosθ/R為毛細壓力，其中，σ為液態(tài)推進劑表面張力，R為大孔的特征半徑，θ為接觸角；Ru為氣體常數(shù)；T為初始溫度；D*P*為氣相推進劑的擴散能力，D*P*=4.487×10-5T1.823；(ρA)vp為推進劑飽和蒸汽密度，(ρA)vp=2.45×108exp(-5 175/T)/T；ks=1010exp(-1 389/T)為推進劑的分解反應(yīng)速率常數(shù)。

1.3.2 催化床建壓過程模型

考慮到催化床內(nèi)復(fù)雜的物理化學(xué)過程，引入燃燒時滯模型，假設(shè)推進劑在一定轉(zhuǎn)化時間之后分解為燃氣(此時間包括了推進劑蒸發(fā)、流動、分解時間的總和)。假設(shè)在該轉(zhuǎn)化時間內(nèi)推進劑的轉(zhuǎn)化速率均勻，并且分解后的氣體在任何瞬時也是均勻分布的，即氨的解離過程也近似認為與DT-3的分解過程同時進行。則液體推進劑及氣體分解產(chǎn)物在催化床中的質(zhì)量變化分別為[16-17]

(7)

(8)

催化床內(nèi)壓力由氣體狀態(tài)方程可得

(9)

式(9)中：Pc為催化床平均壓力；Vg為催化床內(nèi)氣體體積；Rc為催化床內(nèi)氣體常數(shù)；Tc為催化床內(nèi)溫度。

由于催化床內(nèi)裝載著緊密的催化劑，因此推進劑通過催化床時需要克服阻力稱為床流阻，床流阻是催化床設(shè)計中的主要參數(shù)之一，采用文獻[17]中的試驗擬合公式計算。

2 仿真模型的建立及驗證

以典型組件的數(shù)學(xué)模型為基礎(chǔ)，根據(jù)模塊化的建模思想，在MWorks軟件平臺上編制開發(fā)各典型組件的仿真模塊。根據(jù)某次單臺發(fā)動機地面試驗的試車狀態(tài)，對圖1所示單臺發(fā)動機地面試驗系統(tǒng)搭建了相應(yīng)的仿真模型。

考慮到所建立的數(shù)學(xué)模型是一組常微分代數(shù)方程組，在數(shù)值計算方法上選用四階Range-Kutta方法進行求解模型，在初始參數(shù)賦值后，對該地面試驗系統(tǒng)的冷起動過程進行了仿真計算。

考慮到試驗系統(tǒng)中電磁流量計的測量頻率較低，在較短的起動過程中無法準確測量系統(tǒng)的瞬時流量，因此主要對發(fā)動機室壓Pc、閥前壓力Pv的試驗測量值與仿真值進行比較，驗證仿真模型的準確性，計算結(jié)果和試驗測量的對比如圖2所示。對所有狀態(tài)參數(shù)以設(shè)計值為參考都進行無量綱化，時間則以起動加速性要求為參考。

P為測量壓力；Q為測量流量圖1 發(fā)動機地面試驗系統(tǒng)示意圖Fig.1 engine ground test system schematic

圖2 發(fā)動機冷起動過程仿真結(jié)果與試驗結(jié)果對比Fig.2 Simulation results compare with ground test results of engine cold starting process

從圖2可以看出，發(fā)動機閥前壓力的計算值與實測值吻合較好，而發(fā)動機室壓的計算值與實測值之間有著一定的誤差，穩(wěn)態(tài)室壓計算值較試驗值偏高，引起這一誤差的主要原因是試驗測量壓力值為噴管入口處的壓力，而計算得到的壓力是整個催化床內(nèi)平均壓力值；其次發(fā)動機冷起動時催化床初溫較低隨著反應(yīng)進行分解效率隨時間逐漸升高，而仿真中沒有考慮分解效率的變化，將分解效率設(shè)為恒定的穩(wěn)態(tài)分解效率，這兩方面因素導(dǎo)致計算的穩(wěn)態(tài)室壓較高。計算的建壓曲線相比試驗建壓曲線更為陡峭，則是因為催化床模型采用了反映平均效應(yīng)的推進劑分解過程，假設(shè)在延遲時間內(nèi)所積存的推進劑在同一時間發(fā)生分解，而實際中推進劑是由先進入的部分催化劑分解建壓然后引起所積存推進劑發(fā)生鏈式分解，這導(dǎo)致計算建壓過程快于試驗測量的建壓過程。

3 結(jié)果及分析

3.1 貯箱壓力

對于擠壓式發(fā)動機系統(tǒng)而言，貯箱壓力直接決定著起動過程中推力裝置的入口壓力，對起動過程中進入管路及推力裝置的瞬時流量有著至關(guān)重要的作用。對貯箱壓力從設(shè)計值逐步增大進行了仿真計算，圖3給出了在不同貯箱壓力情況下發(fā)動機的冷起動過程曲線。

圖3 貯箱壓力對發(fā)動機冷起動過程的影響Fig.3 Influence of tank pressure on cold starting process of engine

貯箱壓力越高，在起動過程推力室未建壓前的推進劑流路中瞬時壓差就越大，推進劑起動瞬時流量越大，在推力裝置電磁閥打開后充填毛細管以及集液腔的速度就越快，在相同的冷起動時滯內(nèi)催化床中積存推進劑就越多，提供給爆式分解的初始能量就越大，從而導(dǎo)致發(fā)動機冷起動時壓力峰較高。

從圖3可以看出，發(fā)動機冷起動壓力峰隨貯箱壓力的升高略有增大；同時，貯箱壓力的提高也引起了發(fā)動機在穩(wěn)態(tài)工作時壓力的升高，這是由于貯箱壓力升高后系統(tǒng)壓降提高，進入催化床內(nèi)推進劑流量增大而導(dǎo)致了發(fā)動機穩(wěn)態(tài)室壓的提高。進一步得到貯箱壓力變化與壓力峰變化和穩(wěn)態(tài)室壓變化之間的關(guān)系，如圖4所示。

圖4 貯箱壓力變化對冷起動壓力峰及室壓的影響Fig.4 Influence of variable tank pressure on pressure spike and chamber pressure of cold starting

從圖4可以看出，對該型單機系統(tǒng)在所計算的貯箱壓力變化范圍內(nèi)，壓力峰變化和穩(wěn)態(tài)室壓變化與貯箱壓力變化呈線性關(guān)系，冷起動壓力峰的增大值約是貯箱壓力提高值的62%，所造成穩(wěn)態(tài)室壓的增大值約是貯箱壓力提高值的33%。穩(wěn)態(tài)室壓的增大幅度約為冷起動壓力峰增大幅度的1/2。這是由于冷起動時推力室未建壓，而穩(wěn)態(tài)工作下推力室室壓是隨貯箱壓力增大而增大的，因而在相同貯箱壓力變化下，冷起動過程系統(tǒng)壓降增大值大于穩(wěn)態(tài)工作的系統(tǒng)壓降增大值，冷起動過程流量相對變化也就大于穩(wěn)態(tài)工作的流量變化。同時，在穩(wěn)態(tài)流量較大的情況下，床流阻也會相應(yīng)增大，故穩(wěn)態(tài)室壓隨入口壓力的變化幅值較小。

3.2 系統(tǒng)流阻

單臺推力室的動力系統(tǒng)組件較少，在發(fā)動機系統(tǒng)參數(shù)發(fā)生變化時都會改變相應(yīng)組件的流阻特性進而影響發(fā)動機的冷起動過程。因此，對不同毛細管長度、毛細管直徑、電磁閥壓降以及限流圈直徑等典型組件的結(jié)構(gòu)參數(shù)進行了仿真計算。

圖5 不同結(jié)構(gòu)參數(shù)對發(fā)動機冷起動特性的影響Fig.5 Influence of variable structure parameters on characteristic parameters of cold starting

圖5給出了各個因素引起系統(tǒng)流阻增大的情況下，發(fā)動機冷起動特性參數(shù)的變化，如響應(yīng)時間、壓力峰和穩(wěn)態(tài)室壓?？梢钥闯觯毠荛L度的增加顯著增大了發(fā)動機的起動響應(yīng)時間并減小了發(fā)動機的穩(wěn)態(tài)室壓；毛細管通徑、電磁閥壓降以及限流圈喉徑的變化對響應(yīng)時間的影響較??；毛細管通徑和限流圈喉徑對起動壓力峰值的抑制作用相當(dāng)，但是毛細管通徑所引起流阻變化對發(fā)動機穩(wěn)態(tài)室壓減小的作用更強，不利于系統(tǒng)穩(wěn)態(tài)參數(shù)調(diào)節(jié)；而電磁閥壓降所對穩(wěn)態(tài)室壓減小作用較小，相應(yīng)的抑制壓力峰作用也相對較差。因此，增大毛細管長度、電磁閥壓降，減小毛細管直徑、限流圈直徑都可以在不同程度上抑制發(fā)動機的冷起動壓力峰，但同時也影響了發(fā)動機的起動加速性和穩(wěn)態(tài)工作性能。因此需要進一步對上述因素進行對比，綜合考量3種結(jié)構(gòu)參數(shù)對壓力峰的影響。

結(jié)合工程應(yīng)用中需要盡量保證穩(wěn)態(tài)室壓不變化，比較在保證穩(wěn)態(tài)室壓不發(fā)生改變的情況下4種增大流阻方式抑制壓力峰的能力?？紤]到3.1節(jié)貯箱壓力和穩(wěn)態(tài)室壓成線性關(guān)系，因此通過提高貯箱壓力以保證穩(wěn)態(tài)室壓恒定的情況下對上述4種結(jié)構(gòu)參數(shù)變化進行仿真計算，對比不同結(jié)構(gòu)參數(shù)引起相同流阻情況下冷起動壓力峰變化的結(jié)果如圖6所示。

圖6 相同流阻變化下結(jié)構(gòu)參數(shù)對冷起動壓力峰的影響Fig.6 Influence of variable structure parameters on pressure spike of cold starting

從圖6可以看出，限流圈喉徑所表示的曲線斜率最大，這意味著在引起系統(tǒng)流阻相同改變的前提下，限流圈減小壓力峰的作用最強，而毛細管長度增長的抑制作用最弱。因此，從保證起動加速性與穩(wěn)態(tài)室壓恒定的要求出發(fā)，調(diào)節(jié)系統(tǒng)結(jié)構(gòu)參數(shù)中設(shè)置合理的限流圈對抑制冷起動壓力峰的作用最強，同時綜合研制周期與系統(tǒng)復(fù)雜性采取調(diào)節(jié)限流圈參數(shù)的措施也較為合理。

4 結(jié)論

(1)采用所建立的計算模型和計算方法，可以對單組元肼催化分解發(fā)動機冷起動過程進行數(shù)值仿真，得到的冷起動室壓和閥前壓力與試驗結(jié)果吻合較好。

(2)貯箱壓力對發(fā)動機的冷起動響應(yīng)時間影響較小，貯箱壓力升高時，冷起動壓力峰和穩(wěn)態(tài)室壓都隨之增加，冷起動壓力峰增大幅值約為穩(wěn)態(tài)室壓增大幅度的2倍。

(3)增大毛細管長度會顯著降低發(fā)動機的冷起動加速性，且對冷起動壓力峰的抑制作用較小。

(4)對比了3種結(jié)構(gòu)參數(shù)變化對冷起動過程的影響規(guī)律，發(fā)現(xiàn)從減小冷起動過程瞬時流量角度出發(fā)，綜合控制限流圈內(nèi)徑是最為合理有效的抑制冷起動壓力峰的措施。