劉 濤， 劉錦凡， 唐國(guó)安

(1.復(fù)旦大學(xué) 航空航天系，上海 200433； 2.上海宇航系統(tǒng)工程研究所，上海 201109)

POGO是液體火箭飛行過(guò)程中一種常見(jiàn)的不穩(wěn)定自激振動(dòng)，是由于液體火箭推進(jìn)系統(tǒng)與結(jié)構(gòu)系統(tǒng)之間的耦合作用而產(chǎn)生，國(guó)內(nèi)外眾多液體火箭在飛行過(guò)程中都經(jīng)歷了POGO振動(dòng)[1-2]。因此，POGO的抑制是液體火箭設(shè)計(jì)中的一項(xiàng)重要課題。

國(guó)內(nèi)外學(xué)者對(duì)POGO振動(dòng)開(kāi)展了大量的研究工作，Rubin等[3-4]首次提出了POGO振動(dòng)的傳遞矩陣模型，采用臨界阻尼法對(duì)POGO穩(wěn)定性進(jìn)行了研究。由于傳遞矩陣模型中含有超越方程和高階多項(xiàng)式，直接求解特征值十分困難，臨界阻尼法近似認(rèn)為POGO振動(dòng)的頻率為結(jié)構(gòu)系統(tǒng)的振動(dòng)頻率，這種假設(shè)很大程度上降低了POGO穩(wěn)定性分析的計(jì)算量，因而得到廣泛應(yīng)用。譚述君等[5]對(duì)臨界阻尼法的適用性進(jìn)行了研究，指出當(dāng)結(jié)構(gòu)阻尼比小于推進(jìn)系統(tǒng)阻尼比時(shí)，臨界阻尼法是適用的，反之則不適用。Oppenheim等[6]利用類(lèi)似有限元的方法對(duì)推進(jìn)系統(tǒng)的直管、蓄壓器、泵等典型物理部件進(jìn)行了動(dòng)力學(xué)建模，進(jìn)而建立整個(gè)推進(jìn)系統(tǒng)二階振動(dòng)微分方程，實(shí)現(xiàn)了POGO振動(dòng)的時(shí)域分析，該方法亦被稱之為狀態(tài)方程法。Wang等[7]針對(duì)狀態(tài)方程法中矩陣不滿秩給時(shí)域求解帶來(lái)不便的問(wèn)題，對(duì)該方法進(jìn)行了改進(jìn)，建立了全微分形式的POGO振動(dòng)模型，提高了計(jì)算的效率。牛澤雄等[8]指出傳遞矩陣模型比狀態(tài)方程模型更加精確，狀態(tài)方程法因?qū)苈穭?dòng)力學(xué)的簡(jiǎn)化描述，影響了POGO分析的精度。郝雨等[9-10]基于推進(jìn)系統(tǒng)的傳遞矩陣模型，通過(guò)引入中間變量將其轉(zhuǎn)換為時(shí)域動(dòng)力學(xué)模型，實(shí)現(xiàn)了POGO振動(dòng)特征值的直接求解，但該方法為保證計(jì)算精度，需引入大量的中間變量，增加了計(jì)算規(guī)模和難度，此外該方法在處理時(shí)滯問(wèn)題也存在一定的不足。

針對(duì)傳遞矩陣模型精度更高，但由于傳遞函數(shù)中存在超越函數(shù)，導(dǎo)致特征值密集且難以直接求解這一問(wèn)題，本文以某型號(hào)液體火箭為研究對(duì)象，首先建立火箭POGO振動(dòng)閉環(huán)傳遞函數(shù)，基于矢量擬合法獲取閉環(huán)傳遞函數(shù)的極點(diǎn)。進(jìn)而對(duì)不同蓄壓器PV值下的POGO穩(wěn)定性進(jìn)行分析，確定了蓄壓器的設(shè)計(jì)狀態(tài)。通過(guò)與臨界阻尼法的對(duì)比，證明該方法具有更高的計(jì)算精度。該方法對(duì)其他液體火箭的POGO穩(wěn)定性分析同樣適用，具有一定的工程意義。

1 POGO振動(dòng)傳遞函數(shù)

(1)

圖1 POGO閉環(huán)系統(tǒng)原理圖

1.1 結(jié)構(gòu)系統(tǒng)傳遞函數(shù)

在廣義坐標(biāo)系下箭體結(jié)構(gòu)的振動(dòng)方程為

(2)

(3)

