推力作用下變質(zhì)量柔性自旋飛行器穩(wěn)定性分析

潘成龍，榮吉利，徐天富，項大林

(1. 北京理工大學宇航學院，北京 100081；2. 兵器工業(yè)集團航空彈藥研究院，哈爾濱 150030；3. 北京宇航系統(tǒng)工程研究所，北京 100076)

0 引 言

大長徑比和高推重比自旋飛行器彈性效應(yīng)明顯，耦合問題突出，推力穩(wěn)定性問題備受關(guān)注。文獻[1]采用梁結(jié)構(gòu)振動理論建立柔性彈箭運動模型，分析了推力作用對系統(tǒng)振動頻率和系統(tǒng)穩(wěn)定性的影響。文獻[2-4]采用Timoshenko梁為模型，考慮陀螺效應(yīng)和剪切效應(yīng)，分析了隨動推力作用下柔性自旋飛行器穩(wěn)定性問題。文獻[5]研究了變推力作用下系統(tǒng)質(zhì)量偏心因素對自旋飛行器系統(tǒng)動力特性的影響。文獻[6]采用變質(zhì)量系統(tǒng)的方法，分析了旋轉(zhuǎn)固體火箭發(fā)動機工作過程中的章動不穩(wěn)定性問題。除此之外，飛行器在飛行時是一個典型的時變系統(tǒng)，燃料消耗引起變質(zhì)量效應(yīng)，對飛行過程有很大影響，其振動規(guī)律不同于恒定質(zhì)量系統(tǒng)[7]。文獻[8]以Kane運動方程為基礎(chǔ)，建立考慮時變質(zhì)量和幾何剛度的柔性火箭結(jié)構(gòu)動力學方程。文獻[9]運用拉格朗日方程建立時變質(zhì)量柔性火箭動力學方程，同時考慮了推進劑質(zhì)量減少而導致質(zhì)心移動的復雜情況。文獻[10-11]采用改進歐拉中點辛差分格式、模態(tài)疊加方法和自適應(yīng)Newmark法，研究變質(zhì)量系統(tǒng)動態(tài)的動態(tài)響應(yīng)，證明了自適應(yīng)Newmark法對時變系統(tǒng)的適應(yīng)性。

本文在國家自然科學基金項目“自旋飛行器橫向振動的動力學與控制機理研究”工作基礎(chǔ)上，以Timoshenko梁為模型，考慮陀螺效應(yīng)和剪切效應(yīng)，基于有限元法和穩(wěn)定模態(tài)基底法[12]建立變質(zhì)量系統(tǒng)振動方程，分析質(zhì)量變化、推力作用和自旋轉(zhuǎn)速對飛行器穩(wěn)定性的影響，采用自適應(yīng)Newmark法求解系統(tǒng)的振動響應(yīng)。

1 動力學方程

1.1 動能、勢能和外力功

圖1 彈體模型圖

如圖1所示，建立準彈體坐標系Oxyz和隨體坐標系O′ξηζ，忽略軸向變形，系統(tǒng)的動能為[3]：

(1)

式中：uy，uz分別為梁截面y向，z向的橫向位移，θy，θz分別為該截面的轉(zhuǎn)角，Ω為自旋轉(zhuǎn)速，A為截面面積，ρ為軸段微元質(zhì)量密度，I為軸段微元截面的慣性矩，lb為梁長。

考慮剪切影響，系統(tǒng)的彈性勢能為[3]：

(2)

式中：EI為彎曲剛度，κGA為剪切剛度，上標符號“′”表示變量對x的偏導數(shù)。

隨動推力做的功[3]：

(3)

軸向力PN：

(4)

在尾部，隨動推力非保守部分做的虛功[3]：

δWP=Pθz(0,t)δuy(0,t)-Pθy(0,t)δuz(0,t)

(5)

1.2 振動方程

用有限元方法，將彈體模型沿軸向進行離散，劃分為n個梁單元。采用Timoshenko梁單元對橫向位移和轉(zhuǎn)角插值，即：

(6)

由于燃料消耗，質(zhì)量隨時間變化，Φ和Ψ與時間有關(guān)。采用文獻[12]中提出的穩(wěn)定模態(tài)基底法：

(7)

將式(7)代入式(6)中：

(8)

式中：η1(t)=B(t)q1(t)，η2(t)=B(t)q2(t)。

非保守系統(tǒng)Lagrange拉格朗日方程的一般形式如下：

(9)

式中：η和Q分別為廣義坐標和包括非保守力在內(nèi)的廣義力。取η分別為η1和η2,由此得到單元的振動方程：

(10)

K=K1-KPN1-KP1

(11)

(12)

(13)

(14)

(15)

(16)

(17)

(18)

2 數(shù)值仿真及結(jié)果分析

以長徑比為25.3的等截面細長回轉(zhuǎn)梁模型作為仿真模型，整體分為5段，推力作用在軸段1的左端，初始時刻相關(guān)參數(shù)見表1。圖2和3分別給出主動段初至末時刻單位長度質(zhì)量及0～3 s燃料消耗。

表1回轉(zhuǎn)梁模型參數(shù)

圖2 單位長度質(zhì)量

圖3 主動段飛行器質(zhì)量

將式(10)簡化為：

(19)

將式(19)寫成狀態(tài)方程形式：

(20)

2.1 頻率分析

由圖4知，隨著質(zhì)量的消耗，系統(tǒng)的固有頻率逐漸升高；與無推力時相比，推力能夠降低系統(tǒng)固有頻率。

