邱 豐，宋征宇

(北京航天自動控制研究所，北京100854)

采用聯(lián)立法求解大姿態(tài)終端約束的上升段軌跡優(yōu)化

邱 豐，宋征宇

(北京航天自動控制研究所，北京100854)

針對恒定推力和未配置輔助動力系統(tǒng)的運載火箭，提出了一種基于聯(lián)立框架的直接法(“聯(lián)立法”)來求解帶大姿態(tài)終端約束的動態(tài)在線軌跡規(guī)劃問題，解決了傳統(tǒng)迭代制導方法中飛行軌跡線性規(guī)劃的不足，并可根據(jù)姿態(tài)跟蹤精度的需求來控制程序角速率，從而同時滿足各種約束條件。該聯(lián)立法結(jié)合了直接配點法和偽譜法的特點，采用有限元正交配置法將狀態(tài)變量和控制變量進行離散化;結(jié)果表明與偽譜法相比，計算效率得到了大幅提升。在此基礎(chǔ)上，提出了通過合理選擇初值、調(diào)整入軌精度偏差、自適應(yīng)有限元分段等提升計算速度的措施，進一步提高了這一算法的實時性，為未來解決更復雜多樣的制導任務(wù)提供了基礎(chǔ)。

運載火箭;姿態(tài)終端約束;聯(lián)立法;軌跡規(guī)劃

0 引 言

近年來，自主控制的需求在航天應(yīng)用中愈發(fā)增強。上升段飛行中，對于目標變更和推力下降等異常情況，自主規(guī)劃的適應(yīng)性較跟蹤標準軌跡的攝動制導要高很多。自載人航天交會對接任務(wù)起，我國運載火箭首次采用了迭代制導方法［1］。這是一種在線動態(tài)自主軌跡規(guī)劃的方法，但由于飛行軌跡的控制是通過調(diào)整姿態(tài)角來實現(xiàn)的，使得入軌姿態(tài)的散布很大。

解決入軌姿態(tài)可以從總體方案著手，例如在火箭的末級增加輔助動力系統(tǒng)(也稱調(diào)姿系統(tǒng))，或者依靠有效載荷自身進行姿態(tài)調(diào)整，但這些方案均額外增加了系統(tǒng)資源和復雜性。為此，文獻［2］提出了一種預(yù)測修正的迭代制導改進方法。但是，對于一些較為特殊的入軌姿態(tài)角，現(xiàn)有迭代制導中線性規(guī)劃的前提條件不再存在，導致入軌誤差很大。對于此類問題，即使傳統(tǒng)攝動制導也難以滿足精度要求(這將在后文進行分析)。

本文提出了一種基于聯(lián)立框架的直接法(下文簡稱“聯(lián)立法”)，將入軌時間最短作為性能指標，通過控制推力方向變化使得火箭以要求姿態(tài)進入預(yù)定軌道。該方法采用統(tǒng)一化的處理流程，非常適合計算機自動處理，在自主制導領(lǐng)域有較好的應(yīng)用前景。

1 帶終端姿態(tài)約束的最優(yōu)入軌問題

對于推力大小不可調(diào)節(jié)的運載火箭上升段實時軌跡優(yōu)化問題，其本質(zhì)上是帶有狀態(tài)約束和控制約束的最優(yōu)控制問題，間接法和直接法是最常用的兩種方法，但間接法對于末端有姿態(tài)角約束且入軌時間浮動的情況無法完成求解［3］。而直接法將連續(xù)空間最優(yōu)控制問題轉(zhuǎn)化為非線性規(guī)劃(NLP)問題，利用數(shù)值方法求解。文獻［4］利用偽譜法對固體運載火箭上升段軌跡優(yōu)化進行研究，但計算耗時無法滿足實時控制要求;文獻［5］改善了計算效率，但收斂性還需進一步驗證。目前研究成果表明［6－8］，直接法的實時性和收斂性問題并沒有完全解決。

