1、1航天器軌跡規(guī)劃的變分技術(shù)摘要：本文描述了一種基于變分法的軌跡規(guī)劃算法，它可以解決六自由度的航天器對接和接近操作問題。目標(biāo)函 數(shù)綜合考慮燃料消耗、避障距離和到達(dá)時間。以非線性軌道動力學(xué)方程作為動力學(xué)約束條件，推導(dǎo)了歐拉-拉格 朗日方程，依據(jù)龐德里亞金原理推導(dǎo)了最優(yōu)控制的輸入，給出了推力器實際的開關(guān)曲線。采用對初始條件不敏感 的間接配點法解決相應(yīng)的邊值問題，解的開拓可以進一步改進算法的魯棒性。文章討論了對歐拉-拉格朗日方程 和橫截條件的處理，以使其適合目前的配點法的應(yīng)用形式。文章最后通過一個與翻滾衛(wèi)星進行首尾相連的對接機 動對所研究的方法和結(jié)果進行了驗證。關(guān)鍵詞：航天器，軌道，空間探測，算法，航

2、天器任務(wù)規(guī)劃1.1引言一個成功的自主航天器必須能夠?qū)崿F(xiàn)智能導(dǎo)航。智能導(dǎo)航能夠找到并沿著無碰撞的路徑以合 理的時間和燃料穿過預(yù)定空間。很多自主航天任務(wù)目前已經(jīng)發(fā)射，還有一些其他任務(wù)，包括 DAPPA（美國國防高級研究計劃局）的SUMO/FREND服務(wù)航天器，AFRL （美國空軍研究實驗室） 和Lockheed-Martin公司的ANGELS觀察衛(wèi)星，目前都在設(shè)計之中。尤其是，DARPA的軌道快 車任務(wù)，該任務(wù)在技術(shù)方面向我們展示了自主交會、抓捕、燃料加注和重構(gòu)等技術(shù)，并在2007 年7月成功完成了其主要任務(wù)。軌道快車 項目的目的是研發(fā)新一代在軌服務(wù)航天器。 SUMO/FREND目的是通過交會對接

3、和機械臂的抓捕技術(shù)，為失效的客戶衛(wèi)星提供軌道定位服務(wù)。 所有的這些飛行器必須在限定的環(huán)境中運動，SUMO/FREND和軌道快車之后的服務(wù)航天器也需 要在預(yù)定的環(huán)境中與目標(biāo)接觸。為了能夠完成滿足需求的先進的自主任務(wù)，必須找到一種燃料消耗較低同時又能避免碰撞的 軌跡規(guī)劃算法。這是一個很困難的問題，因為這些飛行器相對于它們的對接目標(biāo)一般是從不同的 （初始）軌道開始的；因此，這些航天器的軌跡規(guī)劃必須綜合考慮遠(yuǎn)程軌道機動和近程的臨近操 作，而遠(yuǎn)程機動伴隨有非線性動力學(xué)和燃料消耗問題，避免碰撞對近程臨近操作是極為重要的。軌道機動中的軌跡規(guī)劃技術(shù)已經(jīng)成功使用了很多年，但是這些技術(shù)通常并不能處理避障問 題。另

4、一方面，大部分處理避障的軌跡規(guī)劃技術(shù)都在陸地應(yīng)用，如輪式機器人或機械臂。遺憾的 是，這些算法不能很好地在空間中應(yīng)用，因為它們所基于的假設(shè)在空間并不適用。特別是動力學(xué) 效應(yīng)不能被忽略。而且，在燃料的節(jié)省和保存方面，空間應(yīng)用比地面應(yīng)用更重要。對于給定的軌 道機動問題，到達(dá)時間也是一個很重要的方面，燃料消耗最少的軌跡可能需要無窮大的轉(zhuǎn)移時間。 因此，需要這樣的軌跡規(guī)劃技術(shù)，這些技術(shù)能夠生成滿足軌道和姿態(tài)動力學(xué)限制的路徑，權(quán)衡燃 料消耗和避障距離的約束，滿足寬變的終端范圍限制，能應(yīng)對和躲避移動障礙物，以及能優(yōu)化完 成時間。本文提出了一套基于變分法的算法，能夠滿足以上要求。該算法以懲罰燃料消耗和障礙 物

5、間距離為目標(biāo)函數(shù)，找到一條滿足一系列邊界約束條件且使目標(biāo)函數(shù)最優(yōu)的軌跡，這些邊界約 束條件允許終端位置和時間固定或自由。1.2假設(shè)航天器軌跡規(guī)劃的研究領(lǐng)域極度寬闊，不同的研究者在試圖解決特定的問題時，有時采用不 同的假設(shè)。比如，不同的研究者，要么包括，要么不包括多引力體作為次要影響因素；有些考慮 地球非均勻重力場和牛頓均勻場的作用（稱為J2影響）。類似的，有些研究者試圖直接包含導(dǎo) 航誤差，而另一些人則沒有。一般的軌跡規(guī)劃研究，主要來自移動機器人領(lǐng)域，它們的假設(shè)條件 通常不同于空間飛行器的軌跡規(guī)劃條件。特別是，機器人的研究者們通常（盡管不是一直）忽略 所有系統(tǒng)動力學(xué)，他們假設(shè)機器人能夠精確跟蹤任

