黃育秋,何麟書
(北京航空航天大學宇航學院,100191北京,huangyuqiu13@163.com)
傳統的彈道式再入飛行器和彈道-升力式再入飛行器存在兩大缺點,即存在造成再入器及其有效載荷損傷的著陸沖擊過載和由于各種干擾造成的不易控制的大落點散布[1].升力式再入飛行器的出現有效的克服了這兩種缺點,是航天技術取得的巨大進步.升力式再入飛行器和傳統的再入飛行器相比,升力式再入飛行器升力的增大和可調整,大大增加了飛行器機動飛行的能力.平緩的再入段和大范圍的機動飛行能力,使升力式再入飛行器水平著陸到指定機場跑道和實現全球快速打擊成為可能,正因為如此,升力式再入飛行器具有重要的經濟和軍事意義.
1940年末到1950年初,蘇聯進行了Silbervogel飛行器的風洞試驗,并積累豐富的試驗數據. 1960年中期,Mikyan設計局設計了自己的升力式再入飛行器Mig-105.最近,俄羅斯已秘密研制了新型再入機動飛行的白楊-M導彈[2].
由于升力式再入飛行器的飛行速度較快,飛行空域較大,飛行環境復雜,其彈道優化需要考慮控制、動壓、過載、氣動熱等限制,成為總體設計中的1個難點.且隨著導彈防御系統的日臻完善,再入飛行器的彈道設計還必須考慮實際作戰需求.再入飛行器要順利地飛向目標,必須能夠避開勿入區域,突破敵方攔截區域,且經過設定路徑點,如圖1所示.
勿入區域是指由于地理或政治等原因,再入飛行器不得飛入的領空.路徑點是指為了搜索目標或制導導航等需求而設定的彈道必經地理位置.再入飛行器必須精確飛經此點,且時間可不做約束.攔截區域是指敵方攔截彈等的殺傷范圍.在攔截導彈殺傷區域內,攔截彈可以捕獲再入飛行器,預測其飛行彈道,并實施攔截.在此攔截區域內,再入飛行器可進行一定的規避機動,以降低敵方對彈道的預測精度,減小攔截概率.再入飛行器的規避機動可根據攔截導彈的飛行能力和敵方攔截態勢等選擇適當時機進行[2-3].

圖1 突防彈道示意
由圖1可見,考慮突防約束情況下,再入飛行器的再入機動突防彈道主要分為:飛抵第1個路徑點階段、繞過勿入區域飛抵第2個路徑點階段、滑翔至敵方導彈攔截區階段、攔截區內機動飛行至目標上空階段、下壓攻擊目標階段.在這些限制條件下,升力式再入飛行器的再入機動突防彈道優化設計實際上是1個復雜的多約束多階段彈道優化問題.
本文首先建立了升力式再入飛行器再入機動突防的多約束多階段彈道優化模型,研究了飛行器軌跡優化問題的數值解法.在此基礎上選擇直接法將該軌跡優化問題轉化成參數優化問題,而后采用序列二次規劃法來解該參數優化問題,最后采用C++語言編寫了優化算法,在特定參數下,進行了仿真分析.
考慮地球自轉和地球扁率的情況下,再入飛行器的動力學模型為[1]


其中:V、γ、ψ、r、θ、φ、σ、α、P、X、Y、ωe、g'r、gωe、m分別為飛行器相對地球的速度、飛行路徑角、航向角、地心距、經度、緯度、傾斜角、攻角、推力、氣動阻力、氣動升力、地球自轉角速度、引力加速度在地心矢徑r上的投影、引力加速度在地球自轉角速度ωe上的投影、飛行器質量.航向角ψ是速度矢量V在當地水平面上的投影線順時針與緯度切線正東方向的夾角.
升力式再入飛行器的飛行環境非常復雜,涉及諸多實際問題,彈道設計時必須綜合考慮諸多實際約束,同時其機動突防彈道還需考慮如圖1中的突防約束,下面將分別闡述.
1)終端約束.升力式再入飛行器是一種遠程精確制導武器,用以摧毀敵方高價值目標.終端約束用以表示其攻擊目標時的狀態等,相應的終端約束取為

其中:h為距地高度;tf為飛行器飛行總時間,是一優化參數;C1、C2、C3為根據打擊需要而給定的常數.同時,給定了地面目標的具體地理位置時,相應的終端位置約束為

其中θf表示目標經度,φf表示目標緯度.
2)控制量約束.由于飛行器結構與姿態控制系統設計的限制,飛行器飛行時的攻角不能過大,且變化不能過于劇烈,以免飛行失控.同理,為了滿足飛行器的控制要求,其攻角與傾側角的大小及其變化率也應加以限制.

