侯振東，曹喜濱，張錦繡

(哈爾濱工業(yè)大學(xué) 衛(wèi)星技術(shù)研究所，150080 哈爾濱，hzd1988120@126.com)

衛(wèi)星編隊(duì)飛行技術(shù)是近幾年提出的空間科學(xué)領(lǐng)域的新技術(shù)，具備系統(tǒng)成本低、性能好、可靠性高和適應(yīng)性強(qiáng)等優(yōu)點(diǎn).但由于地球扁率、大氣阻力等影響，通常需要附加作用力與力矩對(duì)成員衛(wèi)星的相對(duì)位置與姿態(tài)進(jìn)行控制.傳統(tǒng)采用的噴氣推進(jìn)由于能量消耗，影響編隊(duì)飛行壽命，且推力器排出的羽流可能污染鄰近衛(wèi)星的光學(xué)器件表面.而電磁編隊(duì)飛行(EMFF)利用通電線圈產(chǎn)生的電磁力來實(shí)現(xiàn)編隊(duì)構(gòu)型的保持與控制，并由太陽能提供機(jī)動(dòng)能量，且不存在羽流污染，在編隊(duì)飛行中將會(huì)有很大的應(yīng)用前景.

自電磁編隊(duì)飛行的概念提出以來，國內(nèi)外學(xué)者進(jìn)行了先期探索和研究，取得了部分成果.其中，Schweighart［1］提出了近場、中場和遠(yuǎn)場三種電磁力模型.Hashimoto［2］與Kaneda［3］等分析了低地球軌道(LEO)的電磁編隊(duì)飛行動(dòng)力學(xué)模型，并對(duì)在同一平面運(yùn)行的相對(duì)距離固定的雙星系統(tǒng)進(jìn)行了詳細(xì)分析.Elias［4］提出了電磁編隊(duì)飛行的通用非線性動(dòng)力學(xué)模型，對(duì)衛(wèi)星星體、飛輪、電磁體以及飛輪和星體的連接部分都建立了數(shù)學(xué)模型.Elias［5］還研究了穩(wěn)定旋轉(zhuǎn)的雙星系統(tǒng)動(dòng)力學(xué)，該雙星系統(tǒng)繞其質(zhì)心旋轉(zhuǎn)，同時(shí)徑向距離保持恒定，維持這種運(yùn)動(dòng)的徑向力由電磁力提供，在論述中作者將非線性動(dòng)力學(xué)模型在理想軌跡附近線性化，并進(jìn)行了穩(wěn)定性分析.

本文假設(shè)編隊(duì)系統(tǒng)中每顆衛(wèi)星都安裝3 個(gè)正交的通電線圈，通過改變每個(gè)線圈的電流大小、方向，可產(chǎn)生任意大小、方向的電磁控制力與力矩，實(shí)現(xiàn)編隊(duì)構(gòu)型的保持與控制.另外星體上還安裝3 個(gè)正交的飛輪，用來控制星體姿態(tài).產(chǎn)生電磁體的激勵(lì)電流和飛輪旋轉(zhuǎn)所需的能源由太陽帆板提供.本文在此基礎(chǔ)上，考慮大氣攝動(dòng)、地球扁率J2攝動(dòng)，建立電磁編隊(duì)飛行的非線性相對(duì)軌道動(dòng)力學(xué)模型，并對(duì)相對(duì)軌道動(dòng)力學(xué)方程進(jìn)行線性化，得到參考坐標(biāo)系下的線性化相對(duì)軌道動(dòng)力學(xué)方程.

1 坐標(biāo)系定義

為敘述方便，建立以下4 種直角坐標(biāo)系.坐標(biāo)系示意圖見圖1.

圖1 坐標(biāo)系定義

地心慣性坐標(biāo)系FI:原點(diǎn)在地心Oe;XI軸指向赤道平面與黃道平面相交接線的升交點(diǎn)(即春分點(diǎn));地球自轉(zhuǎn)軸為ZI;選擇YI使其與XI、ZI滿足右手坐標(biāo).

主星參考坐標(biāo)系Fr:原點(diǎn)在主星質(zhì)心Or;Yr軸沿地心與主星的連線方向，指向主星;Xr軸指向主星速度方向;Zr軸與Xr、Yr構(gòu)成右手坐標(biāo)系，指向軌道面負(fù)法線方向.

