張夢櫻,唐乾剛,韓小軍,張青斌,葛健全
(1. 國防科學技術大學 航天科學與工程學院,湖南 長沙410003;2. 第二炮兵裝備部 高新技術辦公室,北京100085)
隨著對高超聲速滑翔飛行器研究的不斷深入,滿足禁飛區約束的滑翔軌跡優化問題逐漸成為了一個研究熱點。滿足飛行安全的同時,需兼顧復雜航路禁飛區約束,其軌跡優化問題將變得更為復雜。
在禁飛區建模方面,文獻[1]建立了基于地理位置的威脅模型,將暴露危險作為優化目標設計最優軌跡。文獻[2]考慮了威脅和目標可動,建立數學模型并采用遺傳算法在線優化軌跡。文獻[3]在平面地球假設下,建立無限長圓柱的禁飛區模型,設計了滿足航路點和禁飛區約束的最優軌跡。文獻[4]的禁飛區模型建立在球形地球表面,采用球面三角形原理建立無限長圓柱的空間禁飛區模型。文獻[5]建立了半橢球的禁飛區模型,但未改變橢球赤道半長軸方向以適用更一般的情況。現有航路約束或禁飛區約束建模并不完全適用于實際情況,為使最優軌跡能夠滿足任務需求,有必要考慮更為復雜和真實的禁飛區模型。
復雜約束條件下的高超聲速軌跡優化問題可歸結為強非線性、多階段、多約束的最優控制問題[6],對其求解十分困難。因具有收斂速度快和求解精度高的特點,偽譜法近年來發展迅速,成為求解復雜最優控制問題的常見工具。文獻[5,7 -9]用Gauss 偽譜法求解包括端點、航路點和禁飛區約束的通用航空飛行器(CAV)軌跡優化問題。國內學者也研究了包括過載、動壓、駐點熱流等多約束軌跡優化問題[10-14]。在處理復雜路徑約束問題時國內外研究者多采用分段優化的方法[5,8-9,11],有些會采用HP自適應偽譜方法求解[14-17]此類問題。然而當Gauss 節點較多時,處理復雜約束問題會存在設計變量數目龐大和約束眾多,從而影響求解魯棒性的問題[9,11]。同時,采用分段方法,也將增加設定設計變量初值的難度。
針對帶有復雜禁飛區約束的高超聲速滑翔飛行器軌跡優化問題直接采用偽譜法難以求解的問題,本文提出依序迭代求解策略。
為簡化計算,設地球為無旋圓球,在地心坐標系下建立飛行器的三自由度再入動力學模型[18]:


氣動力系數通常是馬赫數Ma 和攻角α 的函數,但因高超聲速飛行器的氣動遵循馬赫數無關原則,可采用(3)式擬合,相應氣動參數詳見文獻[15].

1.2.1 目標函數
對于指定目標點的再入軌跡優化問題,飛行時間是影響突防性能一個重要的優化目標。

式中:t0和tf分別為起始時間和終止時間。
飛行器在大氣中的運動微分方程采用(1)式。
1.2.2 約束模型
考慮飛行器任務需求,需要考慮飛行器的指定端點約束,同時確保熱防護安全、結構可靠性和姿態穩定性。
1)端點約束:

2)氣動熱約束[21]:

3)過載約束
再入段總過載通常較小,主要考慮法向過載ny的約束。

4)動壓約束:

1.2.3 控制量選取
攻角α 和側傾角σ 是影響飛行軌跡的重要變量,選取它們作為優化問題的控制變量,主要限制二者的取值范圍,其表達式為

禁飛區是指飛行過程中不允許經過的空域,如領空、防空識別區、雷達探測區和反導防御空域等。本文按邊界地面投影類型分為規則和一般兩類。規則禁飛區包括半球形、半橢球、橢圓柱和圓柱形,可建立參數化數學模型描述其覆蓋范圍。對一般禁飛區難于數學描述,可轉化為若干規則型禁飛區的疊加組合。因此,典型規則形狀禁飛區的數學建模是基礎。
如圖1所示,將地球球形表面某點(rnfz,λnfz,φnfz)作為禁飛區中心,以其為坐標原點建立禁飛區坐標系Ox1y1z1,其中Ox1軸指向禁飛區在過O 點切面上投影的短半軸方向,Oy1指向天,Oz1指向長半軸方向。αnfz為投影短半軸和正北方向的夾角,以順時針方向為正。OeXeYeZe為地心直角坐標系,(r,λ,φ)為飛行器軌跡上某點坐標。
將以上定義的位置的球坐標轉換到地心坐標系下:

圖1 禁飛區模型示意圖Fig.1 Schematic diagram of no-fly zone model

式中:(x,y,z)和(xnfz,ynfz,znfz)分別為飛行器和禁飛區中心在地心坐標系下的坐標。
從地心系轉化到禁飛區坐標系Ox1y1z1的轉化矩陣為

式中:M2、M1、M3表達式分別為

則飛行器在禁飛區坐標系下的位置坐標為

橢球禁飛區約束可描述為

式中:a 為橢球赤道長半軸;b 為短半軸;c 為橢球的極半徑。
特別地,當a=b=c 時,橢球禁飛區退化成圓球禁飛區。
當c=inf 時,(18)式退化為

表示一個無限高橢圓柱禁飛區約束。
特別地,當{rnfz=re}∩{αnfz=0}∩{a =b}時,禁飛區退化成圓柱禁飛區,亦可用球面三角表示為

式中:rnfz為禁飛區在切平面上投影的半徑,rnfz=a =b.
近年來,偽譜法因具有收斂速度快和求解精度高的特點,成為求解復雜最優控制問題的常見工具[8]。影響偽譜法收斂速度和計算精度的主要因素有參數初值和離散節點的數量。由于偽譜法具有以較少的節點數獲得較精確的解的優點,故參數初值成為主要因素[9]。偽譜法現廣泛應用于軌跡優化問題中,但其算法性能仍受初值的影響[9]。
考慮多約束的高超聲速飛行器軌跡優化問題較為復雜,直接運用偽譜法求解將很難得到優化結果。同時,第2 節所建立的禁飛區模型的復雜程度也給直接優化增添難度。為解決上述問題,彌補算法對初值的敏感性并提高計算效率,考慮分段優化之外的方法,即采用如圖2所示的基于Gauss 偽譜法的最優軌跡依序迭代串行求解流程。
依據各約束的復雜程度,逐步施加約束,并將優化成功得到的彈道作為下一步的初值彈道進行迭代。具體描述如下:

圖2 優化流程Fig.2 Flow chart of optimization
步驟1 生成初值彈道。根據任務設定起點和終點的位置速度參數,擬定初始方向角和速度傾角,在GPOPS 軟件中設定初值猜測范圍、迭代精度和節點數,開始進行求解。若滿足約束,則輸出彈道。若不滿足,則返回調整初值,直到得到優化結果。
步驟2 基于Gauss 偽譜法求滿足端點約束和過程約束的最優解。優化模型中添加包括動壓、駐點熱流以及過載等過程約束,設置迭代精度和節點數,并將上一步所得優化結果作為下一步增加了約束的優化問題的初值彈道進行迭代求解。若結果滿足約束,則輸出彈道,若不滿足,首先調整節點數,不能滿足則返回上一步重新調整計算。
步驟3 求解滿足端點約束、過程約束和禁飛區約束的最優彈道。如圖3所示在優化模型中添加禁飛區約束。由2.3 節可得,所提出的禁飛區模型均可由橢球禁飛區模型退化變形得到,故首先添加半橢球禁飛區模型。為提高求解效率,加快收斂,考慮如下逐步擴大禁飛區大小的策略:

圖3 禁飛區約束處理方法Fig.3 Process of dealing with no-fly zones constraint
設禁飛區加權的長半軸、短半軸及極半徑分別為

將每一步優化得到的可行解迭代到下一步,得到滿足禁飛區約束的軌跡。如在這一步始終不能得到最優解,則返回步驟2 調整初值,重新優化彈道。
步驟4 檢查結果是否滿足所有約束,并判斷收斂性,輸出結果彈道。
根據前述的飛行器運動模型、優化模型和優化算法等,設計如下3 個算例:算例1 通過對經典算例的仿真計算驗證了算法的正確性;算例2 對比不同禁飛區模型對優化結果的影響;算例3 對多組合禁飛區約束優化問題進行仿真,驗證禁飛區模型的正確性和所提出計算方法的有效性。
4.1.1 算例1:算法驗證
文獻[19]所提出的最大橫程算例,是典型最優控制問題的算例。設置初始參數如表1所示,采用所提出優化策略進行計算。