其中:t0為飛行器再入初始時刻,可取為零;˙α表示攻角變化率;˙σ表示傾側角變化率;C4~C7為給定常數.
3)動壓約束.動壓極限值主要取決于熱防護材料強度與氣動控制鉸鏈矩.從防熱系統設計來說,飛行器表面均采用耐高溫絕熱材料,以保證飛行器飛行過程中內部結構所受到的加熱量最小和在高溫加熱時保持應有的氣動外形.這些材料直接面對來流作用,因此,動壓必須限制在一定范圍內,以確保表面絕熱材料結構不受破壞.氣動控制鉸鏈力矩隨動壓的增加而增大,動壓也應保證不超過控制氣動操縱面所要求的最大鉸鏈力矩所允許的動壓.同時對動壓加以限制也可以在一定程度上保證飛行器側向飛行穩定.
因此,為了滿足升力式再入飛行器的結構設計與控制要求,相應的動壓q約束為

式中C8為給定常數.
4)法向過載約束.法向過載最大值主要取決于飛行器的結構強度和彈載設備的承受范圍.為了滿足臨近空間飛行器的結構設計要求,相應的法向過載ny約束為

式中C9為給定常數.
5)氣動熱約束.升力式再入飛行器的長時間大氣層內高速飛行將產生非常高的氣動加熱率與壁溫.由于升力式再入飛行器一般采用小攻角飛行,氣動加熱最嚴重的部位為飛行器的鼻頭與翼前緣.為了便于研究,彈道優化時可以僅考慮鼻頭氣動熱問題.為確保氣動熱不超過壁面材料容忍限度,飛行器鼻頭駐點熱流與壁溫約束取作

其中C10、C11為給定常數.
6)突防約束.如前所述,由于導彈防御系統的出現,升力式再入飛行器的彈道優化必須考慮實際作戰需求.升力式再入飛行器要順利地飛向目標,必須能夠避開勿入區域,突破敵方攔截區域,且經過設定路徑點.
為避免進入敵方領空或攔截范圍,飛行器不得進入勿入區域飛行,將此約束取為

其中0≤t≤tf.攔截導彈殺傷區可取作地理約束,即

其中tM0≤t≤tMf.此殺傷區內的機動形式可以設定為“S”形側向程序機動,相應的傾側角變化規律取為

其中:tM0≤t≤tMf;TM=(tMf-tM0),tM0表示飛行器飛入攔截導彈殺傷區的時刻,是一優化參數,tMf表示程序機動結束時刻,也指俯沖段開始時刻,是一優化參數.式(1)中,等號右端第1項為未加入姿態擾動時的傾側角變化規律;第2項為相應的姿態擾動,a1為機動幅值,即飛行器在導彈攔截區域內進行2個最大幅度為a1的正弦機動.
為了保證升力式再入飛行器準確經過指定路徑點,將路徑點約束取為

其中tpass1、tpass2分別表示經過第1個、第2個路徑點的時刻,是優化參數.以上公式中C12~C21為給定常數[3].
同時,現代防空武器的有效攔截高度多在30 km以下.因此,為了更好地突破敵方防空火力網,飛行器在俯沖攻擊之前,飛行高度也應予以限制,即h≥30 km,t0≤t≤tD0.其中tD0為俯沖前某特定時刻.
多階段優化問題是指,優化模型包含多個狀態階段,各階段之間依靠時間、狀態量銜接,并同時進行優化.多約束多階段優化問題可將多個階段的彈道優化模型統一于1個優化算法下,同時進行解算,從而提高了優化模型的通用性.本文采用直接法將該軌跡優化問題轉化成參數優化問題,而后采用序列二次規劃法來解該參數優化問題,下面將闡述其相關理論.
軌跡優化問題一般可以描述為,確定容許控制u(t)∈Rm和參量p∈Rnp,使得由1個微分方程組確定的系統,從給定的初始狀態過渡到終端狀態,并使性能指標函數J達到最小,同時滿足規定的約束.其數學描述如下[4-8]:

滿足

其中x(t)∈Rn,表示系統的狀態變量,t表示時間變量.標量性能指標函數J,由末值型性能指標函數Φ(x(tf),p)和積分型性能指標函數組成,其被積函數為L(x,u,p,t),并且積分是從t0時刻到tf時刻.方程組(3)表示系統狀態方程,方程組(4)表示狀態變量、控制變量和參量的等式約束,方程組(5)表示狀態變量、控制變量和參量的不等式約束,方程組(6)表示狀態變量的初始條件,方程組(7)表示狀態變量和參量的終端條件.
1)劃分時間區間[t0,tf]為N個子區間,節點為ti,即t0<t1<… <tN-1<tN=tf.
2)控制變量的參數化.在每個子區間t∈[ti,ti+1](i=0,1,…,N-1)里,將控制變量近似為