本文研究更為一般性的帶終端姿態(tài)約束的精確入軌問題，采用有限元正交配置法將狀態(tài)變量和控制變量進行離散化，實現(xiàn)對原命題的分段低階逼近處理，使原問題轉(zhuǎn)化為NLP問題，同時提出了提高計算效率的措施。

1.1 火箭質(zhì)心運動模型

在不考慮氣動力影響的條件下，建立其相對于發(fā)射慣性坐標系(以下簡稱發(fā)慣系)的質(zhì)心運動方程如下［9］:

式中:x，y，z為火箭當前點到發(fā)射點矢徑在發(fā)慣系上的投影，Vx，Vy，Vz為火箭速度矢量在發(fā)慣系上的投影，gx，gy，gz為引力加速度在發(fā)慣系上的投影，m0為飛行器初始質(zhì)量，m·為秒耗量，P為推力，φ為俯仰角，ψ為偏航角。對于定推力發(fā)動機，秒耗量及推力設(shè)為定值。

由于狀態(tài)變量之間數(shù)量級相差較大，為了便于數(shù)值求解，需要對運動方程進行單位歸一化處理:

1.2 需要滿足的約束條件

有效載荷的運行軌道由五個軌道根數(shù)決定(如果考慮其某一時刻在軌道上的具體位置，則需要六個軌道根數(shù))［10］。只有所處位置在軌道上的一點，并且速度大小及方向恰好是軌道上該點所需速度的情況下，才能實現(xiàn)對五個軌道根數(shù)的精確控制。

帶姿態(tài)約束的最優(yōu)入軌問題中，不僅要控制火箭精確進入預(yù)定軌道，而且需滿足一定的入軌姿態(tài)。此外，由于火箭燃料有限，為了確保得到的解是可行的，還需要加入最大入軌時間的約束。

綜上所述，帶姿態(tài)約束的火箭最優(yōu)入軌問題中要滿足的約束如表1所示。

表1 命題約束列表Table 1 Constraints list of the proposition

2 基于有限元正交配置的聯(lián)立法

2.1 有限元正交配置法

假定整個問題的時間區(qū)間為［t0，tf］，將上升段入軌問題描述為一般動態(tài)優(yōu)化命題形式［11］:

式中:Φ(x(tf))是目標函數(shù)，x(t)∈Rnx為狀態(tài)變量組成的nx維向量;u(t)∈Rnu為控制變量組成的nu維向量;p∈Rnp是模型參數(shù);f為微分方程組的右邊約束;gE為代數(shù)方程組;微分代數(shù)方程組中的狀態(tài)變量和控制變量的邊界約束均歸結(jié)到gI中; x0∈Rnx為狀態(tài)變量的初始狀態(tài)，N1(x(tf))為終端目標集約束。

傳統(tǒng)偽譜法采用全局插值多項式擬合，為了保證逼近精度要求多項式具有較高階次，離散化后NLP問題規(guī)模很大，導致了計算難度和求解時間的增加。為了在保證計算精度的同時降低離散化后命題的復雜程度，可將整個時域劃分為N個區(qū)間，每個區(qū)間稱為有限元，在每段有限元［ti－1，ti］上，采用正交多項式對狀態(tài)變量x及控制變量u進行逼近，采用拉格朗日插值多項式可使變量在配置點上的值正好等于其系數(shù)。這種方法稱之為有限元正交配置法［12］。

此外，在有限元端點上需通過連接方程來保證狀態(tài)變量的連續(xù)性;同時應(yīng)滿足初始和終端約束條件。結(jié)合以上離散化后的條件，可以得到有限元正交配置法處理后的NLP問題。

2.2 聯(lián)立法對問題的重構(gòu)

以本問題為例，將整個時域［t0，tf］分為N個有限元，每個有限元內(nèi)用K階拉格朗日插值多項式對變量和約束進行逼近。采用Gauss-Radau點［13］作為配置點進行配置。

目標函數(shù)為入軌時間:

狀態(tài)變量為火箭的速度、位置:

控制變量為火箭的姿態(tài)程序角:

在推力大小不變的情況下，定點入軌是無法保證的。因此在求解過程中將真近點角f作為代數(shù)變量進行處理，通過五個軌道根數(shù)及f，可以確定入軌點速度、位置，從而形成六個等式約束方程。

(1)過程約束

火箭質(zhì)心運動的狀態(tài)變量微分方程表示如下:

式中，i∈［0，N)，j∈［0，K)。

狀態(tài)變量的連續(xù)性方程表示為(以x為例):

實際飛行中，還應(yīng)加入姿態(tài)角及姿態(tài)角速率的不等式約束(以俯仰角為例，ω為角速率上限):

(2)終端約束

狀態(tài)變量的初始約束和終端約束(以x為例):

目標入軌時間最大值的約束:

綜合式(7)～(11)可得帶終端姿態(tài)約束最優(yōu)入軌問題離散化后的大規(guī)模非線性規(guī)劃命題。

由式(9)和(11)可以看出該非線性規(guī)劃命題中含有不等式約束，直接對其求解非常困難，且約束個數(shù)隨著有限元個數(shù)N和插值階數(shù)K增加而成倍增加，故采用內(nèi)點算法［12］對命題進行求解。

至此，聯(lián)立法的計算過程可歸納為圖1所示，分為設(shè)計輸入、離散處理和實時規(guī)劃三部分。通過利用基于內(nèi)點算法的NLP求解器迭代求解，將求得變量代入離散化后的命題得到約束殘差再返回求解器。通過迭代計算使殘差漸漸逼近零，求得的變量也越來越逼近最優(yōu)解。在該計算過程中，終端姿態(tài)角約束被統(tǒng)一處理，并不區(qū)分姿態(tài)角的大小，體現(xiàn)了該方法的通用性。

2.3 聯(lián)立法的性能分析

在［t0，tf］內(nèi)將整個時域分成N段，聯(lián)立法通過有限元將狀態(tài)變量x(t)分別定義到多個時段內(nèi)。在有限元［ti－1，ti］上，由式(7)得到的近似最優(yōu)解xk+1(t)和真實解x(t)之間的誤差邊界可以通過下式來描述［13］:

式中:C為只與階次K有關(guān)的常數(shù)，hi為第i個有限元的長度，為真實解K+1次導數(shù)在區(qū)間［ti－1，ti］上的無窮范數(shù)。因此，當有界時，可以得到:

也就是說當有限元的劃分無限加密的時候，聯(lián)立法得到的近似解無限逼近真實的解，故可以通過控制離散化精度來間接保證解的最優(yōu)性。但工程實際應(yīng)用中不能為了少量精度的提高而增加大量的有限元，在線規(guī)劃時可根據(jù)實際用例選取合適的有限元，以達到求解精度和運算時間的平衡。

3 仿真校驗

3.1 仿真場景

表2～表3以某火箭上升段飛行為例列出了仿真場景的各項參數(shù)。暫不考慮大氣影響，即從火箭二級飛行段開始。表4為終端姿態(tài)要求，其中終端俯仰角為－110°，正是這一特殊約束，導致傳統(tǒng)算法無法同時兼顧入軌精度和終端姿態(tài)的要求。

表2 初始參數(shù)Table 2 Initial parameters

表3 目標軌道參數(shù)Table 3 Target orbit parameters

表4 姿態(tài)角約束Table 4 Attitude constraints

平均取30段有限元和3階Radau配點進行離散化，初值取法如下:

1)各配點速度、位置初值均設(shè)為初始點值;

2)真近點角變量的初值取為0°;

3)目標函數(shù)(總飛行時間)的初值取標準彈道入軌時間376.4 s;

4)由于俯仰角變化大，其初值采用等差選取的方法，而偏航角變化不大，初值取初始狀態(tài)的值。

3.2 仿真結(jié)果與分析

本次數(shù)值實驗是基于 IntelRCoreTMi5－3470 CPU處理器、3.2 GHz主頻的計算機進行的，求解器為IPOPT(3.6.1版)［14］。由于決定運算速度的主要因素是處理器主頻，所以不同硬件平臺下的計算耗時可以根據(jù)CPU的主頻按比例近似折算。