圖4 固有頻率變化

圖5 不同推力作用下進動頻率的變化

圖5分析了質(zhì)量和恒定推力對進動頻率的影響。取恒定推力P=3.2×106N，Ω=102 rad·s-1，推力使正進動頻率減少，負反進動頻率減少，隨著質(zhì)量消耗，一階和二階正進動頻率和負反進動頻率都升高。在恒定推力和質(zhì)量變化共同作用下，一階進動頻率變化斜率越來越小，在t>2.5 s，接近零，進動頻率趨于不變。

圖6 推力變化

圖7 不同工況下進動頻率變化

圖7分析了質(zhì)量和變推力對進動頻率的影響，推力變化如圖6所示。取Ω=102 rad·s-1，由圖7(a)知，工況1推力增大時，一階正進動頻率和負反進動頻率先增大后減小，整體呈增大趨勢，說明這段時間質(zhì)量消耗對進動頻率影響大于推力作用，工況2在推力增大段，一階正進動頻率和負反進動頻率先增大后減小，整體呈減小趨勢，推力作用更明顯；工況1推力不變和減小時，一階正進動頻率和負反進動頻率逐漸正增大，推力減小時，進動頻率變化斜率要大于推力不變時。由圖7(b)知，二階正進動頻率和負反進動頻率逐漸增加，進動頻率變化斜率受推力變化影響。

圖8分析了質(zhì)量和推力對臨界轉(zhuǎn)速ωcr的影響，ωcr采用文獻[14]中方法計算得到。無推力時，隨著質(zhì)量消耗，臨界轉(zhuǎn)速逐漸增大；恒定推力作用時，臨界轉(zhuǎn)速變小，臨界轉(zhuǎn)速變化斜率越來越小，在t>2.5 s，接近零，臨界轉(zhuǎn)速率趨于不變；工況1和2推力增大時，臨界轉(zhuǎn)速先增大后減小，工況1推力不變和減小時，臨界轉(zhuǎn)速逐漸增大。

圖8 臨界轉(zhuǎn)速

圖9 質(zhì)心位置

圖9中，初時刻，質(zhì)心到尾部的距離為2.64 m，隨著質(zhì)量減少，質(zhì)心向頭部移動，移動速度變大，末時刻，到尾部的距離變?yōu)?.5 m。

2.2 振動響應(yīng)分析

經(jīng)典Newmark法的迭代計算公式為：

(21)

(22)

在時變系統(tǒng)的振動問題中，參數(shù)隨時間變化，采用自適應(yīng)Newmark法[10]進行求解，步長h和β為：

(23)

(24)

(25)

(26)

式中：h(t)為自適應(yīng)步長;κ為步長控制因子；ω(t)為最高階振動角頻率；ξ為系統(tǒng)模態(tài)阻尼比。

采用自適應(yīng)Newmark法分析由激振力引起的變質(zhì)量系統(tǒng)橫向振動響應(yīng)問題，激振力作用于飛行器尾部y,z向，大小均為2000 N，激振頻率與自旋轉(zhuǎn)速相同。

圖10 不同轉(zhuǎn)速下尾部位移響應(yīng)

分別計算自旋轉(zhuǎn)速Ω為90 rad·s-1,101 rad·s-1,112 rad·s-1時橫向位移響應(yīng)情況。由圖8知，恒定推力P=3.2×106N時，ωcr在初和末時刻分別為95.8 rad·s-1和107.2 rad·s-1。圖10(a)中自旋轉(zhuǎn)速小于ωcr，受質(zhì)量消耗影響，尾部振幅減少；圖10(b)中振幅逐漸增大，在1.8 s時振幅達到最大，后逐漸減小，由圖8可知，在1.1 s時，自旋轉(zhuǎn)速等于臨界轉(zhuǎn)速，系統(tǒng)發(fā)生短暫共振,共振響應(yīng)峰值最大時刻與自旋轉(zhuǎn)速和臨界轉(zhuǎn)速重合時刻相比延遲0.7 s；圖10(c)中質(zhì)量消耗，臨界轉(zhuǎn)速增加，趨向于自旋轉(zhuǎn)速，導致尾部振幅逐漸增加。

3 結(jié) 論

以柔性自旋飛行器為研究對象，采用變質(zhì)量Timoshenko梁為模型，基于有限元法和穩(wěn)定模態(tài)基底法建立時變系統(tǒng)橫向振動方程，通過算例仿真，分析變質(zhì)量、推力作用和轉(zhuǎn)速對飛行器穩(wěn)定性響，采用自適應(yīng)Newmark法，對系統(tǒng)橫向振動響應(yīng)進行分析，得到以下結(jié)論：

1)隨著質(zhì)量減少，系統(tǒng)固有頻率、進動頻率和臨界轉(zhuǎn)速逐漸增加，推力作用能夠減小系統(tǒng)剛度，不考慮質(zhì)量變化時，正進動和負反進動頻率減小，推力和質(zhì)量共同作用時，正進動和負反進動頻率增大，說明質(zhì)量消耗對進動頻率影響大于推力作用。

2)當自旋轉(zhuǎn)速小于臨界轉(zhuǎn)速時，受質(zhì)量影響振幅減小；靠近臨界轉(zhuǎn)速時，振幅增加，等于臨界轉(zhuǎn)速頻率時發(fā)生短暫共振，共振響應(yīng)峰值最大時刻與轉(zhuǎn)速重合時刻相比稍有延遲。