高床載單組元肼發動機冷起動過程仿真研究

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

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

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

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

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

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

1 數學模型

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

1.1 管路模型

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

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

(1)

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

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

(2)

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

1.2 電磁閥模型

對于系統仿真中所關注的電磁閥響應時間和閥門開度建立相應的數學模型為[14]

(3)

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

當閥門打開時，有

(4)

當閥門關閉時，有

(5)

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

1.3 催化床模型

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

1.3.1 冷起動延遲模型

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

(6)

式(6)中：M為氣體產物的平均摩爾質量；a0為催化劑顆粒直徑；Pm=2σcosθ/R為毛細壓力，其中，σ為液態推進劑表面張力，R為大孔的特征半徑，θ為接觸角；Ru為氣體常數；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)為推進劑的分解反應速率常數。

1.3.2 催化床建壓過程模型

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

(7)

(8)

催化床內壓力由氣體狀態方程可得

(9)

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

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

2 仿真模型的建立及驗證

以典型組件的數學模型為基礎，根據模塊化的建模思想，在MWorks軟件平臺上編制開發各典型組件的仿真模塊。根據某次單臺發動機地面試驗的試車狀態，對圖1所示單臺發動機地面試驗系統搭建了相應的仿真模型。

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

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

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

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

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

3 結果及分析

3.1 貯箱壓力

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

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

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

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

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

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

3.2 系統流阻

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

圖5 不同結構參數對發動機冷起動特性的影響Fig.5 Influence of variable structure parameters on characteristic parameters of cold starting

圖5給出了各個因素引起系統流阻增大的情況下，發動機冷起動特性參數的變化，如響應時間、壓力峰和穩態室壓。可以看出，毛細管長度的增加顯著增大了發動機的起動響應時間并減小了發動機的穩態室壓；毛細管通徑、電磁閥壓降以及限流圈喉徑的變化對響應時間的影響較小；毛細管通徑和限流圈喉徑對起動壓力峰值的抑制作用相當，但是毛細管通徑所引起流阻變化對發動機穩態室壓減小的作用更強，不利于系統穩態參數調節；而電磁閥壓降所對穩態室壓減小作用較小，相應的抑制壓力峰作用也相對較差。因此，增大毛細管長度、電磁閥壓降，減小毛細管直徑、限流圈直徑都可以在不同程度上抑制發動機的冷起動壓力峰，但同時也影響了發動機的起動加速性和穩態工作性能。因此需要進一步對上述因素進行對比，綜合考量3種結構參數對壓力峰的影響。

結合工程應用中需要盡量保證穩態室壓不變化，比較在保證穩態室壓不發生改變的情況下4種增大流阻方式抑制壓力峰的能力。考慮到3.1節貯箱壓力和穩態室壓成線性關系，因此通過提高貯箱壓力以保證穩態室壓恒定的情況下對上述4種結構參數變化進行仿真計算，對比不同結構參數引起相同流阻情況下冷起動壓力峰變化的結果如圖6所示。

圖6 相同流阻變化下結構參數對冷起動壓力峰的影響Fig.6 Influence of variable structure parameters on pressure spike of cold starting

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

4 結論

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

(2)貯箱壓力對發動機的冷起動響應時間影響較小，貯箱壓力升高時，冷起動壓力峰和穩態室壓都隨之增加，冷起動壓力峰增大幅值約為穩態室壓增大幅度的2倍。

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

(4)對比了3種結構參數變化對冷起動過程的影響規律，發現從減小冷起動過程瞬時流量角度出發，綜合控制限流圈內徑是最為合理有效的抑制冷起動壓力峰的措施。