由此可得到箭體結(jié)構(gòu)系統(tǒng)的傳遞函數(shù)

(4)

1.2 推進(jìn)系統(tǒng)傳遞函數(shù)

液體火箭推進(jìn)系統(tǒng)包含氧化劑系統(tǒng)和燃料系統(tǒng)，主要由輸送管(直管、彎管和波紋管等)、泵、阻力元件、蓄壓器和燃燒室等部分組成，如圖2所示。采用模塊化方法建立各部分的傳遞矩陣模型，該方法可以直觀并且精確的描述管路元件輸入輸出端脈動(dòng)壓力與脈動(dòng)流量之間的傳遞關(guān)系。對(duì)于直管、泵和阻力元件等部件的傳遞矩陣可參考文獻(xiàn)[3]。

圖2 某火箭輸送系統(tǒng)結(jié)構(gòu)示意圖[11]

蓄壓器安裝在推進(jìn)劑輸送管路上，是POGO抑制最有效的措施之一，蓄壓器傳遞矩陣為

(5)

(6)

(7)

式中:P1,P2,Q1和Q2分別為蓄壓器入口和出口的脈動(dòng)壓力和脈動(dòng)流量;P,V分別為蓄壓器初始?jí)毫统跏既莘e;Px為發(fā)動(dòng)機(jī)啟動(dòng)活門(mén)入口壓力;K為絕熱指數(shù);Ra,La分別為蓄壓器的等效電阻和等效電感。

文獻(xiàn)[3]中將燃燒室方程的時(shí)滯項(xiàng)進(jìn)行了線性化近似，一定程度上簡(jiǎn)化了計(jì)算，但當(dāng)頻率較高時(shí)，傳遞函數(shù)的精度會(huì)受到明顯影響，本文采用的發(fā)動(dòng)機(jī)燃燒室方程為

(8)

(2) 1981年年平均流量和年最大流量分別作為邊界條件時(shí)，巴塘河河道縱剖面總流速與流速水頭峰值出現(xiàn)位置基本一致。年最大流量下，河道縱剖面最末端最大總流速可達(dá)14 m/s，其他位置為0～5 m/s；年最大流量下，河道縱剖面最末端流速水頭達(dá)到11 m，其他位置為0～1.5 m。年平均流量下，河道縱剖面最末端最大總流速可達(dá)9 m/s，其他位置為0～4 m/s；年平均流量下，河道縱剖面最末端流速水頭達(dá)到4 m，其他位置為0～1 m。

作用在運(yùn)載火箭上典型外力主要有：發(fā)動(dòng)機(jī)燃燒室脈動(dòng)壓力引起的脈動(dòng)推力，貯箱底開(kāi)口處的脈動(dòng)壓力對(duì)結(jié)構(gòu)的作用力，貯箱底開(kāi)口處相對(duì)脈動(dòng)流量引起的液體質(zhì)心脈動(dòng)對(duì)結(jié)構(gòu)的作用力，泵前脈動(dòng)壓力和動(dòng)量變化對(duì)結(jié)構(gòu)的作用力。通過(guò)模態(tài)坐標(biāo)轉(zhuǎn)換，廣義力TN可表示為

TN=NAthCφePc+AtoφtoPto+NAtfφtfPtf-NAo3φpoPso-

(9)

式中:N為發(fā)動(dòng)機(jī)數(shù)量;C為發(fā)動(dòng)機(jī)推力系數(shù);ho和hf分別為氧化劑和燃料的液位高度；Ato和Atf分別為與氧化劑箱和燃料箱底部相連管路的截面積；Apo和Apf分別為泵前管路的截面積；Pso,Psf,Qpo和Qpf分別為氧化劑泵和燃料泵的泵前脈動(dòng)壓力和泵后脈動(dòng)流量；Qto和Qtf分別為氧化劑箱和燃料箱箱底脈動(dòng)流量；φto,φtf,φpo和φpf分別為氧化劑和燃料的箱底、泵前管路對(duì)應(yīng)的振型。

根據(jù)推進(jìn)系統(tǒng)各組件傳遞函數(shù)矩陣之間的關(guān)系，可得

(10)

式中:A和B為系數(shù)矩陣;X為各部件的脈動(dòng)壓力、流量等狀態(tài)變量。由式可計(jì)算推進(jìn)系統(tǒng)反饋傳遞函數(shù)H(s)。由于H(s)中存在時(shí)滯項(xiàng)和超越函數(shù)和高階多項(xiàng)式，給極點(diǎn)的直接求解帶來(lái)了極大的困難。