3.2.1 與標準彈道進行對比分析

圖2給出了標準彈道設(shè)計與聯(lián)立法計算結(jié)果的飛行程序角對比:對于終端俯仰角－110°的約束條件，相當于火箭向下飛行，這對運載能力是有損失的，因此標準彈道設(shè)計時將姿態(tài)調(diào)整段壓縮在飛行的末斷，時間要盡可能短，俯仰角速率最大達到9.35°/s，如圖3所示。六自由度仿真結(jié)果表明，在該角速率條件下，姿態(tài)跟蹤誤差達到了2°，這反過來影響了入軌精度。這就是引言中所提采用攝動制導跟蹤標準軌跡飛行也難以滿足入軌精度的原因。

為此在聯(lián)立法計算中將角速率約束定為4°/s，終端姿態(tài)和入軌精度均能滿足要求，入軌時間為377.5 s，比標準彈道多 1.1 s(多消耗推進劑約118.6 kg，在安全余量范圍內(nèi))，如圖3所示。且軌跡可以實時動態(tài)調(diào)整，及時修正之前的誤差，即誤差不具有積累的效果。

從圖4速度和位置曲線可以看出，該飛行軌跡盡管前后段均近似線性，但整體難以用線性規(guī)劃來統(tǒng)一描述;而大角速率轉(zhuǎn)彎飛行的起點是浮動的，無法提前獲知，使得分段線性化處理也不可能，這是導致現(xiàn)有迭代制導無法滿足該算例入軌精度要求的主要原因。

3.2.2 與偽譜法求解結(jié)果進行對比分析

在同一運算平臺下，分別用全局30和300個Radau配置點的偽譜法對命題進行求解，得到的結(jié)果與30段3階Radau配置點的聯(lián)立法求解情況對比，如表5和圖5所示。

表5 聯(lián)立法與偽譜法結(jié)果對比Table 5 Comparison of the simultaneous method and the pseudo-spectral method

由圖5可見，求解時間近似的情況下，聯(lián)立法求解精度高于使用全局近似的偽譜法;而求解精度相似時，偽譜法的耗時遠遠超過聯(lián)立法。因此聯(lián)立法在求解效率上優(yōu)于全局插值近似的偽譜法。

3.3 提升算法實時性的幾種手段

3.3.1 選取合理的初值

由于每個控制周期均迭代計算，即使在一個周期內(nèi)受到干擾等影響，實際飛行軌跡也不會偏離過大。因此可采用當前點速度位置作為起始點初值，上一次規(guī)劃的相關(guān)參數(shù)作為其他狀態(tài)及控制變量的初值，以縮短優(yōu)化時間，這對實時飛行的迭代計算意義十分重大。

針對本算例，將所有變量的初值均設(shè)置為已有優(yōu)化結(jié)果得到的最優(yōu)值附近，假設(shè)偏差為3%，則各個變量的初值取為:

計算結(jié)果對比如表6所示。

表6 不同初值結(jié)果對比Table 6 Comparison of different initial value

由表6可以看出，該方法可以大大縮短運算時間;且兩次計算之間的時間間隔越短，初值與最優(yōu)解之間的偏差也越小。

3.3.2 選擇合適的入軌點精度偏差

聯(lián)立法是將軌道根數(shù)轉(zhuǎn)化得到的速度、位置作為終端等式約束條件，即要求入軌點的速度、位置與理論值的差值為0。而實際應(yīng)用中，不同任務(wù)對不同軌道根數(shù)的精度要求也不同;并且在合理的偏差范圍內(nèi)，即差值不為0也是可以接受的。因此可以對不同軌道根數(shù)區(qū)別處理以提升計算效率。

以本算例為例，在保持整體容差設(shè)置為1×10－8不變的情況下，將變量值與要求值的差除以10000，則相當于把六個約束放寬到1×10－4，結(jié)果對比如表7和表8所示。

表7 放寬誤差約束后計算效率對比Table 7 Comparison of computational efficiency after the relaxation of error constraints