6、何給定軌跡，即使軌跡是分段線性的，這是由 于距離很近時動力學(xué)作用可以忽略。反而文獻中主要關(guān)注在約束環(huán)境中找到無碰撞的路徑。顯然， 航天器在約束環(huán)境中的系統(tǒng)動力學(xué)是不能忽略的。本文將兩個領(lǐng)域的問題進行結(jié)合。所以，為了使問題既模型簡單又易于運算，這里采取一定 的假設(shè)。特別是，航天器動力學(xué)采用標(biāo)準(zhǔn)的二體模型，J2項攝動、N體動力學(xué)、大氣阻力等都忽 略不計。這些作用是影響航天器軌跡規(guī)劃的重要因素，近幾年的航天器軌跡優(yōu)化在這方面做了很 多工作。但是，本文的目的是提出一種算法，能夠在一定約束條件下找到轉(zhuǎn)移軌跡，并不是在處 理擾動方面開辟一個新領(lǐng)地，因此這里對攝動項不加以考慮?？傊M管這些擾動會使數(shù)學(xué)模型

7、 變復(fù)雜，但是不能從根本上改變算法。在有障礙物限制的軌跡規(guī)劃問題中，尋找全局最優(yōu)解的方法中有我們熟知的NP-hard問題， 它被列為已知的最困難的計算問題之一。因此，很多有障礙物限制的軌跡規(guī)劃算法，并不會試圖 解決全局問題，而是試圖得到局部最優(yōu)的結(jié)果。在很多情況下，尤其障礙物相對較少時，局部的 最優(yōu)解是可以接受的。本文也考慮這種情況，我們并不會試圖尋找全局最優(yōu)解，因為這對現(xiàn)有的 動力學(xué)模型進行計算太困難了。曾用到過的尋找全局最小值的方法的綜述在下一節(jié)給出。1.3研究背景解決特定的軌道機動問題有很多經(jīng)典的方法，這些方法在軌道力學(xué)領(lǐng)域是大家所熟知的。轉(zhuǎn) 移時間自由的軌道轉(zhuǎn)移，如霍曼轉(zhuǎn)移1和雙橢圓轉(zhuǎn)

8、移2已經(jīng)使用很多年了，同樣的情況還有時間 自由或受限制特殊點攔截3,4。對于更加復(fù)雜的軌道機動，如考慮多引力體情況，圓錐曲線拼接 法通常被認(rèn)為是已知的首選方法。盡管這些算法所獲得的結(jié)果并不是最精確的，但是它能夠為更 復(fù)雜最優(yōu)化方案得到一個比較合適的初始值。這些方法一般都能得到燃料最優(yōu)解，而且考慮動力 學(xué)限制。雖然這些方法假設(shè)推力是脈沖形式且大小不受限，但是考慮到化學(xué)推力器的輸出能力， 以及大多數(shù)軌道轉(zhuǎn)移問題的時間尺度，這個假設(shè)通常是有效的。但是這些方法既不能處理避障問 題，也不能處理姿軌耦合的機動問題。有很多研究將線性或者非線性最優(yōu)方法應(yīng)用在航天器軌跡規(guī)劃中。大體來說，線性方法僅僅 在近距離操

9、作中很有效，因為真正的軌道動力學(xué)，甚至是最簡單的二體問題也還是有很強的非線 性。而當(dāng)兩個航天器在相對于軌道周期的短時間內(nèi)近距離相對運動時，它們的相對運動可以被線 性化，這就是著名的Hill方程。通過這樣的簡化，線性化最優(yōu)方法可以有效使用。例如，目前 有很多基于線性或者非線性的創(chuàng)新方法，解決單個或者多個合作航天器的機動問題。在這些方法 當(dāng)中，對于航天器的軌跡規(guī)劃來說，MILP（混合整數(shù)線性規(guī)劃）法在解決有約束線性最優(yōu)問題時 是一種很有效的方法5問。MILP方法已經(jīng)應(yīng)用在包括碰撞規(guī)避及羽流的接近操作和編隊飛行中。 這種方法必須使用樣條或者其它近似方法對解空間（即所有可能的軌道集合區(qū)間）進行近似。定

10、 義一個目標(biāo)函數(shù)，采用MILP法找到滿足給定約束條件使目標(biāo)函數(shù)最小化的樣條曲線。因此，這 種方法是將軌跡規(guī)劃問題轉(zhuǎn)化為限制性最優(yōu)問題，進而用一種已有的優(yōu)化算法解決。MILP法在 計算方面異常高效。遺憾的是，MILP法似乎不能應(yīng)用于明顯非線性的動力學(xué)問題，因此，它不 能使用在軌道注入或者軌道（平移）和姿態(tài)（旋轉(zhuǎn)）耦合的問題中。單元分解法也已經(jīng)被用在最優(yōu)軌跡規(guī)劃中。單元分解法通常是形成一個離散化的網(wǎng)格狀的狀 態(tài)空間，然后使用一些標(biāo)準(zhǔn)的搜索技術(shù)搜索得到的離散空間。單元分解法已經(jīng)用于解決航天器的 軌跡規(guī)劃問題，著名的有文獻7，他使用基于A *的尋找方法，應(yīng)用于航天器近距離操作，得到 一個六自由度的軌跡

11、。路線圖方法也能夠產(chǎn)生一個鄰接圖表示搜索空間，然后用一個標(biāo)準(zhǔn)的圖像 搜索方法從圖像中得到一個路徑。但是，與單元分解法不同的是，圖像法不是狀態(tài)空間離散分解 的結(jié)果，相反圖像法試圖獲得空間的“內(nèi)在連續(xù)性”。Generalized Voronoi diagrams在航天器在 軌跡規(guī)劃中的應(yīng)用被廣泛接受8,9 這些方法主要用來處理避障問題，它們并不是燃料最優(yōu)，盡 管燃料消耗在處理要求相對小間距的問題時似乎很合理。它們不能解決終端時間或者終點位置限 制的問題，并且不能夠使用實際的推力器模型。也不能解決實際的軌道動力學(xué)（附加擾動）和姿 態(tài)耦合的問題。但是，對于大部分近距離操作，相對軌道動力學(xué)特性可以忽略，

