胡靖，高普云

(國防科學(xué)技術(shù)大學(xué)航天科學(xué)與工程學(xué)院，湖南長沙410073)

0 引言

神舟八、九、十號飛船與天宮一號交會對接飛行試驗的順利完成，標(biāo)志著我國交會對接技術(shù)逐漸走向成熟應(yīng)用［1］。交會對接以航天器首次接觸為劃分，可分為交會和對接兩大問題。圍繞這兩類問題，已有大量研究成果和較為成熟的、公認(rèn)的動力學(xué)研究方法和模型。如交會段一般只須考慮航天器的軌道機動，可以將航天器視為質(zhì)點，應(yīng)用形式簡便的CW方程［2］設(shè)計控制方法［3］和分析航天器的相對運動［4］;對接段由于兩航天器對接機構(gòu)已相互接觸，彼此間存在碰撞或接觸力的相互作用，一般應(yīng)用多體動力學(xué)建立碰撞模型［5］研究對接機構(gòu)及航天器的力學(xué)性質(zhì)［6-7］。然而作為交會和對接的過渡階段，航天器對接逼近段又有著不同于交會段和對接段的控制要求和動力學(xué)特點。也有不少文獻(xiàn)［8-10］對逼近段航天器的控制策略和相對位置姿態(tài)估計進(jìn)行了研究。雖然已取得一定的研究成果，并使航天器安全地互相逼近在交會對接實際操作中得以實現(xiàn)，但是針對逼近段的研究仍存在一些理論與應(yīng)用上的局限性。逼近段追蹤器與目標(biāo)器尚未接觸，仍屬于兩個獨立的在軌航天器，需要同時控制追蹤器軌道和姿態(tài)，避免逼近過程中和首次接觸時的碰撞，因而需要建立能夠精確描述相對位置和姿態(tài)的模型。特別是對接航天器之間的連接點(分別處于追蹤器和目標(biāo)器對接機構(gòu)上，在對接段對應(yīng)相連接的兩點)，其相對運動參數(shù)是對接成功與否的關(guān)鍵。而連接點位于對接機構(gòu)的前端位置，也是逼近段最有可能發(fā)生意外碰撞的部位。

本文考慮建立包含連接點相對位置、速度及航天器姿態(tài)、角速度等12個獨立狀態(tài)變量的逼近段動力學(xué)方程，并嘗試在符合逼近段運動參數(shù)的簡化條件下，求解各狀態(tài)變量的解析表達(dá)式。在此基礎(chǔ)上完成逼近段追蹤器軌姿的控制設(shè)計和仿真。最后通過仿真算例，驗證所建立模型和所設(shè)計控制方法的適用性，并給出一些關(guān)于逼近段相對運動和控制的分析結(jié)論。

1 逼近段的動力學(xué)模型

1.1 模型假設(shè)和坐標(biāo)系

假定目標(biāo)器在既定軌道上運行，不進(jìn)行機動。逼近段初始時刻追蹤器到達(dá)目標(biāo)器的軌道高度，接近目標(biāo)器的位置。追蹤器的軌姿控制方式為多次線速度和角速度脈沖，以每次脈沖速度增量或機動時刻為控制變量。引入如下坐標(biāo)系研究航天器的運動:

(1)慣性坐標(biāo)系Oxyz，原點為地心，z軸指向軌道上的一個定點，x軸處于軌道平面內(nèi)并垂直于z軸，y軸與它們構(gòu)成右手系。

(2)體坐標(biāo)系C(k)x(k)y(k)z(k)，原點固連于航天器質(zhì)心，z(k)軸沿航天器軸向，x(k)軸處于垂直于追蹤器軸向的平面內(nèi)，并指向?qū)訖C構(gòu)某一連接點在該平面上的投影，y(k)軸與x(k)，z(k)軸構(gòu)成右手系。文中只考慮它們均為慣量主軸的情況。上標(biāo)k=1，2分別對應(yīng)于追蹤器和目標(biāo)器。

式中:i=1，2，…，n為角速度機動的序數(shù)，n為角速度機動次數(shù)為沿坐標(biāo)軸的單位矢量。

將Hi在中分解(見圖1)，可得:

其中:

式中:Θi，Φi表示由轉(zhuǎn)向的章動角和自轉(zhuǎn)角;為追蹤器繞各軸的轉(zhuǎn)動慣量;為航天器角速度ω(k)在體系下的分量;Δωu1(i)為第i次角速度增量Δω(i)在追蹤器體系下的分量;ti為第i次角速度機動時刻，t0，tf分別表示逼近段初始時刻和航天器首次接觸時刻。