衛(wèi)星軌道坐標(biāo)系Fo:原點(diǎn)在衛(wèi)星質(zhì)心Oo，各軸指向類似于主星參考坐標(biāo)系，只不過衛(wèi)星軌道坐標(biāo)系針對(duì)的是相應(yīng)的衛(wèi)星，而非主星.

線圈安裝坐標(biāo)系Fm:原點(diǎn)在衛(wèi)星質(zhì)心Om(即Oo);Xm、Ym、Zm三軸分別沿星上3 個(gè)正交的通電線圈平面的法線方向.

為敘述方便，記編隊(duì)衛(wèi)星數(shù)目為n，編號(hào)依次為0，1，2…，n-1，其中衛(wèi)星0 為主星.

2 任意兩星間電磁力建模

電磁編隊(duì)飛行的電磁力由衛(wèi)星上的通電線圈產(chǎn)生.如果用fij表示衛(wèi)星B 上第j 個(gè)線圈產(chǎn)生的磁場對(duì)衛(wèi)星A 上第i 個(gè)線圈的作用力，則兩顆衛(wèi)星間的電磁力為

參照文獻(xiàn)［1］，建立遠(yuǎn)場模型.以任意兩顆衛(wèi)星i 和j 為例，記μi為衛(wèi)星i 的磁偶極子矢量，μj為衛(wèi)星j 的磁偶極子矢量，μ0為真空磁導(dǎo)率，Rij=Rj-Ri為兩顆衛(wèi)星的相對(duì)矢徑，Rij為相應(yīng)距離.則衛(wèi)星j 對(duì)衛(wèi)星i 的電磁力與電磁力矩分別為

為便于分析，假設(shè)以反作用飛輪完成衛(wèi)星的姿態(tài)控制，在建模過程中暫不考慮電磁力矩對(duì)相對(duì)軌道運(yùn)動(dòng)的影響.

式(2)中，μi在安裝坐標(biāo)系下可表述為

式中:nix表示安裝坐標(biāo)系Xm軸的線圈匝數(shù);Iix表示Xm軸的電流強(qiáng)度;Aix表示Xm軸的線圈所圍面積;Fm表示分量所在坐標(biāo)系，其他符號(hào)的意義依此類推.

采用zxy 順序的歐拉角轉(zhuǎn)換形式，衛(wèi)星i 的軌道坐標(biāo)系到安裝坐標(biāo)系的坐標(biāo)轉(zhuǎn)換矩陣為

采用zxz 順序的歐拉角變換，地心慣性坐標(biāo)系到衛(wèi)星i 的軌道坐標(biāo)系的坐標(biāo)轉(zhuǎn)換矩陣為

式中Ω 為衛(wèi)星的升交點(diǎn)赤經(jīng)，ˉi 為軌道傾角，f 為緯度幅角.

由參考坐標(biāo)系的定義可知，地心慣性坐標(biāo)系到主星參考坐標(biāo)系的轉(zhuǎn)換矩陣與式(6)相似.

由公式(5)、(6)可得μi在地心慣性坐標(biāo)系下的表達(dá)式為

式中μj在地心慣性坐標(biāo)系下的表述與μi類似.將公式(5)～(7)代入式(2)，可得衛(wèi)星j 對(duì)衛(wèi)星i 的電磁力在地心慣性坐標(biāo)系下的分量形式，將作為仿真程序的理論基礎(chǔ).

由式(2)～(7)可見，兩顆衛(wèi)星間的電磁力與兩星的磁偶極子相耦合，與相對(duì)距離的四次方成反比，具有強(qiáng)烈的非線性.而磁偶極子是兩顆衛(wèi)星姿態(tài)及線圈電流的函數(shù)(一般情況下線圈匝數(shù)與所圍面積是常數(shù)).結(jié)合電磁力產(chǎn)生的機(jī)理與公式分析可知，任意兩星間的電磁力是兩星相對(duì)位置、相對(duì)姿態(tài)及各線圈電流的1 個(gè)非線性函數(shù).

同理，將式(9)代入式(2)中，可得衛(wèi)星j 對(duì)衛(wèi)星i 的電磁力在主星參考坐標(biāo)系下的分量形式.

