劉宇航，楊洪偉，李 爽

(南京航空航天大學(xué)航天學(xué)院，南京 211106)

0 引 言

電推進因其高比沖的特點可以大幅降低燃料消耗，在深空探測任務(wù)中具有重要的作用。電推進發(fā)動機推力小，持續(xù)開機時間長，相應(yīng)的小推力軌跡優(yōu)化問題具有強非線性。間接法是求解小推力軌跡優(yōu)化問題的主要方法之一。在間接法中，因引入的協(xié)態(tài)變量沒有實際物理含義，且取值范圍不定，所以協(xié)態(tài)變量初值猜測問題是間接法主要困難之一。對于燃料最優(yōu)問題而言，最優(yōu)控制律為bang-bang控制，這種控制律的不連續(xù)性則進一步加劇了協(xié)態(tài)變量猜測的敏感性,提升了間接法打靶求解時的收斂困難。針對協(xié)態(tài)變量初值猜測困難問題，Bertrand等提出了平滑技術(shù)，即構(gòu)造優(yōu)化指標(biāo)同倫函數(shù)，先求解較易的能量最優(yōu)問題，然后通過同倫參數(shù)的改變，逐步同倫到較難求解的燃料最優(yōu)問題。Jiang等進一步提出了協(xié)態(tài)變量歸一化技術(shù)和開關(guān)函數(shù)檢測技術(shù)，發(fā)展了求解高效的小推力軌跡優(yōu)化同倫法。潘迅等以及沈紅新等利用同倫法分別解決了由中途飛越約束和多圈轉(zhuǎn)移引起的小推力軌跡優(yōu)化困難，豐富了同倫法的使用場景。為了避免協(xié)態(tài)初值猜測， Wu等、Yang及其團隊分別基于線性化、引力場同倫、推力同倫等方法推導(dǎo)了協(xié)態(tài)變量初值的近似解析解。此外，也有李鑒等利用無損卡爾曼濾波(UKF)將問題轉(zhuǎn)化，避免了初值猜測的困難。但上述文獻研究的小推力軌跡優(yōu)化方法針對的是定比沖情形。

對于深空探測任務(wù)而言，由于電推進發(fā)動機功率會隨著航天器相對太陽的距離增大而降低，調(diào)節(jié)比沖可以提供所需的小推力幅值，因此變比沖電推進更具有工程實用性。但是使用變比沖模型會使變量增加，使得軌跡優(yōu)化產(chǎn)生了新的困難。Taheri則將含雙曲正切函數(shù)的復(fù)合平滑控制方法應(yīng)用于變比沖發(fā)動機模型中。Chi等提出了新的同倫指標(biāo)函數(shù)，解決了傳統(tǒng)同倫法求解時存在的控制量耦合困難。Li等進一步研究了雙同倫方法。雖然同倫法求解變比沖小推力軌跡優(yōu)化問題效率很高，但仍需協(xié)態(tài)變量初值猜測，并且需要同倫迭代，無法直接求解燃料最優(yōu)問題，不適用于在線軌跡優(yōu)化。對于深空探測小推力轉(zhuǎn)移過程而言，由于存在各種擾動，探測器會偏離標(biāo)稱軌跡，存在軌跡在線重規(guī)劃的需求。本文將針對燃料最優(yōu)變比沖小推力軌跡優(yōu)化問題，研究無需猜測的協(xié)態(tài)變量初值快速確定技術(shù)，發(fā)展具備在線軌跡規(guī)劃能力的軌跡優(yōu)化方法。

近年來，人工智能在航天領(lǐng)域的應(yīng)用備受關(guān)注，為解決航天工程問題提供了新思路，在航天動力學(xué)與控制領(lǐng)域也已取得了較為豐碩成果。在深空探測軌跡優(yōu)化與設(shè)計方面，基于機器學(xué)習(xí)的方法可以用于估計轉(zhuǎn)移可達性快速評估、小推力轉(zhuǎn)移剩余質(zhì)量快速預(yù)測、最優(yōu)軌跡快速規(guī)劃等。針對深空探測任務(wù)中航天器自主性的需求迫切，基于人工智能技術(shù)的小推力最優(yōu)實時制導(dǎo)與控制方法也逐漸被提出。結(jié)合間接法理論和人工智能技術(shù)，發(fā)展小推力軌跡快速優(yōu)化或者實時制導(dǎo)方法是當(dāng)前研究的一個主要技術(shù)途徑。本文研究小推力軌跡快速優(yōu)化也將以間接法為基礎(chǔ)，圍繞協(xié)態(tài)變量初值敏感性這一制約軌跡優(yōu)化效率的因素，采用機器學(xué)習(xí)方法建立協(xié)態(tài)變量與狀態(tài)量的映射關(guān)系，實現(xiàn)協(xié)態(tài)變量初值的高精度高效估計。