2 矢量擬合法

矢量擬合法是由B.Gustavsen提出的一種高效的有理分式擬合方法[12-14]，能夠?qū)㈩l域離散數(shù)據(jù)按極點(diǎn)-留數(shù)的形式進(jìn)行擬合，由于避免了在數(shù)值擬合過(guò)程中引起的不平衡加權(quán)及病態(tài)缺陷等問(wèn)題，被廣泛應(yīng)用于電力、控制等領(lǐng)域的頻率特性擬合。

對(duì)于某一函數(shù)f(s)，其有理分式逼近可表示為

(11)

(12)

將式中σ(s)乘以f(s)可得

(13)

式(13)為超定線性矩陣方程，將復(fù)頻域離散采樣值[s,f(s)]代入其中，通過(guò)最小二乘法可確定an,cn,d和h等未知量。

進(jìn)一步，可將式有理分式采用零極點(diǎn)形式表示

(14)

(15)

通過(guò)式(14)～式(15)可得到

(16)

式(16)將求解f(s)的極點(diǎn)問(wèn)題轉(zhuǎn)化為求解σ(s)的零點(diǎn)問(wèn)題，即從非線性問(wèn)題轉(zhuǎn)變?yōu)榫€性求解問(wèn)題。

將求解得到新的極點(diǎn)作為初始極點(diǎn)再進(jìn)行迭代計(jì)算，經(jīng)過(guò)多個(gè)迭代步計(jì)算后可得到較為精確的極點(diǎn)。

3 POGO穩(wěn)定性分析

將上述建模方法應(yīng)用于某液體火箭的POGO穩(wěn)定性分析中，選取火箭一級(jí)飛行0, 30 s, 60 s, 90 s, 120 s和150 s等6個(gè)典型秒態(tài)，研究各典型秒態(tài)下推進(jìn)系統(tǒng)與箭體一階縱向振動(dòng)下的POGO穩(wěn)定性。

圖3為火箭一級(jí)飛行150 s時(shí)未安裝蓄壓器狀態(tài)下的推進(jìn)系統(tǒng)和推進(jìn)-結(jié)構(gòu)閉環(huán)系統(tǒng)傳遞函數(shù)曲線。從圖中可以看出，推進(jìn)系統(tǒng)傳遞函數(shù)在0～100 Hz有多個(gè)共振頻率點(diǎn)，而推進(jìn)-結(jié)構(gòu)閉環(huán)傳遞函數(shù)有一個(gè)較為明顯主共振點(diǎn)。

基于矢量擬合法對(duì)圖3中的傳遞函數(shù)曲線進(jìn)行有理分式擬合，由于有理分式擬合精度與式(18)的N值有關(guān)，通常N值越大，擬合精度越高，但隨著N值增加，將產(chǎn)生大量不穩(wěn)定的數(shù)學(xué)極點(diǎn)，影響傳遞函數(shù)物理極點(diǎn)的識(shí)別，本文引入模態(tài)辨識(shí)中的穩(wěn)態(tài)圖法[15]來(lái)確定傳遞函數(shù)中穩(wěn)定的極點(diǎn)。穩(wěn)態(tài)圖法基本思路是：假定系統(tǒng)的階次從Nmin到Nmax，把各階次下的特征值計(jì)算結(jié)果畫(huà)在橫軸為頻率縱軸為階次的坐標(biāo)圖中，可得到穩(wěn)定圖。若每次增加階次后，得到的極點(diǎn)和留數(shù)在容差范圍內(nèi)，則認(rèn)為是相同的極點(diǎn)，且隨著階次的增加在穩(wěn)態(tài)圖中形成多個(gè)穩(wěn)定軸，穩(wěn)定軸所對(duì)應(yīng)的極點(diǎn)便是系統(tǒng)的特征值。

(a) 推進(jìn)系統(tǒng)

(b) 推進(jìn)-結(jié)構(gòu)閉環(huán)系統(tǒng)