表1 仿真參數設置Tab.1 The values of simulation parameters
4.1.2 算例2:4 種禁飛區模型優化結構比較
在相同初始運動參數和端點約束條件下,以相同地理坐標為中心,規避某一特定區域,設置如下4 組禁飛區,以最小時間為優化目標,研究不同禁飛區模型對優化彈道結果的影響。為了嚴格比較4 種禁飛區模型優化結果,應用GPOPS 求解時,設置彈道分段數為10 段,數值計算精度為10-3,每段節點數最小為4,最大為8.
端點和禁飛區參數設置如表2、表3所示。

表2 端點參數設置Tab.2 The values of start point and end point

表3 不同禁飛區約束分組設置Tab.3 The values of different group of no-fly zone constraints
4.1.3 算例3:多組合禁飛區約束優化問題
考慮禁飛區約束的滑翔彈道優化問題,低速飛行器將其考慮成一個低空突防的航路規劃問題[21]。文獻[22]采用了遺傳算法、航路極坐標編碼方式、航跡適應度函數和遺傳操作算子法進行仿真,卻不適合高超聲速條件下的軌跡優化問題。為滿足禁飛區約束,采用Gauss 偽譜法進行彈道優化設計。
設計一條固定端點的時間最短飛行軌跡,在飛行過程中需繞過兩個禁飛區。第1 個禁飛區模擬行政區域的領空,第2 個禁飛區模擬航母戰斗群的防御半徑。具體仿真數值設置如表4所示。

表4 多約束條件設置Tab.4 Multi-constraints
4.2.1 算例1
本例仿真運算中,設置彈道分段數為10 段,精度為10-7,每段節點數最小為4,最大為8. 將所得=∞優化結果作為初值,迭代到添加熱流密度約束的問題中,所得優化結果與文獻[19]比較如表5所示。

表5 優化結果比較Tab.5 Comparison of optimal results
本文優化所得結果中,無約束結果較好復現了文獻[19]的結果。比較添加熱流約束的結果,仿真得到全程駐點最大熱流為795.14 kW/m2,考察熱流極值約束滿足情況,其誤差為0.020 8%,滿足精度要求。與文獻[19]結果相比,在滿足相同約束的情況下,優化目標最大橫程更大,因而結果更優。同時,若采用直接添加約束進行優化,其仿真時間為38.70 s,相比迭代結果要慢。綜上,可驗證本文所提出方法的正確性。
4.2.2 算例2
在用偽譜法難以直接求解原問題的情況下,采用本文策略得到4 組優化時間分別為2 021.3 s、1 944.8 s、1 833.8 s 和1 828.5 s,對應圖4第2 ~5 條軌跡結果。比較可得,球形禁飛區模型所得優化時間短于柱形禁飛區模型所得優化時間,且在4 種類型禁飛區均能滿足規避如圖陰影區域的情況下,橢球禁飛區對應時間最短,即結果最優。

圖4 優化軌跡局部投影Fig.4 Part ground projection view of optimal trajectories
4.2.3 算例3
類似地,由于本例約束復雜,偽譜法直接求解無法得到結果,采用3.2 節提出策略。運算中,設置彈道分段數為10 段,精度為10-3,每段節點數最小為4,最大為8. 取攻角α 和傾側角σ 為優化變量,得到優化軌跡及其地面投影和結果如圖5~圖7所示。

圖5 三維優化軌跡Fig.5 3-D optimal trajectory