基于機器學(xué)習(xí)估計協(xié)態(tài)變量初值時，首先需要生成大量最優(yōu)軌跡，建立數(shù)據(jù)集。Yin等提出通過擾動標(biāo)稱最優(yōu)軌跡參數(shù)產(chǎn)生數(shù)據(jù)集的方法。該方法將標(biāo)稱軌跡的協(xié)態(tài)變量初值代入擾動后的燃料最優(yōu)控制問題求解，因考慮的狀態(tài)量擾動小，故易收斂。但是，當(dāng)擾動量增大時，標(biāo)稱軌跡與受擾軌跡的協(xié)態(tài)變量初值偏差將增大，導(dǎo)致求解難以收斂。本文研究擬基于最優(yōu)軌跡延拓，發(fā)展適用于大擾動情形的數(shù)據(jù)集高效生成方法。同時，也將分析軌跡延拓策略與擾動上限和數(shù)據(jù)集生成效率的關(guān)系。其次，基于機器學(xué)習(xí)的方法需要設(shè)計合理的人工神經(jīng)網(wǎng)絡(luò)映射關(guān)系。現(xiàn)有的研究中映射輸入量往往僅考慮位置速度或者軌道根數(shù)的一到兩種狀態(tài)量形式，應(yīng)用于小推力軌跡優(yōu)化時收斂率還存在一定的不足。本文將研究多形式狀態(tài)量組合輸入的方法，進一步提升收斂率。此外，本文也將分析人工神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)對求解效率的影響并對其進行優(yōu)化設(shè)計。

1 變比沖小推力軌跡優(yōu)化問題

1.1 變比沖小推力轉(zhuǎn)移動力學(xué)模型

本文假設(shè)航天器采用離子電推進發(fā)動機，其能量全部來源于太陽能。在太陽引力的作用下，航天器的動力學(xué)模型為：

(1)

≤≤

(2)

推力幅值由下式計算得到：

(3)

式中:為發(fā)動機輸入功率利用效率；為最大輸入功率。

在實際航天器工作中，為保證系統(tǒng)功能正常運行，所以由太陽能電池板產(chǎn)生的電能首先供給除發(fā)動機外的設(shè)備，其次再供給發(fā)動機產(chǎn)生推力。而太陽能電池板的輸出功率受航天器與太陽距離的影響，具體關(guān)系如下:

(4)

式中:為太陽能電池板的輸出功率；為1個AU距離下太陽能電池板輸出功率；括號中的部分代表太陽能板效率隨航天器相對太陽的距離改變的經(jīng)驗值。進一步，可以列出發(fā)動機的輸入功率：

(5)

式中:為發(fā)動機輸入功率；為發(fā)動機功率處理器最大輸入功率；為占空比，代表發(fā)動機功率處理器真實工作時間與總開機時間的比值；為除發(fā)動機系統(tǒng)外其它系統(tǒng)所需總功率。

1.2 小推力軌跡燃料最優(yōu)控制

首先建立燃料最優(yōu)控制問題性能指標(biāo)：

(6)

在求解燃料最優(yōu)控制問題中,由于最后推導(dǎo)出的最優(yōu)控制律為bang-bang控制(只取0或1)，所以采用式(6)的性能指標(biāo)在開關(guān)切換點處對狀態(tài)和協(xié)態(tài)微分方程積分時不連續(xù)，且增加了對于協(xié)態(tài)變量初值猜測的難度。于是在構(gòu)造性能指標(biāo)時采用同倫方法，引入同倫參數(shù),由1逐漸變?yōu)?，=0時為燃料最優(yōu)控制問題。