圖4為有理分式擬合結(jié)果，基于穩(wěn)態(tài)圖可以確定在0～100 Hz內(nèi)推進(jìn)系統(tǒng)有5對(duì)共軛極點(diǎn)，如表1所示(表中只列出正虛部的極點(diǎn)，表2同)。推進(jìn)-結(jié)構(gòu)閉環(huán)系統(tǒng)有7對(duì)共軛極點(diǎn)，前5階極點(diǎn)如表2所示。表1中推進(jìn)系統(tǒng)的一階頻率為12.57 Hz，與箭體縱向一階共振頻率12.8 Hz較為接近。由表2可知，推進(jìn)-結(jié)構(gòu)閉環(huán)系統(tǒng)一階共振頻率為12.44 Hz，對(duì)應(yīng)極點(diǎn)實(shí)部為0.97，說(shuō)明閉環(huán)系統(tǒng)是不穩(wěn)定的，即出現(xiàn)POGO振動(dòng)。表2中其他極點(diǎn)的實(shí)部均小于0，說(shuō)明當(dāng)箭體結(jié)構(gòu)頻率與推進(jìn)系統(tǒng)頻率較為接近時(shí)，容易出現(xiàn)POGO振動(dòng)，而且POGO振動(dòng)最有可能發(fā)生在箭體結(jié)構(gòu)頻率點(diǎn)附近。

(a) 推進(jìn)系統(tǒng)

(b) 推進(jìn)-結(jié)構(gòu)閉環(huán)系統(tǒng)

表1 推進(jìn)系統(tǒng)極點(diǎn)

表2 推進(jìn)-結(jié)構(gòu)閉環(huán)系統(tǒng)極點(diǎn)

設(shè)計(jì)上，通過(guò)改變蓄壓器膜盒充氣壓力P和膜盒容積V來(lái)達(dá)到POGO抑制的效果。本文對(duì)蓄壓器PV值在0(未安裝蓄壓器)，0.003 92 MPaL, 0.039 0 MPaL, 0.078 4 MPaL, …,0.588 MPaL等17個(gè)狀態(tài)下的POGO穩(wěn)定性分別采用矢量擬合法和臨界阻尼法進(jìn)行分析。

圖5為各秒態(tài)下的臨界阻尼值。圖6為閉環(huán)傳遞函數(shù)在結(jié)構(gòu)縱向振動(dòng)頻率附近的極點(diǎn)分布。兩種方法都說(shuō)明通過(guò)改變蓄壓器PV值可以達(dá)到POGO抑制的效果。而且兩種方法對(duì)POGO系統(tǒng)的穩(wěn)定性的判斷是一致的。但矢量擬合法可有效解決臨界阻尼法的適用性問(wèn)題，同時(shí)比臨界阻尼法具有更高的精度。究其原因是臨界阻尼法求解時(shí)近似認(rèn)為出現(xiàn)POGO振動(dòng)的點(diǎn)為箭體結(jié)構(gòu)的共振頻率點(diǎn)，而實(shí)際上推進(jìn)-結(jié)構(gòu)耦合后POGO振動(dòng)的頻率點(diǎn)有一定的改變，因此矢量擬合法的分析精度更高。

圖5 各秒態(tài)下POGO系統(tǒng)臨界阻尼值

圖6表明并不是蓄壓器PV值越高越好，除0外，其他秒態(tài)閉環(huán)極點(diǎn)實(shí)部隨蓄壓器的PV值均呈現(xiàn)先減小，后增大趨勢(shì)。從圖6可知，當(dāng)蓄壓器PV值在0.157～0.196 MPaL左右時(shí)，閉環(huán)極點(diǎn)的實(shí)部最小，POGO抑制效果最好。

(a)

(b)

(c)

(d)

(e)

(f)

4 結(jié) 論

本文采用矢量擬合法實(shí)現(xiàn)對(duì)火箭推進(jìn)系統(tǒng)-結(jié)構(gòu)系統(tǒng)傳遞函數(shù)的高精度擬合，借助穩(wěn)態(tài)圖法可確定傳遞函數(shù)穩(wěn)定的極點(diǎn)，進(jìn)而通過(guò)極點(diǎn)實(shí)部值的正負(fù)判斷POGO穩(wěn)定性。通過(guò)與臨界阻尼法的分析結(jié)果對(duì)比，表明矢量擬合法和臨界阻尼法結(jié)論一致，但矢量擬合法具有更高的分析精度。

當(dāng)推進(jìn)系統(tǒng)頻率與箭體結(jié)構(gòu)頻率較為接近時(shí)，易出現(xiàn)POGO振動(dòng)，通過(guò)改變蓄壓器PV值可實(shí)現(xiàn)對(duì)POGO振動(dòng)的抑制，針對(duì)本文所研究的液體火箭，當(dāng)蓄壓器PV值在0.157～0.196 MPaL，POGO振動(dòng)抑制的效果最好。