圖6 優化軌跡地面投影Fig.6 Ground projection of optimal trajectories
圖6為無約束和僅添加過程約束的兩條軌跡穿過禁飛區,而添加禁飛區約束的優化結果滿足所設置的路徑約束,即從起點出發,先后繞過兩個禁飛區的側面實現規避。結果反映了4 條軌跡逐步逼近原問題最優解的過程。由圖7(e)可見,控制量的變化滿足了過程約束要求,且在滑翔段,攻角保持在17°附近,即對應最大升阻比攻角。由于所需規避區域距離較遠,軌跡緯度變化平緩(見圖7(c)),其地面投影呈不甚明顯的“S”型,因而傾側角在飛行過程中并未出現較大翻轉,且對應規避動作時,角度變化率較為穩定,在工程上易于實現。由圖7(h)~圖7(j)可知,飛行器的駐點熱流密度、法向過載及動壓均滿足了約束條件。同時,由圖7可知飛行器終端參數λf=62.42°,φf=42.28°,hf=20 km,vf=712.62 m/s,γf= -10°,均滿足設計要求。
在Windows XP 系統下,Matlab2009 計算平臺下,采用第3 節描述的優化策略,從初值計算到優化結果的得到,仿真計算總耗時為4 036.1 s,約1 h,易于仿真實現。
本文研究了復雜約束條件下的再入軌跡快速優化問題,從復雜約束的解構入手,建立了常用過程約束模型和更符合實際情況的典型禁飛區約束模型;針對采用Gauss 偽譜法所面臨的初值敏感問題,提出依據復雜程度逐步添加約束的迭代求解策略;通過與文獻[19]算例進行對比,驗證了方法的正確性;分析了不同禁飛區模型適用的實際情況,并以仿真算例驗證所建立模型的準確性;以時間最小為例設計了固定端點并考慮復雜約束條件下的最優軌跡。結果表明,依序逐步添加約束的迭代求解策略,不但能方便解決考慮多種禁飛區約束的問題,并且避免了Gauss 偽譜法在初值彈道生成上的難題,可作為分段求解軌跡優化問題方法的補充和參照。