理論上編隊(duì)衛(wèi)星所受的電磁力應(yīng)包括地磁場的作用.但由參考文獻(xiàn)［6］可知，在軌道高度不低于200 km，編隊(duì)距離幾十米的情況下，地磁場對(duì)衛(wèi)星的電磁力至少比兩星相互作用的電磁力小1 個(gè)數(shù)量級(jí)，因此在電磁力建模時(shí)可以忽略地磁場對(duì)衛(wèi)星的影響.

3 線性化相對(duì)軌道動(dòng)力學(xué)模型建立

相對(duì)軌道動(dòng)力學(xué)方程建立在參考坐標(biāo)系下.參照希爾方程，可得衛(wèi)星i 的相對(duì)軌道動(dòng)力學(xué)方程為

式中:ω 為主星軌道角速度;μe為引力常數(shù);m 為衛(wèi)星質(zhì)量;R0為主星地心距;且

3.1 電磁力線性化

在多顆衛(wèi)星(多于兩顆)的電磁編隊(duì)飛行中，每兩顆衛(wèi)星之間都會(huì)產(chǎn)生電磁力.多顆衛(wèi)星編隊(duì)的電磁力實(shí)質(zhì)上是任意兩顆衛(wèi)星間電磁力的疊加.為方便起見，本文先研究雙星編隊(duì)系統(tǒng)，兩顆衛(wèi)星標(biāo)號(hào)為i、j.由坐標(biāo)系定義可知，衛(wèi)星0 為主星.這里以衛(wèi)星i 為主星來分析衛(wèi)星i、j 之間的電磁力線性化方法.由于衛(wèi)星姿態(tài)由飛輪協(xié)同控制，理論上可將衛(wèi)星姿態(tài)穩(wěn)定在任意方位.為方便討論，假設(shè)在飛輪作用下，各衛(wèi)星的線圈安裝坐標(biāo)系與主星參考坐標(biāo)系指向重合.在參考坐標(biāo)系下記Rij=［x y z］T.由式(9)可得:

由式(14)、(15)及式(2)可得到在上述假設(shè)條件下的電磁力表達(dá)式為

利用泰勒定理對(duì)Fmij在理想點(diǎn)

展開得

由式(19)可知，根據(jù)某一理想點(diǎn)S0可求出任意S 處的電磁力，電磁力簡化為S 的線性函數(shù).但這一表述基于泰勒展開定理，舍去二階以上表達(dá)式，只在S-S0較小時(shí)精確成立.因此為達(dá)到精度要求，實(shí)際應(yīng)用中需要實(shí)時(shí)校正理想點(diǎn)S0.

下面推導(dǎo)n 顆衛(wèi)星相互作用時(shí)電磁力的線性化.先假定有多顆主星，表述時(shí)以衛(wèi)星i為主星，按上述方法可表述如下:

記

為線性化的常數(shù)項(xiàng).

由式(21)可知，多顆衛(wèi)星的電磁編隊(duì)飛行中電磁力的線性化表述還有坐標(biāo)轉(zhuǎn)換矩陣的影響，增加了電磁力表達(dá)式的耦合性.這是由電磁力的作用機(jī)理決定的，給實(shí)際應(yīng)用帶來困難.對(duì)于多星電磁編隊(duì)飛行的電磁力線性化，本文只做初步的分析總結(jié).

3.2 氣動(dòng)力線性化

假設(shè)n 顆衛(wèi)星氣動(dòng)參數(shù)的CD和密度ρ 均相同，有一定的面質(zhì)比偏差.記主星迎流面面積為A，衛(wèi)星i 迎流面面積為A+δA.則衛(wèi)星i 的氣動(dòng)阻力在主星參考坐標(biāo)系下表述為

分別對(duì)式(22)中的［˙x ˙y ˙z］T、［x y z］T、δA求梯度，可導(dǎo)出

其中:

3.3 J2 攝動(dòng)線性化

結(jié)合本文的坐標(biāo)系定義，對(duì)參考文獻(xiàn)［7］的公式略作修改，可得

式中:

3.4 相對(duì)軌道動(dòng)力學(xué)方程最終表述