由于協(xié)態(tài)變量沒有具體的物理含義，在選取協(xié)態(tài)變量初值時沒有明確的猜測范圍，極大地影響了求解效率，本文采用文獻[5]中提出協(xié)態(tài)變量歸一化的方法，引入歸一化協(xié)態(tài)變量，將包括在內(nèi)的8個協(xié)態(tài)變量初值限制到一個8維的單位球面，可極大地提高猜測效率。

由此燃料最優(yōu)控制問題性能指標(biāo)改寫為：

(7)

進一步，哈密頓函數(shù)為：

(8)

式中：，和分別表示位置、速度和質(zhì)量的協(xié)態(tài)變量。

根據(jù)龐德里亞金極小值原理，若使哈密頓函數(shù)取極小值，則·取極小值，此時與方向相反，最優(yōu)推力方向：

(9)

推導(dǎo)出歐拉-拉格朗日方程為：

(10)

(11)

(12)

由于比沖在區(qū)間主動變化，則最優(yōu)比沖大小計算方法為：

(13)

式中:如下：

(14)

對于采用改進對數(shù)同倫函數(shù)的式(7)，最優(yōu)發(fā)動機節(jié)流度寫為：

(15)

開關(guān)函數(shù)：

(16)

考慮航天器在任務(wù)中始末位置速度固定，任務(wù)時間固定。末態(tài)質(zhì)量自由，則對應(yīng)的末態(tài)質(zhì)量的協(xié)態(tài)變量()=0。打靶方程為：

(17)

式中:?[,,,]。

至此，將燃料最優(yōu)控制問題轉(zhuǎn)化為兩點邊值問題，滿足打靶方程的情況下，將狀態(tài)初值和協(xié)態(tài)變量初值代入運動微分方程(式(1))、協(xié)態(tài)微分方程(式(10)～(12))中即可得到燃料最優(yōu)控制問題的最優(yōu)解。所以，協(xié)態(tài)變量初值的猜測對于問題求解十分重要，本文介紹一種基于機器學(xué)習(xí)進行協(xié)態(tài)變量初值估計的方法。

1.3 求解方法簡述

本文采用間接法解決小推力軌跡優(yōu)化問題。間接法是將小推力問題作為最優(yōu)控制問題，首先給定一個性能指標(biāo),引入?yún)f(xié)態(tài)變量,構(gòu)造哈密頓函數(shù)，根據(jù)最優(yōu)性的一階必要條件推導(dǎo)出狀態(tài)變量和協(xié)態(tài)變量的微分方程和最優(yōu)控制律。此時將最優(yōu)控制問題轉(zhuǎn)換為兩點邊值問題，通過打靶的形式在滿足邊界條件下對狀態(tài)和協(xié)態(tài)微分方程進行積分，得出小推力控制問題的最優(yōu)解。

2 基于最優(yōu)軌跡延拓的數(shù)據(jù)集高效生成方法

在使用機器學(xué)習(xí)方法前，需生成大量燃料最優(yōu)軌跡的數(shù)據(jù)集，而大量求解燃料最優(yōu)控制問題會花費大量時間成本。文獻[21]中提出了一種最優(yōu)軌跡快速生成的方法，其原理是給標(biāo)稱軌跡末端狀態(tài)一個很小的擾動值，將標(biāo)稱軌跡的協(xié)態(tài)變量初值當(dāng)作猜測值代入擾動后的燃料最優(yōu)控制問題中，因擾動量很小，所以標(biāo)稱軌跡的協(xié)態(tài)變量值接近真實值，猜測很容易收斂。

本文中提出了一種基于最優(yōu)軌跡延拓的數(shù)據(jù)集高效生成方法，相比文獻[21]的優(yōu)勢在于：擾動量大的情形下也可高效地生成數(shù)據(jù)，有利于生成空間范圍更大的數(shù)據(jù)集。

以軌道根數(shù)作為狀態(tài)量的軌跡生成為例，基于最優(yōu)軌跡延拓生成數(shù)據(jù)集的步驟是：

1)確定新軌跡與標(biāo)稱最優(yōu)軌跡不同的狀態(tài)量，可以是末狀態(tài)不同，也可以是初狀態(tài)不同。以初狀態(tài)不同為例，新軌跡的初始狀態(tài)記為(+Δ,+Δ,+Δ,,,)，其中(Δ,Δ,Δ)為擾動量。

2)設(shè)定擾動上限(Δ,Δ,Δ)，使擾動量(Δ,Δ,Δ)在不超過擾動上限的范圍內(nèi)隨機選取。

