亞軌道飛行器返回軌跡快速優化

王文虎

中北大學機電工程學院，太原 030051

亞軌道飛行器(Suborbital Launch Vehicle，SLV)是指從地面發射，在亞軌道空間運行作業，然后重返地面的跨大氣層可重復使用飛行器[1]。SLV返回段具有飛行高度低、任務多樣性、飛行環境多變的特點。

快速、準確、魯棒的SLV返回軌跡優化可以實現快速任務規劃從而極大地降低成本。返回軌跡優化本質是求解一類高度復雜的最優控制問題。常見的數值方法主要包括間接法和直接法。近年來，直接法中的偽譜法因其求解精度高、收斂速度快、具有很好的魯棒性而備受關注[2]。偽譜法主要有：勒讓德偽譜法(Legendre Pseudospectral Method, LPM)、拉道偽譜法(Radau Pseudospectral Method, RPM)和高斯偽譜法(Gauss Pseudospectral Method, GPM)。

國內外關于GPM及LPM研究較多，而RPM國內未見相關文獻。Huntington[3]對上述3種偽譜法在計算效率和精度方面進行了比較，但所用的RPM是向后Radau偽譜法，并且算例較為簡單，沒有比較各種偽譜法在處理復雜問題時的能力。目前，國內外還沒有文獻對向前RPM與GPM進行過比較。

本文針對亞軌道飛行器，分別采用GPM和向前RPM進行了SLV返回軌跡快速優化研究，同時對2種偽譜法進行了比較分析。結果表明，Gauss偽譜法在處理含控制量約束的問題時比向前Radau偽譜法更具優勢。

1 Gauss偽譜法與向前Radau偽譜法

1.1 Gauss偽譜法與向前Radau偽譜法基本原理

偽譜法是一種正交配點法，其基本思路為：通過選擇合適的正交配點與節點，構造全局插值多項式來逼近狀態和控制變量，配點處狀態量的導數可由全局插值多項式求導來近似，從而將微分方程約束轉換為一組代數約束。性能指標中的積分項可由相應的Gauss求積公式計算得到。經上述離散化方法，原最優控制問題轉化為非線性規劃問題，而后通過有效的大規模稀疏NLP求解器即可求解。

Gauss偽譜法與向前Radau偽譜法基本原理如圖1所示，其主要區別在于：1)選用的配點不同。GPM選勒讓德-高斯(Legendre-Gauss, LG)點∈(-1,1)，而向前RPM選標準勒讓德-高斯-拉道(Legendre-Gauss-Radau, LGR)點∈[-1,1)；2)近似狀態量所用的方式不同。GPM用LG點與初始時刻點為節點構造插值多項式，終端狀態通過高斯求積公式得到，而向前RPM用全部離散點為節點構造插值多項式近似所有狀態量。關于Gauss偽譜法文獻較多，在此不再贅述。具體GPM離散化方法可參見文獻[4]。

圖1 Gauss偽譜法與向前Radau偽譜法基本原理圖

1.2 向前Radau偽譜法離散化最優控制問題[5]

1)時域變換

為將原問題由t∈[t0,tf]轉化到τ∈[-1,1]，對時間變量t作變換：

(1)

時域變換后的Bolza問題如下：

性能指標

J=Φ(x(-1),t0,x(1),tf)+

(2)

動力學微分方程約束

(τ),u(τ),τ;t0,tf)

(3)

邊界條件約束

φ(x(-1),t0,x(1),tf)=0

(4)

路徑約束

C(x(τ),u(τ),τ;t0,tf)≤0

(5)

2)配點選取及狀態與控制變量的近似

選N次與N-1次Legendre正交多項式之和PN(τ)+PN-1(τ)的根即標準LGR點作為配點，配點數為N。以τN+1=1和N個標準LGR點為節點構造Lagrange插值多項式來近似狀態變量：

(6)

Li(τ)(i=1,…,N+1)為Lagrange基函數。

以N個LGR點為節點構造Lagrange插值多項式來近似控制變量：