圖7 優化結果Fig.7 Optimal results
References)
[1]Vian J L,Moore J R. Trajectory optimization with risk minimization for military aircraft[J]. Journal of Guidance,Control,and Dynamics,1989,12(3):311 -317.
[2]Twigg S,Calise A,Johnson E. On-line trajectory optimization including moving threats and targets[C]∥AIAA Guidance,Navigation,and Control Conference and Exhibit. Providence,Rhode Island:American Institute of Aeronautics and Astronautics,2004:2004 -5139.
[3]Jorris T R,Cobb R G. Three-dimensional trajectory optimization satisfying waypoint and no-fly zone constraints[J]. Journal of Guidance,Control,and Dynamics,2009,32(2):551 -572.
[4]陳小慶,侯中喜,劉建霞. 高超聲速滑翔式飛行器再入軌跡多目標多約束優化[J]. 國防科技大學學報,2009,31(6):77-83.CHEN Xiao-qing,HOU Zhong-xi,LIU Jian-xia. Multi-objective optimization of reentry trajectory for hypersonic glide vehicle with multi-constraints[J]. Journal of National University of Defense Technology,2009,31 (6):77 -83. (in Chinese)
[5]謝愈,劉魯華,湯國建. 多約束條件下高超聲速滑翔飛行器軌跡優化[J]. 宇航學報,2011,32 (12):2499 -2504.XIE Yu,LIU Lu-hua,TANG Guo-jian. Trajectory optimization for hypersonic glide vehicle with multi-constraints[J]. Journal of Astronautics,2011,32(12):2499 -2504. (in Chinese)
[6]李惠峰. 高超聲速飛行器制導與控制技術[M]. 北京:中國宇航出版社,2012:242 -243.LI Hui-feng. Guidance and control technology of hypersonic vehicles[M]. Beijing:China Astronautic Publishing House,2012:242 -243.(in Chinese)
[7]Yaple D E. Simulation and application of GPOPS for a trajectory optimization and mission planning tool[D]. Ohio,US:Air Force Institute of Technology,2010.
[8]Jorris T. Common aero vehicle autonomous reentry trajectory optimization satisfying waypoint and no-fly zone constraints[D].Ohio,US:Graduate School of Engineering and Management,Air Force Institute of Technology,2007.
[9]雍恩米. 高超聲速滑翔式再入飛行器軌跡優化與制導方法研究[D]. 長沙:國防科學技術大學,2008.YONG En-mi. Study on trajectory optimization and guidance approach for hypersonic glide-reentry vehicle[D]. Changsha:National University of Defense Technology,2008.(in Chinese)
[10]姚寅偉,李華濱. 基于Gauss 偽譜法的高超聲速飛行器多約束三維再入軌跡優化[J].航天控制,2012,30(2):33 -45.YAO Yin-wei,LI Hua-bin. The generation of three-dimensional constrained entry trajectories for hypersonic vehicle based on the gauss pseudospectral method [J]. Aerospace Control,2012,30(2):33 -45.(in Chinese)
[11]Zhao J,Zhou R. Reentry trajectory optimization for hypersonic vehicle satisfying complex constraints[J]. Chinese Journal of Aeronautics,2013,26(6):1544 -1553.
[12]Sun Y,Zhang M. Optimal reentry range trajectory of hypersonic vehicle by Gauss pseudospectral method[C]∥International Conference on Intelligent Control and Information Processing. Harbin:IEEE,2011:545 -549.
[13]Liu P,Zhao J S,Gu L X. Rapid approach for three-dimensional trajectory optimization of hypersonic reentry vehicles[J]. Flight Dynamics,2012,30(3):263 -268.
[14]Darby C L,Hager W W,Rao A V. An HP-adaptive pseudospectral method for solving optimal control problems[J]. Optimal Control Applications and Methods,2011,32(4):476 -502.
[15]國海峰,黃長強,丁達理,等. 多約束條件下高超聲速導彈再入軌跡優化[J].彈道學報,2013,25(1):10 -15.GUO Hai-feng,HUANG Chang-qiang,DING Da-li,et al. Reentry trajectory optimization for supersonic missile considering multiple constraints [J]. Journal of Ballistics,2013,25(1):10 -15.(in Chinese)
[16]王麗英,張友安,趙國榮. 改進的HP 自適應網格細化算法及應用[J].彈道學報,2013,25(1):16 -21.WANG Li-ying,ZHANG You-an,ZHAO Guo-rong. Improved HP-adaptive mesh refinement algorithm and its application[J].Journal of Ballistics,2013,25(1):16 -21.(in Chinese)
[17]洪蓓,辛萬青. HP 自適應偽譜法在滑翔彈道快速優化中的應用[J].計算機測量與控制,2012,20(5):1283 -1286.HONG Bei,XIN Wan-qing. Application of HP-adaptive pseudospectral method in rapid gliding trajectory optimization[J]. Computer Measurement & Control,2012,20(5):1283 -1286.(in Chinese)
[18]趙漢元. 飛行器再入動力學和制導[M]. 長沙:國防科學技術大學出版社,1997:84 -85.ZHAO Han-yuan. Reentry dynamics and guidance of flight vehicle[M]. Changsha:National University of Defense Technology Press,1997:84 -85.(in Chinese)
[19]Betts J T. Practical methods for optimal control and estimation using nonlinear programming[M]. Philadelphia:Society for Industrial and Applied Mathematics,2009:247 -251.
[20]賈沛然,陳克俊,何力. 遠程火箭彈道學[M]. 長沙:國防科學技術大學出版社,2009:155 -156.JIA Pei-ran,CHEN Ke-jun,HE Li. Long-distance missile ballistics[M]. Changsha:National University of Defense Technology Press,2009:155 -156.(in Chinese)
[21]文葉,姜文志,馬登武. 飛機低空突防中的航路規劃技術研究[J]. 飛行力學,2006,24 (2):22 -26.WEN Ye,JIANG Wen-zhi,MA Deng-wu. Study on the plane low-altitude penetration route planning technique[J]. Flight Dynamics,2006,24 (2):22 -26.(in Chinese)
[22]熊丹君,蔡滿意,劉宇坤. 多約束條件下飛行器航路規劃[J]. 彈箭與制導學報,2009,29 (2):289 -292.XIONG Dan-jun,CAI Man-yi,LIU Yu-kun. Air vehicle route planning for multi-constraints[J]. Journal of Projectiles,Rockets,Missiles and Guidance,2009,29 (2):289-292.(in Chinese)