3)將擾動量劃分為次達到，相鄰兩次擾動量之差為(Δ,Δ,Δ)。每一次增大擾動量時，都計算當(dāng)前擾動量下的最優(yōu)軌跡，逐漸增加擾動量，直至達到設(shè)定的擾動量(Δ,Δ,Δ)。此時，計算了條不同軌跡，且他們初始狀態(tài)之間的差值都是線性關(guān)系。

值得注意的是，在每一次迭代中協(xié)態(tài)變量的猜測值都采用上一次得到的協(xié)態(tài)變量初值，當(dāng)上一次迭代中沒有收斂，則使用上一次的猜測值。雖然總擾動量(Δ,Δ,Δ)很大，但是每兩代間的擾動量(Δ,Δ,Δ)較小，所以猜測很容易收斂?；谧顑?yōu)軌跡延拓數(shù)據(jù)集生成方法的算法邏輯如圖1。產(chǎn)生數(shù)據(jù)集的效率及數(shù)據(jù)集的空間范圍受到總擾動量(Δ,Δ,Δ)和延拓次數(shù)的影響，需要選擇合適值。

3 基于機器學(xué)習(xí)的協(xié)態(tài)變量初值估計

3.1 人工神經(jīng)網(wǎng)絡(luò)

人工神經(jīng)網(wǎng)絡(luò)(ANN)由輸入層、隱含層和輸出層構(gòu)成，其原理是通過許多人工神經(jīng)元組合構(gòu)造出輸入和輸出的映射關(guān)系，人工神經(jīng)網(wǎng)絡(luò)可以用于具有非線性的問題中，原理上來說可以構(gòu)造任意輸入輸出的映射關(guān)系。人工神經(jīng)元機理類似于生物神經(jīng)元，是將輸入加以權(quán)重和偏置，通過激活函數(shù)將此加入權(quán)重和偏置的輸入值映射為輸出值。訓(xùn)練人工神經(jīng)網(wǎng)絡(luò)的過程就是不斷迭代優(yōu)化人工神經(jīng)元中權(quán)重和偏置的值，使得指標(biāo)函數(shù)最小，如均方差(MSE)函數(shù)。本文中神經(jīng)網(wǎng)絡(luò)的訓(xùn)練采用反向傳播算法(BP算法)。

圖1 基于最優(yōu)軌跡延拓數(shù)據(jù)生成算法Fig.1 Data generation algorithm based on optimal trajectory continuation

3.2 超參數(shù)的選擇

超參數(shù)的選擇會直接影響訓(xùn)練效果，選擇合適的超參數(shù)是必要的。主要確定的超參數(shù)有：激活函數(shù)、神經(jīng)元數(shù)量、隱含層層數(shù)、學(xué)習(xí)率等。在本文中，激活函數(shù)使用雙曲正切S型函數(shù)，學(xué)習(xí)率設(shè)置為0.001，神經(jīng)元數(shù)量、隱含層層數(shù)的選取將在4.3節(jié)中給出。

4 神經(jīng)網(wǎng)絡(luò)訓(xùn)練與分析

4.1 仿真參數(shù)設(shè)定

使用Fortran語言編寫最優(yōu)軌跡生成程序并通過MATLAB神經(jīng)網(wǎng)絡(luò)工具箱進行訓(xùn)練，硬件使用Intel i7-8700K的CPU，主頻3.70 GHz。

本文以結(jié)束任務(wù)前500天的時刻進行重優(yōu)化為例，校驗最優(yōu)軌跡延拓數(shù)據(jù)生成法以及人工神經(jīng)網(wǎng)絡(luò)在估計協(xié)態(tài)變量初值的效果。

在本文中擾動軌道根數(shù)以產(chǎn)生數(shù)據(jù)集，動力學(xué)方程仍在笛卡爾坐標(biāo)系下建立。標(biāo)稱軌跡的初末狀態(tài)及協(xié)態(tài)變量初值見表1～2。

表1 初末狀態(tài)Table 1 Initial and final state values

表2 協(xié)態(tài)變量初值Table 2 Initial values of co-state variables

4.2 擾動上限對軌跡生成效率影響