(7)

3)動力學微分方程轉化

對(6)式求導，結合(3)式可以得出：

,Uk,τk;t0,tf)=0

(k=1,…,N)

(8)

其中Xk=X(τk)∈Rn，Uk=U(τk)∈Rm為配點處狀態和控制量，Dki∈RN×(N+1)為微分矩陣，可通過式(9)離線確定：

(9)

其中φ(τi)=(1-τi)[PN(τi)+PN-1(τi)]。

4)性能指標

性能指標中的積分項均通過高斯拉道求積公式來近似：

J=Φ(X0,t0,Xf,tf)+

(10)

其中，wk為高斯拉道求積公式權系數。

5)邊界條件約束與路徑約束

邊界與路徑約束只需將離散狀態和控制量代入式(4)，(5)即可。

至此，由式(8)，(10)以及離散化后的式(4)，(5)將原最優控制問題轉化為非線性規劃問題。

2 SLV返回軌跡優化

2.1 動力學模型

考慮地球旋轉影響，假定飛行器側滑角為0，SLV無動力返回動力學方程為[6]：

γ

2ω(tanγcosφsinψ-sinφ)-

(11)

從實際系統考慮，控制舵面偏轉限制必然使得攻角α、傾側角σ存在帶寬和速率限制，在模型中考慮α，σ角速率限制。令

，

(12)

由式(11)與(12)聯立，可得最優控制問題狀態變量X=[rθφVγψασ]T∈R8，而控制變量為U=[uαuσ]T∈R2。uα,uσ稱為“偽控制量”，可以對其加以限制來滿足實際控制能力約束。

為加快優化速度、有效地利用自動微分技術，密度ρ(r)、當地音速a(r)以及升阻力系數CL,CD均采用擬合模型。

當李莉把交了兩個月租金的出租屋讓給許峰，自己投奔梅子的時候，梅子指著她的腦袋痛心疾首：“世界上最傻的女人就是你，許峰那個白眼狼，一看就是個攀高枝的。他騙你騙得高段啊，你走得那個順溜啊。”

2.2 優化性能指標與各類約束

1)性能指標與終端約束

為描述性能指標與終端約束，引入“末端進場走廊”(FAC)概念(見圖2)，FAC也可以理解為“進場著陸窗口”，是在返回段終端速度、航跡角與航向角、著陸場位置及方向給定情況下，能夠保證飛行器安全返回的終端位置(經緯高)的范圍。

圖2 末端進場走廊(FAC)[6]

終端約束要求在滿足其它各類約束條件下飛行器能夠準確的捕獲FAC，即

≤r(tf)≤

V(tf)=Vf

γ(tf)=γf

ψ(tf)=ψf

其中rFAC，θFAC，φFAC分別為FAC中心的地心距、經度、緯度，FACL，FACW，FACH分別為三維FAC的長、寬、高。

由于SLV返回的主要目標是能夠安全返回指定的著陸場，因此選擇如下性能指標：

w1,2,3為權重系數，這里選w1,2,3=1，即使得終端位置與末端進場走廊中心位置距離最小。

2)其它約束

3 優化算例與結果分析

3.1 優化算例

本文采用X-33總體參數及氣動模型進行仿真計算[6]。SLV質量35828kg，參考面積149.3881m2，地球半徑Re取6378.1km。

返回初始狀態與終端條件約束：

[r0,θ0,φ0,V0,γ0,ψ0]=[6429.1km,-85°,30°,2.6km/s,-1.3°,0°]

[Vf,γf,ψf]=[91.44m/s,-6°,-60°]

[rFAC,θFAC,φFAC]=[6378.7km,-80.7112°,30°]

[FACH/2,FACW/2,FACL/2]=[122m,0.00137°,0.001097°]

返回過程中狀態量約束：

≤≤

偽控制量約束：

≤≤

過載、動壓、熱流等路徑約束：

≤≤

3.2 優化結果