略去未知干擾項(xiàng)，由式(10)～(25)可得衛(wèi)星i 在主星參考坐標(biāo)系下的線性化相對(duì)軌道動(dòng)力學(xué)方程為

其中:

針對(duì)兩顆衛(wèi)星的電磁編隊(duì)飛行，略去未知干擾因素，公式可簡化為

式(26)是相對(duì)位置［x y z］T及控制量［IixIiyIizIjxIjyIjz］T的線性微分方程組，對(duì)控制算法設(shè)計(jì)及其他理論分析帶來便利.

經(jīng)分析，式(26)還可整理為狀態(tài)空間的表達(dá)形式.選狀態(tài)變量為控制量I =［IixIiyIizIjxIjyIjz］T.式(26)的狀態(tài)空間描述如下:

其中:

式中:E 為3×3 階單位矩陣;M1M2為矩陣-2M/m 的子式，且-2M/m=［M1M2］，M1為3×3 階矩陣，M2為3×6 階矩陣;F 為常數(shù)向量.

4 仿真及結(jié)果分析

為考核線性化相對(duì)動(dòng)力學(xué)模型的可用性及精度，對(duì)線性化相對(duì)軌道動(dòng)力學(xué)模型與未線性化軌道動(dòng)力學(xué)模型分別進(jìn)行仿真分析.仿真中設(shè)定主星與從星各軸線圈匝數(shù)、所圍面積相同.大氣密度查表可得.表1 給出初始時(shí)刻兩顆衛(wèi)星的軌道參數(shù).兩顆衛(wèi)星的其他設(shè)計(jì)參數(shù)見表2.軌道參數(shù)的選擇以構(gòu)成串行編隊(duì)為約束，衛(wèi)星間距約為20 m.不失一般性，從星的三軸電流均設(shè)為I =I0+δI=［2000+100sin(πt)］A，主星的三軸電流均設(shè)為I0=2 000 A.

表1 衛(wèi)星軌道參數(shù)

表2 衛(wèi)星其他設(shè)計(jì)參數(shù)

在上述條件下對(duì)線性化與非線性的動(dòng)力學(xué)模型進(jìn)行仿真，電磁力線性模型與非線性模型的變化曲線見圖2.線性模型與非線性模型的運(yùn)動(dòng)情況見圖3.

圖2 電磁編隊(duì)線性模型與非線性模型的電磁力對(duì)比

圖3 電磁編隊(duì)線性模型與非線性模型精度對(duì)比

由圖2 可見，在20 000 s 的仿真時(shí)間內(nèi)，電磁力線性模型與非線性模型計(jì)算值幾乎一致，誤差不超過10%，說明電磁力線性化是合理的.

由圖3 可見，在20 000 s 的仿真時(shí)間(衛(wèi)星軌道周期約5 900 s，即衛(wèi)星運(yùn)行3 個(gè)周期多)內(nèi)，線性模型與非線性模型運(yùn)行軌跡幾乎重合，線性模型與非線性模型的主從星相對(duì)運(yùn)動(dòng)在主星參考坐標(biāo)系的y、z 兩軸誤差不會(huì)發(fā)散，x 軸誤差逐漸發(fā)散，在仿真結(jié)束時(shí)有0.3 m 的誤差.單獨(dú)分析電磁力、J2攝動(dòng)等每一項(xiàng)單獨(dú)作用力線性化對(duì)模型精度的影響.仿真結(jié)束時(shí)，電磁力線性化可使模型誤差達(dá)0.2 m，其他作用(氣動(dòng)阻力線性化、J2攝動(dòng)線性化及Hill 方程誤差)對(duì)模型精度影響不超過0.1 m.

為考察不同軌道高度與不同衛(wèi)星質(zhì)量對(duì)模型精度的影響，設(shè)置軌道高度分別為6 600 km、6 882 km、7 400 km，衛(wèi)星質(zhì)量分別為100 kg、300 kg、500 kg，做9 組仿真實(shí)驗(yàn).鑒于x 軸誤差為主要誤差，圖4 給出各種情況下的線性化模型與非線性模型的x 軸誤差.由圖4 可見，誤差隨衛(wèi)星質(zhì)量的減小而增大，這是因?yàn)樵谝欢ǖ碾姶帕€性化誤差下，小質(zhì)量衛(wèi)星受到的電磁加速度誤差更大.減小線圈電流可以減小輕質(zhì)量引起的誤差.圖4 還表明軌道越高，誤差越小.分析認(rèn)為高軌道衛(wèi)星受大氣阻力、J2攝動(dòng)的影響較小，電磁力對(duì)相對(duì)位置的影響變大，使得模型對(duì)電磁力線性化的誤差更為敏感.