本文中擾動后的軌跡將以與標(biāo)稱軌跡不同的初始狀態(tài)到達與標(biāo)稱軌跡相同的結(jié)束狀態(tài)。以計算10000次為例，比較延拓次數(shù)為1、10和100時所需時間和成功率，=1時對應(yīng)為直接擾動方法，即文獻[21]中所述方法。在表3中，擾動編號從A到F擾動上限依次增大，具體擾動大小從0至擾動上限值中隨機選取。

由圖2可以看出，隨著擾動上限的增大，三種方法生成最優(yōu)解的個數(shù)都會降低。在擾動上限大小相

表3 擾動上限設(shè)定Table 3 Set of different upper limits of the perturbation

圖2 生成軌跡個數(shù)隨擾動上限變化曲線Fig.2 The number of generated trajectories with respect to the upper limit of the perturbation

圖3 計算時間隨擾動上限變化曲線Fig.3 Calculation time varies with the upper limit of the perturbation

同的情況下，生成最優(yōu)解個數(shù)由高到低依次為：延拓100次、延拓10次、直接擾動。由圖3可以看出,延拓100次的方法比直接擾動方法計算10000次所需時間更短。當(dāng)擾動上限為F時，直接擾動用時大量減少，其原因是成功生成解的個數(shù)特別少，所以少了很多積分求解軌跡的時間，總時間相較E時下降較多。

綜上，可以看出基于最優(yōu)軌跡延拓相比于直接擾動成功率更高，計算速度更快，尤其是在擾動量較大時更為明顯。而且，延拓次數(shù)的增大會提高生成最優(yōu)軌跡的效率。經(jīng)計算，在擾動編號D條件下，延拓100次比直接擾動成功生成最優(yōu)解速度快4.9倍；在擾動編號E條件下，快7.7倍；在擾動編號F條件下，快15.8倍。

4.3 神經(jīng)網(wǎng)絡(luò)訓(xùn)練

本文使用50萬的數(shù)據(jù)集進行訓(xùn)練，另取2萬組數(shù)據(jù)作為樣本用于誤差分析?？紤]到設(shè)定的擾動上限盡可能與實際可達到擾動上限接近，并盡量使擾動上限取較大值。最終，數(shù)據(jù)集生成時設(shè)定擾動上限為擾動編號E對應(yīng)的值，并使用延拓次數(shù)為100次的基于最優(yōu)軌跡延拓方法生成數(shù)據(jù)集。最后求得50萬組數(shù)據(jù)的實際擾動上限見表4。

表4 實際擾動上限Table 4 Actual upper limit of perturbation

為對比不同輸入對于協(xié)態(tài)變量初值猜測效果的影響，本文分別將速度位置、軌道根數(shù)、改進春分點軌道根數(shù)作為ANN的輸入；將位置協(xié)態(tài)初值、速度協(xié)態(tài)初值、質(zhì)量協(xié)態(tài)初值、歸一化協(xié)態(tài)初值作為ANN的輸出，共8維輸出。三種ANN不同輸入的方法中，只將由半長軸、偏心率和軌道傾角引起改變的量作為ANN的輸入。在進行神經(jīng)網(wǎng)絡(luò)訓(xùn)練前，對所有的輸入量進行歸一化處理。

1)方法一：速度位置作輸入

采用速度位置作為輸入時，速度和位置各有三個維度，共有六個量作為ANN的輸入。

2)方法二：軌道根數(shù)作輸入

采用軌道根數(shù)作為ANN的輸入時，輸入為半長軸、偏心率、軌道傾角三個量。

3)方法三：改進春分點軌道根數(shù)作輸入

采用春分點軌道根數(shù)作為ANN的輸入時，輸入為除外的其余五個改進春分點軌道根數(shù)元素。其中，=++，，和分別表示航天器的真近點角、近地點幅角和升交點赤經(jīng)。

式(18)和式(19)分別為平均相對誤差和標(biāo)準(zhǔn)偏差的計算公式。表5～10中給出了三種輸入下不同神經(jīng)元數(shù)量和隱含層層數(shù)不同時平均相對誤差(%)和標(biāo)準(zhǔn)偏差(×10)大小。

(18)

(19)

表5 平均相對誤差(方法一)Table 5 Average relative error (Method 1)

表6 標(biāo)準(zhǔn)偏差(方法一)Table 6 Standard deviation (Method 1)

表7 平均相對誤差(方法二)Table 7 Average relative error (Method 2)