圖1 坐標(biāo)系Fig.1 Coordinate systems

1.2 動力學(xué)方程

地心引力下，航天器質(zhì)心運動由如下方程確定:

式中:rCk為航天器質(zhì)心相對于地心的位置矢量;μ為地球引力常數(shù)。

根據(jù)運動學(xué)知識，可由式(3)得到如下相對運動方程:

式中:ΔrA為 A(1)點(追蹤器上的連接點)相對于A(2)點(目標(biāo)器上與A(1)對應(yīng)的連接點)的位置矢量;rCkAk為A(k)點相對于航天器質(zhì)心的位置矢量。

無外力矩時，航天器的角速度由如下方程確定:

式中:J(1)為追蹤器的轉(zhuǎn)動慣量張量。

2 動力學(xué)方程的求解

2.1 追蹤器姿態(tài)和角速度的求解

形如式(5)的歐拉方程的解可以用橢圓函數(shù)表示，文獻(xiàn)［11］給出了以時為初始時刻解的表達(dá)式?？梢詮钠渲型茖?dǎo)出以機動時刻為初始時刻的解為:

其中:

積分常數(shù)c1，c2，c3由角速度的初始條件確定:

理論上，將上述角速度的解代入歐拉角和角速度關(guān)系式中［11］，可以解得航天器的姿態(tài)角。然而該微分方程的求解比較困難。為避免復(fù)雜的積分運算，這里并不直接計算航天器體系與地心慣性系之間的歐拉角，而是采用以下方法求解航天器姿態(tài)。

首先可導(dǎo)出如下表達(dá)式:

式中:ψi，θi，φi為 ti－1＜ t≤ti時由 C(1)xiyizi轉(zhuǎn)向的歐拉角。其余參數(shù)如下:

積分常數(shù)c4，c5，c6由歐拉角初始條件確定:

記歐拉角表示的坐標(biāo)轉(zhuǎn)換矩陣［11］為 A(ψ，θ，φ)。根據(jù)各坐標(biāo)系間的相互關(guān)系，可得到與的坐標(biāo)轉(zhuǎn)換矩陣(當(dāng)ti－1＜t≤ti時)為:

并得到Oxyz與C(1)x(1)y(1)z(1)的坐標(biāo)轉(zhuǎn)換矩陣為:

式(9)給出了航天器體系與慣性系之間的坐標(biāo)轉(zhuǎn)換矩陣，和式(6)一起為追蹤器姿態(tài)和角速度的解。

2.2 連接點相對運動的求解

其中:

式中:E 為單位矩陣;ΔxA，ΔyA，ΔzA為 ΔrA在Oxyz下的分量;為A(k)點在體系中位置坐標(biāo)。

式(10)的解為:

兩邊對t求導(dǎo)得:

式中:ΔvAx，ΔvAy，ΔvAz為 A(k)點相對速度在 Oxyz下的分量;c1u，c2u為積分常數(shù)，若已知相對位置和速度在Tj－1時刻初始值，代入式中即可確定。從t=0開始依次遞推至Tj＜t≤Tj+1，可整理得:

其中:

式中:下標(biāo)j=1，2，…，N為線速度機動的序數(shù)，N為線速度機動次數(shù);Tj為第j次線速度機動時刻，T0=t0;Δvru(j)為第j次線速度增量Δvr(j)在Oxyz下的分量。

式(13)給出了連接點相對運動參數(shù)的解析解。這兩個方程中的 F1，F(xiàn)2，F(xiàn)3，G1，G2，G3有明確的物理意義，它們分別反映了初始條件、速度增量以及航天器的轉(zhuǎn)動對連接點相對位置、速度的影響。

3 逼近段的終端約束和邊值問題

為避免發(fā)生碰撞，連接點之間首次接觸時應(yīng)達(dá)到相對靜止。而隨著追蹤器逼近目標(biāo)器，兩航天器應(yīng)逐漸達(dá)到相對應(yīng)的對接姿態(tài)，避免在預(yù)定首次接觸時刻前發(fā)生意外碰撞;并在首次接觸時刻達(dá)到相同的姿態(tài)和角速度，使整個對接面重合，對接面上所有對應(yīng)接觸的點均達(dá)到相對靜止。這就避免了首次接觸后航天器的姿態(tài)調(diào)整和接觸部位之間的相對運動，從而消除對接段可能產(chǎn)生的碰撞。上述條件反映為逼近段的終端約束，即:

式(14)共提供了12個獨立的標(biāo)量約束，唯一確定了首次接觸時刻追蹤器相對于目標(biāo)器的位置、姿態(tài)、速度和角速度。

補充以上參數(shù)的初始條件如下:

逼近段邊值問題可描述為:求滿足兩端約束式(14)、式(15)，以及過程約束式(6)、式(9)、式(13)的控制量 t1，…，tn，T1，…，TN或 Δω(1)，…Δω(n)，Δvr(1)，…，Δvr(N)及狀態(tài)量 A(1)(t)，ω(1)(t)，ΔrA(t)，ΔvA(t)。當(dāng) n=N=6時為定解問題;當(dāng) n，N其中之一大于6時，可以選取逼近段時間tf－t0最小為優(yōu)化目標(biāo)，轉(zhuǎn)化為最優(yōu)控制問題。以上所有約束均為代數(shù)方程，不難應(yīng)用序列規(guī)劃方法求解。

4 仿真算例

設(shè)目標(biāo)器在距地面1 500 km高度的圓軌道上運行，姿態(tài)保持為軸向始終沿軌道切向且不繞軸向轉(zhuǎn)動。追蹤器到達(dá)目標(biāo)器的軌道上，初始時刻A(k)點的相對距離50 m，相對速度在軌道面內(nèi)兩個方向分量為4.51 ×10－2m/s，2.58 ×10－3m/s，兩航天器軸向夾角 5.63 rad。追蹤器角速度為(0.200，0.201，0.200)rad˙s－1，繞體系各軸的轉(zhuǎn)動慣量為(5 000，5 000，500)kg˙m2，A(k)(k=1，2)點在各自體系下位置坐標(biāo)分別為(0，0，4)m 和(0，0，－4)m。預(yù)定逼近段角速度機動12次，其中第1，6，7，10，12次角速度增量預(yù)設(shè)為沿體系三個方向分量均為 0.1 rad˙s－1;其余各次為 － 0.1 rad˙s－1。線速度機動6次，其中第2和第3次速度增量1 m/s，沿x軸方向;第1和第4次速度增量－2.5 m/s，沿y軸方向;第5和第6次速度增量－0.935 m/s，沿 z軸方向。以機動時刻為控制變量。

仿真結(jié)果表明，逼近段所需最短時間為67.578 6 s。角速度機動時刻分別為:16.336 7 s，25.192 9 s，28.681 2 s，39.579 2 s，50.584 8 s，51.428 2 s，51.428 2 s，61.270 2 s，61.270 2 s，67.023 5 s，67.023 5 s，67.578 6 s。圖 2 為追蹤器角速度的時間歷程。圖3為矩陣A(1)和A(2)的主對角元素之差隨時間的變化曲線。線速度機動時刻分別 為:3.958 78 s，10.817 4 s，16.397 9 s，41.232 9 s，49.130 1 s，67.338 2 s。A(k)點的相對位置、速度的時間歷程分別如圖4、圖5所示。

圖2 追蹤器角速度時間歷程Fig.2 Time history of chaser’s angular velocity

圖3 坐標(biāo)轉(zhuǎn)換矩陣主對角元素之差隨時間變化曲線Fig.3 Time history of the differences between the corresponding elements of the two coordinate transformation matrixes

圖4 連接點相對位置時間歷程Fig.4 Time history of relative positions

圖5 連接點相對速度時間歷程Fig.5 Time history of relative velocity

從圖2～圖5可以看出，首次接觸時刻追蹤器達(dá)到與目標(biāo)器相等的姿態(tài)角和角速度，并且連接點以零相對速度接觸，無碰撞對接得以實現(xiàn)。其中從圖4可以看出，連接點相對位置各方向分量總是沿著一個方向趨近于零的，特別是沿對接方向的分量在整個逼近過程中始終非正。這表明在逼近過程中，該連接點未與其對應(yīng)連接點發(fā)生碰撞。由于得到準(zhǔn)確的連接點相對位置時間歷程，可以根據(jù)該曲線初步判斷逼近段是否發(fā)生意外碰撞。