12、而且軌道和姿態(tài) 運動可以解耦，那么，這種方法可以得到很好的結(jié)果。人工勢函數(shù)制導(dǎo)方法是一種最簡單并且最流行的軌跡規(guī)劃方法。在它的基本形式中，建立了 一個叫做勢函數(shù)的數(shù)學(xué)函數(shù)，它在目標(biāo)點上取得全局最小值，在每個障礙物的中心取得局部最大 值（這是一種緊約束）Ji。該路徑的產(chǎn)生是通過一個基于勢函數(shù)的梯度下降算法完成的，然后 將得到的軌跡投影到基本的歐式空間中。人工勢函數(shù)制導(dǎo)方法在航天器軌跡規(guī)劃中的應(yīng)用也很廣 泛1214。以上的這些方法都考慮是實際的推力器，固定的障礙物，自由的到達(dá)時間。它們并沒有直接 最小化燃料消耗，但是他們使用的是Hill方程并且應(yīng)用在近程接近操作中，所以燃料消耗很少。 另一個方面，

13、它們也不適用于遠(yuǎn)程軌道機動，不能處理自由終端位置和移動的障礙物問題。它們 沒有考慮姿態(tài)規(guī)劃，因此就不能解決與旋轉(zhuǎn)目標(biāo)對接的姿態(tài)與軌道耦合問題。最近，有人提出了用于解決多顆衛(wèi)星星座重構(gòu)的軌跡規(guī)劃問題的一種兩階段的優(yōu)化方法。在 這種方法中，可以利用快速擴展隨機樹法得到一個可行但不是最優(yōu)的軌道。該軌道可以作為平滑 或者優(yōu)化算法的初始估計。該算法或者是由針對線性系統(tǒng)動力學(xué)的線性優(yōu)化方法Pl，或者是由 非線性優(yōu)化方法如偽光譜法16組成。這種方法可以被成功應(yīng)用在考慮碰撞和定向約束，為重構(gòu) 衛(wèi)星星座尋找合理的軌道上。與以上方法不同的是，該算法可以尋找到燃料最優(yōu)軌跡，可以處理 姿態(tài)規(guī)劃和移動目標(biāo)的問題。但是該

14、方法僅僅能夠應(yīng)用在軌道動力學(xué)并不重要的深空星座當(dāng)中。也有將非線性優(yōu)化方法應(yīng)用在航天器軌道優(yōu)化問題中。一般的，非線性方法通常使用在單個 航天器在單一引力體的不同軌道間或者在不同引力體間機動多個軌道周期的情況。這些方法可以 分為直接法和間接法?；旧?，直接法是將最優(yōu)軌跡規(guī)劃問題轉(zhuǎn)化為NLP（非線性規(guī)劃問題），通 過采用一種非線性規(guī)劃算法求解。最常用的方法是直接打靶法和直接配點法。直接打靶法主要通 過參數(shù)化受設(shè)計器控制的變量（通常是控制輸入）實現(xiàn)。系統(tǒng)動力學(xué)模型從初始時刻至終端時刻 計算，進而計算結(jié)果的優(yōu)化值；然后反復(fù)調(diào)整參數(shù)對優(yōu)化結(jié)果進行改進。直接配點法通常被認(rèn)為 是另一種可供選擇的方法，軌跡和控

15、制輸入都通過樣條函數(shù)來參數(shù)化，進而調(diào)整樣條函數(shù)的系數(shù) 來改進優(yōu)化結(jié)果。但是，Betts指出，配點法也可以看作是打靶方法的一個特例，就是將打靶法 的積分器采用Hermite-Simpson算法，并且參數(shù)化的參數(shù)數(shù)量與積分器時間步長數(shù)量相一致。間接法則是推導(dǎo)一組微分方程：Euler-Lagrange方程，方程的解即是最優(yōu)軌跡規(guī)劃問題的極 值曲線。簡單的來說，必須解決邊值問題去得到這些解。相比直接法是一個優(yōu)化問題，間接法是 一個尋根問題。像直接法一樣，有很多解決邊值問題的數(shù)學(xué)方法，包括間接打靶法，間接多重打 靶法，間接配點法。直接法漸漸的成為軌跡優(yōu)化問題中的首選方法是有原因的。首先，直接法比 間接法

16、更具有魯棒性，有更廣泛的收斂范圍。此外，由于間接法以伴隨變量的形式引入其他狀態(tài) 變量，一定要保證對這些變量進行合理估計才會使迭代收斂。但是由于伴隨變量并沒有直接的物 理意義，所以提供一個好的初始估計往往很困難。并且使用間接法要求對每個問題推導(dǎo) Euler-Lagrange方程，這就使得這種方法在完成一般的軌道優(yōu)化中更具有挑戰(zhàn)性。但是，在某些情況下，間接法比直接法更有優(yōu)勢。特別是它能直接（至少在理論上）包含幾 乎所有假設(shè)約束條件，不用事先用數(shù)值解對約束條件進行修正。使用直接法，通過數(shù)值解直接滿 足約束條件，因此，它的解要么很籠統(tǒng)，要么僅適合于特定的問題。特別是在處理較強約束條件 的問題上，例如使