從表7中看出，放寬誤差約束可以縮短優(yōu)化時間，但對入軌精度的影響在可接受的范圍之內(nèi)，如表8所示。因此，未來可以根據(jù)不同算例的特點對約束進行區(qū)別調(diào)整，在不影響精度的情況下進一步提高實時性。

表8 放寬誤差約束后入軌精度對比Table 8 Comparison of orbit injection accuracy after the relaxation of error constraints

3.3.3 自適應(yīng)有限元分段

減少有限元的分段能有效減小離散后NLP命題的規(guī)模，但會降低解的精度和最優(yōu)性。通過上述飛行曲線的觀察可以發(fā)現(xiàn)，控制變量的大幅變化基本集中在初始和終端位置，中間很長的時間變化平緩。根據(jù)這一特點，在變量快速變化的區(qū)間內(nèi)加密有限元劃分，而在變化平緩的區(qū)域減少有限元分段，這種方法可稱作自適應(yīng)有限元分段法。

為此，將整個飛行過程分為三個時間段:［0，20 s］、［20 s，350 s］和［350，tf］，三段內(nèi)各取5個有限元，每個有限元內(nèi)同樣取3階Radau點進行配置，結(jié)果如表9所示。

表9 不同有限元分段結(jié)果對比Table 9 Comparison of different finite element configurations

在實際使用中，可以利用起飛前離線規(guī)劃的軌跡，通過姿態(tài)角變化率的求解，調(diào)整有限元分段;在飛行過程中，也可以利用這一規(guī)律自適應(yīng)地動態(tài)調(diào)整分段。

4 結(jié) 論

本文針對傳統(tǒng)迭代制導方法無法保證終端大姿態(tài)約束的問題，利用有限元正交配置進行低階逼近和聯(lián)立求解，與偽譜法等直接法相比計算效率得到大幅提高。并以某火箭上升段最優(yōu)入軌任務(wù)為例，對算法的效果進行了仿真驗證。結(jié)果表明，入軌精度和終端大姿態(tài)約束均能滿足，求解時間也縮短到百毫秒級，為實時在線規(guī)劃打下了基礎(chǔ)。

［1］ 宋征宇.從準確、精確到精益求精——載人航天推動運載火箭制導方法的發(fā)展［J］.航天控制，2013，31(1):4－10.［Song Zheng-yu.From accurate，precise to perfect-manned space promotes the development of guidance method on launch vehicle［J］.Aerospace Control，2013，31(1):4－10.］

［2］ 呂新廣，宋征宇，鞏慶海.運載火箭軌跡預(yù)測制導方法研究［J］.載人航天，2013，19(6):24－28.［Lv Xin-guang;Song Zheng-yu;Gong Qing-hai.Research on trajectory prediction guidance for launch vehicle［J］.Manned Spaceflight，2013，19 (6):24－28.］

［3］ 宋征宇.信息集成在航天運輸控制系統(tǒng)中的應(yīng)用研究［J］.宇航學報，2015，36(4)，365－374.［Song Zheng-yu.Application study ofinformation integration in aerospace transportation control system［J］.Journal of Astronautics，2015，36(4):365－374.］

［4］ 楊希祥，張為華.基于Gauss偽譜法的固體運載火箭上升段軌跡快速優(yōu)化研究［J］.宇航學報，2011，32(1)，15－21.［Yang Xi-xiang，Zhang Wei-hua.Rapid optimization of ascent trajectory for solid launch vehicles based on gauss pseudo-spectral method［J］.Journal of Astronautics，2011，32(1):15－21.］

［5］ 閆循良，廖守億，張金生，等.基于節(jié)點改善策略的偽譜軌跡優(yōu)化［J］.宇航學報，2013，34(7):891－900.［Yan Xun-liang，Liao Shou-yi，Zhang Jin-sheng，et al.Trajectory optimization using pseudospectral method based on a grid node refinement strategy［J］.Journal of Astronautics，2013，34(7):891－900］