觀察圖2可以發(fā)現(xiàn)，沿z1軸角速度分量終端值只與角速度增量的總量有關(guān)，而與機動時刻無關(guān)。這是由于該算例中，追蹤器繞x1，y1軸兩方向的轉(zhuǎn)動慣量相等導(dǎo)致的。實際情況下，對于繞軸旋轉(zhuǎn)對稱的航天器，沿該軸方向的角速度分量都存在這樣的特性。對于此類航天器，所施加的角速度增量的總量不能任意選擇，只能由確定。

仿真結(jié)果中存在第6、7次，第8、9次，第10、11次角速度機動時刻兩兩重合的情況。這可以理解為該時刻只進(jìn)行一次機動，而所需的角速度增量為預(yù)先設(shè)定的兩次增量之和。該算例實際只需進(jìn)行8次角速度機動。比較不同機動次數(shù)設(shè)置的算例仿真結(jié)果，可以總結(jié)出這樣的規(guī)律:如果預(yù)先設(shè)定的機動次數(shù)越多，那么計算所得的機動時刻重合的現(xiàn)象就越多。也就是說，時間最短的優(yōu)化結(jié)果有使多余次數(shù)的機動時刻重疊的趨勢，因此逼近段無須設(shè)置過多的機動次數(shù)。然而針對軌姿控制發(fā)動機推力(力矩)大小和方向難于實現(xiàn)連續(xù)可控變化的情況，卻可以通過這種算法對每次施加的速度增量大小在發(fā)動機可實現(xiàn)的范圍內(nèi)進(jìn)行優(yōu)化選擇。

5 結(jié)束語

本文所建立的動力學(xué)方程以連接點的相對運動參數(shù)為未知量。通過所得的連接點相對距離時間歷程曲線，可以更直觀地判斷逼近段有無碰撞產(chǎn)生。并且在推導(dǎo)過程中，本文給出了多自由度軌姿耦合控制的兩航天器相對運動方程的解析解，逼近段的邊值問題由此得到解決。由此所得計算結(jié)果更有利于分析這些參數(shù)隨初始條件、控制量、航天器幾何參數(shù)等影響因素變化的解析規(guī)律。與離散動力學(xué)方程的數(shù)值求解方法只能在離散點上滿足方程相比，本文所推導(dǎo)的解析解整個過程滿足方程，提高了計算精度;與打靶法相比，避免了猜測初值的困難，顯著提高了計算效率。上述解析解的推導(dǎo)，是基于航天器質(zhì)心的相對距離遠(yuǎn)小于航天器軌道高度假設(shè)的，這完全適用于逼近段，同時也可以應(yīng)用于其他一些滿足該假設(shè)前提的相對運動問題的求解。

［1］ 周建平.天宮一號/神舟八號交會對接任務(wù)總體評述［J］.載人航天，2012，18(1):1-5.

［2］ Clohessy W H，Wiltshire R S.Terminal guidance system for satellite rendezvous［J］.Journal of the Aerospace Sciences，1960，27(6):653-658.

［3］ Lou Y Z，Tang G J，Lei Y J.Optimal multi-objective linearized impulsive rendezvous［J］.Journal of Guidance，Control，and Dynamics，2007，30(2):383-389.

［4］ Kasdin N J，Gurfil P，Kolemen E.Canonical modelling of relative spacecraft motion via epicyclic orbit elements［J］.Celestial Mechanics and Dynamical Astronomy，2005，92(4):337-370.

［5］ Klisch T.Contact mechanics in multibody systems［J］.Mechanism and Machine Theory，1999，34(5):665-675.

［6］ Xavier C，Sunwook K，Michel I，et al.Post-impact dynamics of two multi-body systems attempting docking/berthing［J］.Acta Astronautica，1997，40(11):759-769.

［7］ Bolotnik N N，Shmatkov A M.Control of Impact Interaction of Docking Spacecraft［J］.Journal of Computer and systems Sciences International，2007，46(4):129-136.

［8］ Philip N K，Ananthasayanam M R.Relative position and attitude estimation and control schemes for the final phase of an autonomous docking mission of spacecraft［J］.Acta Astronautica，2003，52(7):511-522.

［9］ 陳統(tǒng)，徐世杰.非合作式自主交會對接的終端接近模糊控制［J］.宇航學(xué)報，2006，27(3):416-421.

［10］楊少東，汶德勝，馬軍勇，等.基于預(yù)測濾波的衛(wèi)星交會對接相對位姿確定方法線［J］.宇航學(xué)報，2011，32(8):1688-1693.

［11］ Landau L D，Lifshitz E M.Mechanics［M］.London:Pergamon Press，1976.