其中:ui∈Rm,表示在ti時刻的控制變量值;ui+1∈Rm,表示在ti+1時刻的控制變量值.利用分段線性近似方程(8),未知的控制變量u(t)被m(N+1)個未知的控制參數u0,u1,…,uN-1,uN代替,因此,所有的未知參數可以組成1個向量~u∈Rm(N+1)+np.
3)價值函數、約束函數的參數化.假設給定了1個猜測的控制輸入量和參量,即~u,在初始條件方程(6)下,從t0時刻到tf時刻積分狀態方程組(3),得到的狀態變量隨時間的變化歷程可以表示為x(~u,t),也就是說,利用控制輸入量和參量可以唯一地確定狀態量,進而得到價值函數J以及約束c、d、ψ,據此約束可以離散成1個等式約束向量g∈Rnc(N+1)+nf和不等式約束向量h∈Rnd(N+1).
用上述的參數化策略,有限維的最優控制問題,(即方程組(2)~(7))被近似化為有限維的非線性規劃問題[9-10],即

方程(9)所表示的非線性規劃問題可以采用如下序列二次規劃方法求解:


式中:


矩陣Bk的初值B0一般取為單位陣,即B0= I[11-14].其算法流程如圖2所示.

圖2 算法流程
取突防約束如表1所示.

表1 突防約束
利用1節中的模型和2節中的優化方法得到仿真結果如圖3~8所示.

圖3 攻角度變化曲線
圖7給出了升力式再入飛行器的最優三維突防彈道,圖8給出了最優彈道的地表投影.可見,升力式再入飛行器從發射點出發,精確地通過了第1個路徑點,順利繞過勿入區域,到達了設定的第2個路徑點,突入了攔截導彈殺傷區,以“S”形機動彈道在攔截導彈殺傷區內飛行以提高突防概率,最終精確命中了目標.同時可以發現,滑翔段飛行高度保持在30 km以上,滿足了突防高度要求.

圖4 傾側角變化曲線

圖5 速度變化曲線

圖6 飛行路徑角變化曲線

圖7 三維機動突防彈道

圖8 經緯度曲線
本文以升力式再入飛行器機動突防彈道優化設計為研究目的,首先給出了該問題的多約束多階段彈道優化模型,研究了彈道優化數值解法理論,將該多約束多階段優化問題的多個階段彈道優化模型統一于1個優化算法.采用直接法將該軌跡優化問題轉化成參數優化問題,而后采用序列二次規劃法來解該參數優化問題,得到了滿足相應約束的再入機動突防彈道.通過本文的研究可以看出采用本文的方法能夠進行升力式再入飛行器的再入機動突防彈道的優化設計,并具有較好的效果.
[1]趙漢元.飛行器再入動力學和制導[M].長沙:國防科技大學出版社,1997:382-384.
[2]李瑜,楊志紅,崔乃剛.助推-滑翔導彈彈道優化研究[J].宇航學報,2008,29(1):67-69.
[3]李瑜.助推-滑翔導彈彈道優化與制導方法研究[D].哈爾濱:哈爾濱工業大學,2009.
[4]BRAIN C F.Some tools for the direct solution of optimal control problems[J].Advances in Engineering Software,1998,29(1):45-61.
[5]LEWALLEN J M,TAPLEY B D,WILLIAMS S D.Iteration procedures for indirect trajectory optimization methods[J].Journal of Spacecraft and Rockets,1968,5 (3):321-327.
[6]VINH N X,CHERN J S,LIN C F.Phugoid oscillations in optimal reentry trajectories[J].Acta Astronautica,1981,8:311-324.
[7]JACKSON M C,STRAUBE T M,Fill T J,et al.Onboard determination of vehicle glide capability for the shuttle abort flight manager(SAFM)[C]//IEEE Aerospace and E-lectronic Systems Society,Core Technologies for Space Systems Conference.Colorado Springs:[s.n.],2002.
[8]Gath P F.CAMTOS-a software suite combining direct and indirect trajectory optimization methods[D].Stuttgart:University of Stuttgart,2002:10-28.
[9]de O PANTOJA J F O,MAYNE D Q.A sequential quadratic programming algorithm for discrete optimal controlproblemswith controlinequality constraints[C]//Proceedings of the 28th conference on decision and control.Tampa,FL:[s.n.],1989,1:353-357.
[10]BARRON R L,CHICK C M.Improved indirect method for air-vehicle trajectory optimization[J].Journal of Guidance,Control,and Dynamics,2006,29(3):643-652.
[11]BETTS J T.Survey of numerical methods for trajectory optimization[J].Journal of Guidance,Control,and Dynamics,1998,21(2):193-207.
[12]SHIPPEY B M.Trajectory optimization using collocation and evolutionary programming for constrained nonlinear dynamical systems[D].Arlington:University of Texas,2008.
[13]VINH N X,LU P.Chebyshev minim ax problems for skip trajectories[J].Journal of the Astronautical Sciences,1988,36(1):179-197.
[14]YONG E,TANG G J,CHEN L.Three-dimensional optimal trajectory for global range of CAV[C]//1st International Symposium on Systems and Control in Aerospace and Astronautics.Piscataway:IEEE,2006:1396-1400.