圖4 不同軌道高度與衛(wèi)星質(zhì)量情況下的x 軸誤差

仿真結(jié)果還表明這9 組實(shí)驗(yàn)中x 軸誤差主要由電磁力線性化產(chǎn)生，但篇幅有限，沒有列出仿真圖.

綜上所述，本文建立的電磁編隊(duì)線性化軌道動(dòng)力學(xué)模型在低軌道(200 ～1 000 km)、小衛(wèi)星(100 ～500 kg)、編隊(duì)距離20 m 左右、激勵(lì)電流幾千安的條件下是合理的，仿真結(jié)束時(shí)誤差在0.05 ～5 m 范圍內(nèi).且誤差主要來自于電磁力線性化.為電磁編隊(duì)飛行動(dòng)力學(xué)的進(jìn)一步研究提供參考，也可作為控制研究的理論基礎(chǔ).

5 結(jié) 論

本文以基于磁偶極子假設(shè)的遠(yuǎn)場模型為基礎(chǔ)，在考慮大氣攝動(dòng)、地球扁率J2攝動(dòng)的條件下推導(dǎo)了電磁編隊(duì)飛行的非線性相對(duì)軌道動(dòng)力學(xué)模型，并在一定條件下，對(duì)動(dòng)力學(xué)方程線性化，最終得到基于Hill 方程的線性化相對(duì)軌道動(dòng)力學(xué)方程.利用Matlab/Simulink 對(duì)電磁編隊(duì)線性模型與非線性模型進(jìn)行建模仿真，結(jié)果表明本文建立的線性化相對(duì)軌道動(dòng)力學(xué)與精確的非線性模型很好吻合，x 軸相對(duì)運(yùn)動(dòng)有累積誤差，誤差主要源于電磁力線性化，且誤差在一定的運(yùn)行時(shí)間內(nèi)很小.

由本文分析可知，線性化模型結(jié)構(gòu)清晰，形式簡單，便于獨(dú)立分析各控制變量作用，有利于更成熟的線性控制方法的應(yīng)用.非線性模型精度高，但結(jié)構(gòu)復(fù)雜，控制算法設(shè)計(jì)上有困難.本文的線性化模型合理可行，結(jié)果將為電磁編隊(duì)飛行的控制研究奠定理論基礎(chǔ).

［1］SCHWEIGHART S A，SEDWICK R J.Propellantless formation flight operations using electromagnetic formation flight［C］//SpaceOps 2006 Conference.Rome，Italy:［s.n.］，2006:1-38.

［2］HASHIMOTO T，SAKAI S，NINOMIYA K，et al.Formation flight using super-conducting magnets［C］//International Symposium of Formation Flying Missions and Technologies.Toulouse，F(xiàn)rance:［s.n.］，2002:3-10.

［3］KANEDA R，YAZAKI F，SAKAI S，et al.The relative position control in formation flying satellites using superconducting magnets［C］//2nd International Symposium on Formation Flying Missions and Technologies.Washington，DC:［s.n.］，2004:3-5.

［4］ELIAS L M.Dynamics of multi-body space interferometers including reaction wheel gyroscopic stiffening effects:structurally connected and electromagnetic formation flying architectures［D］.Combridge，MA:Massachusetts Institute of Technology，2004:1-23.

［5］ELIAS L M，KWON D W，SEDWICK R J，et al.Electromagnetic formation flight dynamics including reaction wheel gyroscopic stiffening effects［J］.Journal of Guidance，Control，and Dynamics，2007，30(2):499-511.

［6］AHSUN U，MILLER D W.Dynamics and control of electromagnetic satellite formations［D］.Cambridge，MA:Massachusetts Institute of Technology，2007:1-58.

［7］SCHWEIGHART S A.Development and analysis of a high fidelity linearized J2model for satellite formation flying［D］.Cambridge，MA:Massachusetts Institute of Technology，2001.