17、用自由開關(guān)推力器時，就很難使用直接法，因為在解決這個問題之前需要知道 推力器點火的時刻，這樣才能在整個過程中處理變化中的系統(tǒng)動力學(xué)。間接法在解決非常困難的 軌道優(yōu)化問題中相當(dāng)成功，比如小推力的多體轉(zhuǎn)移問題17,18。與直接法類似，這些方法通常是基于牛頓法的變形形式。遺憾的是，這些方法的普遍特點是 其局部最優(yōu)性；它們?nèi)啃枰锰荻刃畔ⅲ瑢μ荻鹊哪；蛑笜?biāo)函數(shù)值進行迭代改進。結(jié)果都是 趨于收斂到局部最小值，在不對目標(biāo)函數(shù)形式作附加要求或者滿足微分、代數(shù)約束時，結(jié)果不能 保證收斂到全局最小值。對與這些優(yōu)化算法在航天器軌跡優(yōu)化中應(yīng)用的綜述，見文獻19。直接 法和間接法在“六自由度軌道機動”章節(jié)有詳細(xì)討

18、論。以上提到的方法中，沒有哪個方法是普遍全局最優(yōu)的。然而，有很多人試圖解決全局最優(yōu)問 題，至少是在統(tǒng)計的意義上。模擬退火法就是其中之一2。21。在模擬退火法中，隨機擾動可以應(yīng) 用于最優(yōu)解最好的估計。擾動大小的變化是依據(jù)“退火過程”，擾動大小本質(zhì)上是隨著時間減小的一個函數(shù)。此外，擾動的大小還是隨著當(dāng)前目標(biāo)函數(shù)的值而減小的函數(shù)。這樣的話，這種算法 可以用統(tǒng)計學(xué)的方式搜索到很大部分的解空間。遺憾的是，很難嚴(yán)格保證軌跡規(guī)劃問題的隨機退 火的全局收斂性，而且選擇合理的退火過程通常需要反復(fù)試驗。尋找全局最優(yōu)的另一種方法是遺傳算法22對于不同類型的遺傳算法通常使用不同的變量， 但是一般來說首先需要為可能的問

19、題解建立一個編碼表。在軌跡規(guī)劃中，編碼表由一些基向量的 權(quán)重或神經(jīng)網(wǎng)絡(luò)的權(quán)重組成。一組候選解的產(chǎn)生，通常是從編碼表中隨機選出。這些候選解通過 一個適應(yīng)性函數(shù)進行評估，選擇一部分子集作為新的候選解。對這些留存的候選解隨機的擾動或 交叉，會產(chǎn)生新的候選解集合。反復(fù)迭代，直到達(dá)到最大設(shè)定迭代次數(shù)或候選解達(dá)到滿意的品質(zhì)。 遺傳方法很普遍，可以應(yīng)用到任何其他優(yōu)化方法可以應(yīng)用到的地方。但是，遺憾的是，遺傳方法 不能很好利用梯度信息，因此對問題的計算量很大。因此還沒有被廣泛應(yīng)用到軌跡規(guī)劃問題中。1.4數(shù)學(xué)基礎(chǔ)1.4.1直接法形式。本文使用Lagrange表述，目標(biāo)函數(shù)取成積分形式，并對末端時刻f的狀態(tài)進行懲

20、罰:J L(),u(),= j七 L L(),u(),山o直接法解決最優(yōu)問題開始是構(gòu)造一個目標(biāo)函數(shù)和一些約束條件。目標(biāo)函數(shù)有很多不同的規(guī)范動態(tài)約束條件：z(z(),u；根據(jù)待解決的問題本身的特點，狀態(tài)變量雙)、控制量 u也需要代數(shù)約束：g - gL,uo,k為。最優(yōu)問題轉(zhuǎn)化為求解最優(yōu)控制量u*()和相應(yīng) 狀態(tài)量x*()，通過滿足一些約束條件使得目標(biāo)函數(shù)最小。有很多特定的情況，其中一種典型的情況是狀態(tài)變量和輸入控制量都是二次的，微分約束條 件是線性的，沒有代數(shù)約束條件。這樣的優(yōu)化問題可以求得解析解。對于不滿足這些標(biāo)準(zhǔn)的問題， 可以通過數(shù)值計算得到數(shù)值解。事實上這是相當(dāng)復(fù)雜的，除了相當(dāng)少的情況下，

21、允許的控制輸入 是一個無限維數(shù)的函數(shù)空間。直接法通過限制函數(shù)控制量輸入搜索空間使得最優(yōu)化問題容易計 算。這個近似函數(shù)通常是B樣條。因此，軌道最優(yōu)問題可重新表示為：u* w w*b(t)% L z(t), wtb(t), t dto切 g V g r x (t), wb(t), t ! g和j L滿足約束：” 2 M(t), wt,u。因此，問題變?yōu)镹LP(非線其中，b(t)是基函數(shù)向量，w是常數(shù)權(quán)重向量。尋找權(quán)重向量w*最小化下式: J r Z (), wtb(t), t = j性規(guī)劃問題)。有很多解決非線性規(guī)劃問題的方法，包括：直接打靶法、多次直接打靶法和直接配點法。1.4.2間接法Eule

22、r-Lagrange 方程本文采用基于間接配點法的變分方法。變分法的提出，是為了解決找到使特定約束最小化(或 最大化)的函數(shù)。變分問題描述如下：23,24J z(t),u(t),t= jlf Lz(t),u(t),tdt問題一：給定指標(biāo)函數(shù)t0，找到最優(yōu)軌跡z*(t)和輸入控制量u*(t)，使得指標(biāo)函數(shù)加最小(或最大)，并且滿足動態(tài)約束條件)=f z,(t),t。變分法的軌跡規(guī)劃允許對燃料消耗和避障進行折中。使用間接優(yōu)化方法允許非線性動態(tài)約束 和開關(guān)控制，滿足推力器閾值約束。此外，到達(dá)時刻和終端位置都是最優(yōu)的。這里，數(shù)學(xué)基礎(chǔ)要 求總結(jié)出一個可以完成6自由度軌跡規(guī)劃的間接配點法方法?！伴g接”這個