［6］ Calise A J，Brqandt N.Generation of launch vehicle abort trajectories using a hybrid optimization method［J］.Journal of Guidance，Control，and Dynamics［J］.2004，27(6)，930－937.

［7］ Mease K D，Chen D T，Teufel P，et al.Reduced-order entry trajectory planning for acceleration guidance［J］.Journal of Guidance，Control，and Dynamics，2002，25(2):257－266.

［8］ 雍恩米，陳磊，唐國金，等.飛行器軌跡優(yōu)化數(shù)值方法綜述［J］.宇航學報，2008，29(2):397－406.［Yong En-mi，Chen Lei，Tang Guo-jin，et al.A survey of numerical methods for trajectory optimization of spacecraft［J］.Journal of Astronautics，2008，29(2):397－406.］

［9］ 徐延萬.彈道導彈、運載火箭控制系統(tǒng)設(shè)計與分析［M］.北京:宇航出版社，1999:41－53.

［10］ 張蒿.固體運載火箭彈道設(shè)計與優(yōu)化［J］.固體火箭技術(shù)，1997，20(1):1－5.［Zhang Song.Calculation and optimization of solid launch vehicle trajectory［J］.Journal of Solid Rocket Technology，1997，20(1):1－5.］

［11］ Betts J T.Survey of numerical methods for trajectory optimization［J］.Journal of Guidance，Control and Dynamics，1998，21 (2):193－207.

［12］ 邱豐，宋征宇.基于聯(lián)立框架的軌跡優(yōu)化算法［J］.航天控制，2015，33(2):9－15.［Qiu Feng，Song Zheng-yu.A trajectory optimization algorithm based on simultaneous strategies［J］.Aerospace Control，2015，33(2):9－15.］

［13］ 陳偉鋒.大規(guī)模復雜過程系統(tǒng)的高性能優(yōu)化理論與方法研究［D］.杭州:浙江大學，2011.［Chen Wei-feng.Research on high-performance optimization for large-scale complex process systems［D］.Hangzhou:Zhejiang University，2011.］

［14］ Wachter A，Biegler L T.On the implementation of an interiorpoint filterline search algorithm forlarge-scale nonlinear programming［J］.Mathematical Programming，2006，106(1):25－57.

通信地址:北京142信箱402分箱(100854)

電話:(010)88525917

E-mail:qfeng1991@126.com

宋征宇(1970－)，男，碩士，研究員，主要從事飛行器控制、制導與仿真，系統(tǒng)集成、智能自主控制技術(shù)等研究。本文通信作者。

通信地址:北京142信箱402分箱(100854)

電話:(010)68389311

E-mail:zycalt12@sina.com

(編輯:張宇平)

Large Terminal Attitude Constrained Trajectory Optimization of Ascent Stage Via Simultaneous Method

QIU Feng，SONG Zheng-yu
(Beijing Aerospace Automatic Control Institute，Beijing 100854，China)

For launch vehicles with constant thrust and without auxiliary propulsion system，a direct method on simultaneous framework(simultaneous method)is presented to solve the dynamic online trajectory planning problem with large terminal attitude constraint，making up for the deficiency of linear programming in the traditional iterative guidance method，meanwhile，the angular rate of the trajectory can be limited according to the requirement of the attitude tracking error，so as to satisfy all constraint conditions.By combing the characteristics of the direct collocation method and pseudo-spectral method，and based on orthogonal collocation on finite element，this paper adopts the simultaneous method and discretized state variables and control variables to solve this problem，which is much more efficient than the pseudo-spectral method.On this basis，the paper puts forward some measures to increase the calculation efficiency，such as initial value selection，orbit precision error adjustment，and adaptive finite element configuration，further improving the real-time performance and providing feasible solutions to more complex and diverse guidance tasks in the future.

Launch Vehicles;Terminal attitude constraint;Simultaneous method;Trajectory design

V4

A

1000-1328(2017)01-0018-08

10.3873/j.issn.1000-1328.2017.01.003

邱 豐(1991－)，女，碩士，主要從事導航、制導與控制。

2016-05-23;

2016-10-16