分別采用Gauss偽譜法和向前Radau偽譜法對不采用偽控制量和采用偽控制量2種情況進行了優化，采用國外大規模稀疏NLP軟件包SNOPT求解，優化結果如表1及圖3～7所示。兩種方法均選用30個配點，優化時間都在10s以內，考慮到本文算法在Matlab平臺下實現，雅克比矩陣通過自動微分方法計算，如采用C代碼、解析雅克比矩陣、無量綱化動力學方程等手段，優化速度有望進一步提高。優化性能指標表明終端時刻飛行器非常接近FAC中心。由圖4可以看出，采用偽控制量后GPM和向前RPM的攻角、傾側角曲線更為平滑。圖5表明采用GPM，控制量變化率均在要求的±5(°)/s范圍內，而采用向前RPM傾側角變化率在終端時刻達到了-5.21738(°)/s，超出了要求的范圍，這是由于標準LGR點不包括終端時刻，導致向前RPM在處理終端控制約束時受限。圖6表明返回過程滿足熱流、法向過載及動壓等路徑約束。兩種方法結果趨勢基本吻合，由于配點個數較少，也導致了一些差異。

表1 優化結果比較

圖3 狀態量變化曲線

圖4 控制量變化曲線

圖5 偽控制量變化曲線

圖6 路徑約束變化曲線

圖7 哈密爾敦函數變化曲線

3.3 可行性與最優性驗證

對于采用離散化方法求解連續最優控制問題而言，驗證解的可行性與最優性是非常必要的。

可行性通過數值積分結果與優化結果相比較來驗證，積分時的控制量由最優離散控制量插值得到。從圖3中可以看出GPM及向前RPM優化與數值積分結果非常吻合，表明所得結果是可行的。

與一般直接法不同，偽譜法可以根據協態映射原理估算配點處的協態信息，從而求得哈密爾敦函數值以檢驗解的最優性。對于存在路徑約束、終端時刻自由的自治系統最優控制問題，最優軌跡對應的哈密爾敦函數應恒為0。圖7中哈密爾敦函數接近于0，實際數據在10-11量級，可以認為所得結果接近最優解。

4 結論

本文利用Gauss偽譜法和向前Radau偽譜法進行了亞軌道飛行器返回軌跡優化，比較了2種方法在處理復雜問題時的能力。仿真結果表明，在滿足控制能力約束、路徑約束等條件下，GPM能夠快速準確地生成捕獲FAC中心的SLV返回軌跡。這2種方法計算效率相當，但向前RPM在處理存在控制量約束問題時效果不佳，不能滿足終端控制能力約束。

參 考 文 獻

[1] Martin J C and LAW G W.Suborbital Reusable Launch Vehicles and Applicable Market[R].The Aerospace Corporation Report, 2002.

[2] Fahroo F, Ross I M.Advances in Pseudospectral Methods for Optimal Control[C].AIAA Guidance, Navigation and Control Conference and Exhibit.Honolulu, Hawaii, 2008:1-23.AIAA 2008-7309.

[3] Huntington G T, Benson D, Rao A V.A Comparison of Accuracy and Computational Efficiency of Three Pseudospectral Methods[C].AIAA Guidance, Navigation, and Control Conference and Exhibit.Hilton Head, South Carolina, 2007:1-25.AIAA 2007-6405.

[4] Huntington G T.Advancement and Analysis of a Gauss Pseudospectral Transcription for Optimal Control Problems[D].Cambridge, MA: Massachusetts Institute of technology, 2007.

[5] Garg D, Patterson M A, Darby C L, et al.Direct Trajectory Optimization and Costate Estimation of General Optimal Control Problems Using a Radau Pseudospectral Method[C].AIAA Guidance, Navigation, and Control Conference.Chicago, Illinois, 2009:1-29.AIAA 2009-5989.

[6] Singh B, Bhattacharya R.Optimal Guidance of Hypersonic Vehicles Using B-Splines and Galerkin Projection[C].AIAA Guidance, Navigation, and Control Conference and Exhibit.Honolulu, Hawaii, 2008:1-11.AIAA 2008-7263.