23、詞是用來區(qū)分這兩種方法，其中一種試圖通過直接在適當(dāng)?shù)暮瘮?shù)近似空間中近 似最小指標(biāo)函數(shù)J日，如順序二次規(guī)劃法(“直接”法)；另一種試圖通過解一組微分方程來 最小化指標(biāo)函數(shù)(間接法)。間接變分優(yōu)化方法的目的是要推導(dǎo)一組微分方程，它的解可使J日 最小。首先，通過在指標(biāo)函數(shù)中加入時變Lagrange乘子項來強制系統(tǒng)動態(tài)約束：J z(t),v(t),t=j Lz(t),v(t),t+Xt(t)f z(t),v(t),tldt jf L z(t),v(t),tdtt0其中，e = f 如,知扇a,at0Lagrange乘數(shù)* (t)叫做協(xié)態(tài)向量，它們可以以Euler-Lagrange方程的解的形式表示:X

24、 (t)=牛z *(t), v* (t), X (t), t(2)J d!z (t ) = 7 .7 (t ) = z最小化 ，需滿足邊界約束條件z(0) z0,z(f) zf和動態(tài)約束條件 z = (t) = fz(t),v(t),t, 一,24,25。注意到，Euler-Lagrange方程是與邊界條件耦合的二階微分方程組。除了一些特殊情況，即當(dāng)系統(tǒng)動態(tài)函數(shù)f (D是線性的，并且指標(biāo)函數(shù)僅包含二次型項，Euler-Lagrange方程通常是非線性的。如果推力器不受限制，可以通過下式解得最優(yōu)控制輸入少(t) 250 = a z(t), V(t), tav如果推力器是有飽和的或者有其他限制(如

25、bang-bang推力器)，最優(yōu)控制輸入V*(t)根據(jù)龐德 里亞金極小值原理25 來選擇：L z(t), X(t), v*(t), t L z(t), X(t), v(t), t 其中，許可控制輸入v(t)是滿足約束條件的任意輸入值。另一種情況，通?？赡軙ㄟ^狀態(tài)向量 和協(xié)態(tài)向量得到v*(t)的解v*(t)=nz *(t ),X*(t)。解決任何微分方程的問題，都需要有邊界條件。對于變分問題，邊界條件通常是分開的，部分邊界條件在t=0時是已知的，部分邊界條件.t=t在 f時是已知的。這就被稱為BVP (邊值問題)。確切的邊界條件取決于問題本身。例如，如果到達(dá)時間、起始點位置和終點位置給定，則邊

26、界條件可以通過以下形式表示：z (t0) = Z0z (t f ) = Zf如果終端時刻自由，終端位置取決于終端時刻，即z =0，則邊界條件是I-1dOL rz(t ),v*(t ),X(t ),t -Xt(t )(t ) = 0Lf f f f f dt fz (t f ) =O (t f) z (t0) = z01.4.3配點法求解邊值問題有很多方法。最流行的是打靶法和它的變體25。打靶法必須為初始條件的未 知部分確定一個初始估計，對微分方程從0到tf積分，用理想值與實際值的差，用Newton法或 者其他相關(guān)方法更新估值。打靶法相當(dāng)精確，但是對初始條件估計的準(zhǔn)確性很敏感。遺憾的是， 協(xié)態(tài)方

27、程沒有一個明顯的物理意義，因此確定一個精確的初始估計可能很困難。配點法雖然不能和打靶方法一樣精確，但是它對初始解估計的準(zhǔn)確性沒有那么敏感。配點法 使用三次樣條或者其他的逼近方法去逼近微分方程的解2629。z(t) w cTb(t)zX (t) w cTb(t)這里向量Cz和向量CX是權(quán)重矢量，b(t)為基函數(shù)向量。每個基函數(shù)都以離散的時間節(jié)點t t t ., 一 cc .(t0，t1，.，m-1)為中心。配點法選擇權(quán)重矢量z和cX，使得插值解滿足邊界條件和每個節(jié)點上 的微分方程。cb(t) = fb(t), 丫 cTb(t), cTb(t), t , tz i L z i X i i i/J

28、- r-|cTb(t) = - 一 Ub(t), cTb(t), 丫 cTb(t), cTb(t), t , tX i6zz i X iz i X iiicTb（t0） = z0crb（t 1） = Zf其中，前兩個方程來自配點法要求，后兩個方程確保滿足邊界條件。邊界值問題由此轉(zhuǎn)化為求解 權(quán)重矢量Cz和土的一個AME（代數(shù)矩陣方程）。對于一個4 -m的代數(shù)方程，其中n是狀態(tài)空 間的大小，m是節(jié)點的數(shù)量，它的結(jié)果是：0 =中（cz, zJ = W （提該方程必須解出系數(shù)&。為了完成轉(zhuǎn)化，假設(shè)微分方程的右邊是利普希茨條件。對本文來說， AME通過Levenberg-Marquardt方法進行求解，

29、該方法同時具有梯度下降法的魯棒性和Newton 法的快速收斂性。為了能夠進一步加強配點法的魯棒性，AME求解器可以融合延拓算法。這種算法可以用來 針對高非線性問題改進配點法的收斂性2729。用延拓修正方程，再去求解，例如，可以通過在 微分方程中為非線性項乘以一個增益。該增益起初為零，成為一個線性問題，這樣可以很容易求 解。在每一步的迭代中，增益不斷增加，前一步的解作為當(dāng)前步的初始估計使用。當(dāng)增益達(dá)到一 致之后，這個方法就可以找到原問題的解。1.4.5姿態(tài)描述為了規(guī)劃六自由度軌跡，必須在系統(tǒng)狀態(tài)空間中增加姿態(tài)的描述。所有的最小的（三參數(shù)） 姿態(tài)描述都在取值范圍內(nèi)某處有奇點。這一事實部分解釋了非最

30、小描述，如四元數(shù)法，普及應(yīng)用 的原因，因為它沒有奇點。但是遺憾的是，對于軌跡規(guī)劃問題，使用非最小姿態(tài)描述的方法會帶 來一個很嚴(yán)重的問題：配點法的AME雅克比行列式是奇異的，這樣就會給基于梯度的如Levenberg-Marquardt的AME求解器帶來問題。克服這種問題的一個解決辦法就是在有利的地方 使用具有動力學(xué)奇異的最小姿態(tài)描述。一種描述方法就是修正的Rodrigues向量，一個相對新的 三參數(shù)描述方法，僅僅是在從原姿態(tài)旋轉(zhuǎn)360時才會出現(xiàn)奇點30。修正Rodrigues向量的元素=+ 2氣卞，用四元數(shù)q=m，河定義為:8 b =1 +門2b8 =相反的:1 + b Tb1 -b Tb1 +

31、 b T b修正Rodrigues向量與Euler軸k和角度之間的關(guān)系如下，,/中、b = k tan()4其中，定義Euler角的范圍0 V中 . . . U其中，p代表在非旋轉(zhuǎn)坐標(biāo)系中的航天器位置，此坐標(biāo)系的原點固連于對接目標(biāo)航天器；v為相 關(guān)的平動速度向量；at為對接目標(biāo)的加速度，b為修正Rodrigues向量，為在體坐標(biāo)系中 的角速度向量；m為航天器的質(zhì)量，H為轉(zhuǎn)動慣量向量。矩陣*(b)為旋轉(zhuǎn)矩陣，將體坐標(biāo)系轉(zhuǎn) 化為慣性坐標(biāo)系。該矩陣可方便地用等價的四元數(shù)法表示為R(q)=(中-t8)1 + 2跖t 2qS()1.5六自由度軌道機動本節(jié)推導(dǎo)了有障礙物情況下姿軌耦合軌道機動問題的指標(biāo)函數(shù)

32、。算例中將用到的指標(biāo)函數(shù)具有以下形式J=七Lz,v, T + Lz + L+ 人tf (z,v, T)dtt controlobstacletime,，V、一 1=ftf |u(t)| + t(t)tRt(t) + aEg (v) h (p) + y +人 tf z,呼(z, X) I dt TOC o 1-5 h z t0 L 1i=1 (7)式中f ()表示航天器動力學(xué)方程(6) ； Lcontrol =U(t)1 + T(t)T Wt為控制力的罰項；L =ag( v)h (p)T =obStaCl(；i=1i i 為障礙間距罰項；time-為完成時間罰項。1.5.1燃料消耗當(dāng)軌跡指標(biāo)函數(shù)

33、具有以下形式時，軌跡燃料消耗最小25J = f tf |u(t )| + t1(8)然而，由該指標(biāo)函數(shù)得到的控制律:u X -1/ mu * (t) = 0-1/ m X 1/ mI maxv ,i是不連續(xù)的，因此導(dǎo)致歐拉一拉格朗日方程右端也不連續(xù)。為克服該問題，需使用與方程(9 )近似的連續(xù)控制律(見圖1)圖1控制律10 max 1。和m 1_ mu (e -1)1 2( e 1)e + 2m 2u(u m2 一 1)e + um2 一 122m(um2 + 1)e + u m2 一 1 TOC o 1-5 h z maxmax._ (2u2 m2 一 u)e1 (um2 一 1)e + u

34、m2maxmax._ (um2 一 1)e 2 一 (3u m2 一 1)eI = maxmax22m2(um2 + 1)e + u m2 一 1maxmaxu (t)=umas X+i1 v,i1s X+i2 v ,i21 f2m 2v ,iX ,i - i2S2X ,i - i1umaxX -2m - u-2m-u X , -(1 E)m-u +e /m-(1 -詞-u+e /m X . -1/m-1/ m X 1/ m1/ m X . (1 -e )m - u +e/m (1-e)m-u+e/m X 2m - u(10)當(dāng)式中e等于1時，上述控制律與最小燃耗律完全相同。本文中大部分例子將

35、使用該連續(xù)方法，e取值范圍從0開始增長到e約等于1。姿態(tài)執(zhí)行機構(gòu)最小能量指標(biāo)函數(shù)如下:J = j tf p T(t)t t(t) +t0式中，p為正常數(shù)。由此可得控制律為TmaxT*(t) = (人tH-1) /2p i廣，max(人tH-1) 2p T-2 p T (人tH-1) 2 p t(11)1.5.2障礙間隔作為燃料消耗懲罰項的補充，在指標(biāo)函數(shù)中引入障礙距離懲罰項:(12), P t , 一 j . 一一 ，一式中，Oj為與第j個障礙固聯(lián)的坐標(biāo)系原點位置O ，b 一 ,，一j與Oj表示障礙坐標(biāo)系的方向。d(.)為航天器到障礙的距離，是障礙物外形的函數(shù)；5j為常值標(biāo)量增益。障礙罰函數(shù)g

36、()值隨著到障礙物邊緣的距離的減小而增大，且當(dāng)與障礙物的距離為某值時，罰函數(shù)值很小或為0。為使歐拉一拉格朗日方程滿足利普希茨光滑條件，英嚴(yán)七須滿足利普希茨條件；例如當(dāng)g 為C 2時，該條件成立。一種滿足這樣要求的方程由分段三次函數(shù)組成:d = d (p, p , c ) sn rOjOjk (a d 3 + a d 2 + a d + a )1 n 2 n 3 n 4式中，S通過對相對距離dr進行縮放決定了 g ()非零區(qū)間的大??；k決定了罰函數(shù)的最大值;數(shù)3，a4的選擇標(biāo)準(zhǔn)是使g ( 40 k)，g (s) = 0島(d )/dd (s) = 0a2g(d)/dd2(s) = 0胞?r rr

37、 r；見見圖2。圖2分段三次障礙罰函數(shù)，光滑化的速度標(biāo)準(zhǔn)化函數(shù)有兩種策略可以使方程(12)最小。一種是通過使軌跡離開障礙物以減小障礙物間距。另一種 是通過增大接近速度以減少接近障礙物消耗的時間(事實上，只需令-0就能使消耗任意小)。顯然，第二種方法是不適用的。作為對相對于障礙物速度增加的補償，需要引入：0j=1(13)J()丁 +a jgdr (p，pq .)叫(p，pq ，v，氣.)+j jj jj j,“，_-,V 式中，0.為第j個障礙坐標(biāo)系的速度，0.為其角速度。航天器與障礙物的相對速度r可由相 對距離對時間求導(dǎo)得出，且速度標(biāo)準(zhǔn)化函數(shù)”隨七的增大而增大。很容易看出，要對速度改變進行精確

38、校正，h(Vr )應(yīng)為七；但是絕對值的偏導(dǎo)在V廣0時 無意義，因此，歐拉一拉格朗日方程在七=0處也無意義。為解決此問題，可以在= 0附近對絕對值函數(shù)進行三次插值以使其平滑:h(v )=rb v 3 + b v 2 + b v + b1 r 2 r 3 r 4V 其它其中，系數(shù)b1,，b4需滿足h( )邛dh(v )/dv () = 152h(v )Mv2() = 0，r r，r rJ且枷匕)/叩(0) = 0 (見圖2)。由于仞牧,)/存在且有界，該方程滿足利普希茨條件。1.5.3歐拉一拉格朗日方程與性能函數(shù)對應(yīng)的歐拉一拉格朗日方程為:V = bl |3口 p + R(b)u(t)(.a.=扁

39、F入=X +E vp j=iq = G (q)wq、ppt入1+ E ah(v )j=iOh(v) Ova t.jOvr, jOg(d ) OdOh(v ) OvFTrg (d )OdOpOvOpr ,jW = H-1S(Hw)w + T(t)O哄wOq0S(Hw)wOwTH-1,t XG (q)人 q q(14)g,h,d的大小如前文定義。式中，U(t)由方程(10)計算，T(t)由方程(11)計算。1.6結(jié)果翻滾目標(biāo)對接下面的例子描述了一個對自主衛(wèi)星服務(wù)來說很現(xiàn)實的任務(wù)，但超出了現(xiàn)有航天器軌跡規(guī)劃程 序的能力。一個與ISS同軌道且在其后方約63km處的衛(wèi)星發(fā)生故障并開始翻滾。將衛(wèi)星視作直

40、徑3m的球體，且兩側(cè)各有一沿y軸伸出體外7m的太陽能電池板。電池板寬3m。衛(wèi)星上對接 點位于兩太陽能電池板正中的星體上。衛(wèi)星相對于慣性系的俯仰角速率為0，即相對ISS的俯仰 角速率約為1.14X 10-3rad/s （每92min旋轉(zhuǎn)一周）。其偏航角速率約為5.2X 10-2rad/s （每2min 旋轉(zhuǎn)一周）。要求設(shè)計一條服務(wù)星軌跡，軌跡的起點位于服務(wù)星在ISS上對接位置，距ISS中心 的距離為10m，且面向地球，終點是服務(wù)星要與翻滾衛(wèi)星對接。翻滾衛(wèi)星如圖3所示。圖3翻滾衛(wèi)星模型（衛(wèi)星俯仰角速率為1.14X10-3rad/s,偏航角速率為5.2X 10-2rad/s。初始軌道速度方向沿z軸正

41、向）假設(shè)機器人服務(wù)航天器的質(zhì)量為2000kg，轉(zhuǎn)動慣量1=32。 ml S 2。沿各軸 的最大作用力為1000N，最大作用力矩為2000Nm。力矩罰函數(shù)增益是10，時間增益Y是 0。如前例，采用千米、分和服務(wù)衛(wèi)星質(zhì)量單位構(gòu)成非標(biāo)準(zhǔn)單位制對問題進行描述。而在表示結(jié) 果時，須轉(zhuǎn)化回標(biāo)準(zhǔn)單位制。結(jié)果如圖4所示。該軌跡需要81611際的推力，66.5Nms的力矩?？倷C動時間為2728s （45.46 min ）。整個軌道攔截軌跡見圖4。衛(wèi)星最后逼近段的時序顯示在圖5中。近距離觀 察顯示服務(wù)衛(wèi)星在對接前“起旋”并持續(xù)了數(shù)分鐘，以應(yīng)對目標(biāo)衛(wèi)星的偏航轉(zhuǎn)動。圖4服務(wù)衛(wèi)星對接問題軌跡。每隔約30s顯示一次服務(wù)衛(wèi)

42、星體坐標(biāo)軸圖5服務(wù)衛(wèi)星逼近軌跡，每2s顯示一次服務(wù)衛(wèi)星體坐標(biāo)軸該算例證明此軌跡規(guī)劃程序能夠解決需要質(zhì)心運動與姿態(tài)運動高度協(xié)調(diào)的問題。在前面的例 子中，也可以獲得滿足bang-coast-bang （開-關(guān)-開）推力器約束的軌跡。最后，這個例子也說明 了該軌跡規(guī)劃器能解決航天器飛行足夠遠(yuǎn)距離的問題，這種情況下不能使用Hill方程或其他線性 化相對運動軌道方程。因此該軌跡規(guī)劃器能夠解決那些現(xiàn)有的軌跡規(guī)劃器所不能解決的問題。1.7結(jié)論由以上算例可以看出，文中論述的變分軌跡規(guī)劃算法能夠滿足引言中提出的要求。該算法具 有下列特點：能夠在存在靜止或運動障礙的情況下獲得自由路徑；能夠在有引力場和其它環(huán)境力作

43、用的情況下尋找得到最小燃料消耗軌跡；能夠處理大范圍端點約束與時間尺度問題；考慮了現(xiàn)實的航天器動力學(xué)模型；考慮了現(xiàn)實的推力器，包括bang-coast-bang(開-關(guān)-開)推力器。同時，本文使用的動力學(xué)模型并非其它大部分軌跡規(guī)劃法所采用的線性化模型。因此，相比 其它軌跡規(guī)劃法，該算法能獲得更大距離范圍的最小燃料消耗軌跡。另外，變分軌跡規(guī)劃法不局 限于質(zhì)心運動軌跡規(guī)劃；如最后一個例子所示，該方法能夠解決姿軌緊耦合問題。綜合以上特點, 可以看出，本文介紹的變分軌跡規(guī)劃法對當(dāng)前的軌跡規(guī)劃有著深遠(yuǎn)的意義。參考文獻Hohmann, W.(1925) Die erreichbarkeit der himm

44、elskarper (the attainability of celestial bodies), Munich, Germany.Hoelker, R. F., and Silber, P. S.(1961) “The bi-elliptical transfer between coplanar circular orbits.” Proc., Fourth AFBMD/STL Symp., Vol. 3, Pergamon, New York, 164-175.Tempelman, W.(1961) “Minimum-energy intercepts originating from

45、 a circular orbit.” J. Aerosp. Sci., 28, 924-929.Gedeon, G. S. (1958). “Orbital mechanics of satellites. Proc., American Astronautical Society, American Astronautical Society.Richards, A., Schouwenaars, T., How, J. P., and Feron, E.(2002). “Spacecraft trajectory planning with avoidance constraints u

46、sing mixedinteger linear programming. J. Guid. Control Dyn., 25_4_.Breger, L. S., and How, J. P.(2006). “Safe trajectories for autonomous rendezvous of spacecraft. Proc., 2006AIAA Guidance, Navigation and Control Conf., American Institute of Aeronautics and Astronautics.Jackson, M. C. (1994). “A six

47、 degree of freedom, plume-fuel optimal trajectory planner for spacecraft proximity operations using an anode search. MS thesis, Massachusetts Institute of Technology, Cambridge, Mass.Choset, H., Fernandez, J. J., and Kortenkamp, D.(1999). “Path planning and control for aercam, a free-flyinginspectio

48、n robot in space. Proc., 1999 IEEE Int. Conf, on Robotics and Automation (ICRA99), Institute of Electrical and Electronics Engineers _IEEE_, Detroit, 1396-1403.Choset, H., La Civita, M., and Park, J. C. (1999). “Path planning between two points for a robot experiencing localization error in known an

49、d unknown environments. Proc., Int. Conf, on Field and Service Robotics, Springer, Pittsburgh.Khatib, O. (1980). “Commande dynamique dans lespace operationnel des robots manipulateurs en presence dobstacles. Ph.D. thesis, Ecole Nationale Superieure de lAeronautique et de lEspace, Toulouse, France.Kh

50、atib, O. (1986). “Real time obstacle avoidance for manipulators and mobile robots. Int. J. Robot. Res., 5_1_.McInnes, C. R. (1999). “Autonomous path planning for on-orbit servicing vehicles.” Reducing space mission cost, American Institute of Aeronautics and Astronautics.McQuade, F., and McInnes, C.

51、 R. (1997). “Autonomous control for onorbit assembly using potential function methods.” Aeronaut. J., 101_1006_, 255-262.Roger, A. B., and McInnes, C. R. (2000). “Safety constrained free-flyer path planning at the international space station” J. Guid. Control Dyn., 23(6): 971-979.Garcia, I. M., and

52、How, J. P. (2005). “Trajectory optimization for satellite reconfiguration maneuvers with position and attitude constraints. Proc., 2005 IEEE American Control Conf., Institute of Electrical and Electronics Engineers _IEEE_, 889-895.Aoude, G. S., Garcia, I. M., and How, J. P. (2007). “Two-stage path p

53、lanning approach for designing multiple spacecraft reconfiguration maneuvers. Proc., 20th Int. Symp. on Space Flight Dynamics.Nah, R. S., and Vadali, S. R. (2001). “Fuel-optimal, low-thrust, threedimensional Earth-Mars trajectories. J. Guid. Control Dyn., 24.Raniera, C., and Ocampa, C. (2005). “Optimization of roundtrip, timeconstrained, finite burn trajectories via an indirect method. J. Guid. Control Dyn., 28.Betts, J. T. (1998). “Survey of numerical methods for trajectory optimization. J. Guid. Control Dyn., 21.Kirkpatrick