表8 標(biāo)準(zhǔn)偏差(方法二)Table 8 Standard deviation (Method 2)

表9 平均相對誤差(方法三)Table 9 Average relative error (Method 3)

表10 標(biāo)準(zhǔn)偏差(方法三)Table 10 Standard deviation (Method 3)

式中:為樣本數(shù)；為神經(jīng)網(wǎng)絡(luò)預(yù)測的協(xié)態(tài)變量初值；為真實的協(xié)態(tài)變量初值。

在列舉的五種神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)中選取平均相對誤差和標(biāo)準(zhǔn)偏差總體較小的網(wǎng)絡(luò)結(jié)構(gòu)。最終確定：當(dāng)輸入為位置速度時，采用4層隱含層，每層15個神經(jīng)元的結(jié)構(gòu)；當(dāng)輸入為軌道根數(shù)和改進春分點軌道根數(shù)時，采用3層隱含層，每層30個神經(jīng)元的結(jié)構(gòu)。

4.4 協(xié)態(tài)估計效果分析

在50萬的訓(xùn)練數(shù)據(jù)外選取2萬數(shù)據(jù)，用于測試人工神經(jīng)網(wǎng)絡(luò)初值估計的效果。

表11中給出了三種不同輸入形式訓(xùn)練的ANN用于初值估計時最優(yōu)控制問題的求解收斂率。神經(jīng)網(wǎng)絡(luò)輸入為改進春分點軌道根數(shù)時收斂率最高，但是三種ANN模型求解收斂率相差不大，均在95.00%左右。經(jīng)計算，同倫法隨機猜測一次的計算收斂率為70.00%，遠低于本文中方法。為了進一步提高收斂率，將三種不同輸入形式的ANN結(jié)合，選定其中一個神經(jīng)網(wǎng)絡(luò)作為主輸出，當(dāng)此神經(jīng)網(wǎng)絡(luò)所估計值求解不收斂時，使用其他兩種輸入的ANN重新估計，這樣可以將求解收斂率提高到99.78%，算法流程見圖4。圖5為對于2萬測試數(shù)據(jù)，三種不同輸入形式的ANN結(jié)合所成功求解出的軌跡。

對2萬個樣本中成功收斂部分的計算時間進行統(tǒng)計。采用ANN估計的協(xié)態(tài)初值進行軌跡優(yōu)化的時間消耗分為兩部分：神經(jīng)網(wǎng)絡(luò)協(xié)態(tài)估計用時和軌跡優(yōu)化用時。經(jīng)計算，采用方法一至方法三的平均時間分別為0.02794 s、0.02722 s和0.02678 s；使用傳統(tǒng)同倫方法的求解平均用時約為0.38770 s，本文方法比傳統(tǒng)同倫方法約快13.88倍。

表11 ANN協(xié)態(tài)初值估計求解的收斂率Table 11 The convergence rate of estimating the initial values of co-state variables using ANN

圖4 基于組合輸入的ANN求解最優(yōu)軌跡算法Fig.4 ANN algorithm for solving the optimal trajectory based on the combined inputs

圖5 求解成功軌跡圖Fig.5 Successfully solved trajectories

5 結(jié) 論

本文針對變比沖小推力軌跡間接優(yōu)化中協(xié)態(tài)變量初值猜測困難的問題，提出一種基于機器學(xué)習(xí)進行協(xié)態(tài)初值高效高精度估計的方法。在最優(yōu)軌跡生成中，隨著擾動上限增大，求解效率會明顯下降。為了提高求解效率，提出了基于最優(yōu)軌跡延拓的數(shù)據(jù)集生成方法。本文提出的數(shù)據(jù)集生成方法在高擾動上限情況下可以將數(shù)據(jù)集生成速度提升幾倍甚至十幾倍。并構(gòu)建了基于位置速度、軌道根數(shù)和改進春分點軌道根數(shù)多形式狀態(tài)量組合輸入的人工神經(jīng)網(wǎng)絡(luò)。多形式狀態(tài)量組合輸入可以將求解收斂率提高到99.78%。此外，軌跡求解的平均總用時僅為0.02700 s左右。仿真結(jié)果表明，采用本文方法可以高精度、高效地估計協(xié)態(tài)變量初值，本文方法可用于燃料最優(yōu)變比沖小推力在線軌跡優(yōu)化設(shè)計。