基于優化的羽流撞擊消旋技術瞄準策略設計

，，

中國空間技術研究院 通信衛星事業部，北京 100094

在軌服務技術是空間應用技術領域的研究熱點。從國內外研究現狀來看，在軌服務的對象正逐漸由合作目標向非合作目標擴展。非合作目標泛指既不能提供有效合作信息，也不具備專門的在軌服務接口的航天器[1]，例如失效航天器和空間碎片。失效的GEO和MEO衛星所受的大氣和地磁阻尼微弱，在太陽光壓的長期作用下往往進入自由翻滾的運動狀態。這類目標給在軌捕獲帶來很大挑戰，因此國內外學者提出了多種消旋技術構想[2]，希望先為目標消旋再進行捕獲。消旋技術的實現方式又可以籠統分為接觸式與非接觸式，非接觸式消旋技術包括電磁力消旋和羽流撞擊消旋[3]等。

羽流撞擊消旋技術由于系統復雜度低、適用范圍廣、任務安全性高等優點而受到關注。羽流撞擊消旋的主要構想是在服務平臺部署化學推力器，令噴管指向目標衛星的太陽翼，利用羽流對目標星的撞擊力形成與目標星角動量方向相反的力矩，實現對目標星的姿態調整和控制，最終為目標星消旋。相比姿態同步交會對接和機械臂抓捕等技術方案，消旋任務模式下服務平臺可以避免進入目標星翻滾運動的包絡區域，具有較高的安全性，且目標星上不需配置專用結構或機構，降低了系統復雜度，具有良好的適用性。

國外學者從動力學與控制角度出發開展了相關研究。Ferrari和Lavagna提出了一種羽流撞擊模型和動力學仿真框架，研究了目標自旋運動情況下的消旋動力學過程，證明了技術的可實現性[4]；Thomas開展了羽流撞擊消旋技術的工程探索[5]，包括服務平臺構型以及軌道和姿態控制系統的設計，對化學推力器羽流模型進行了比較細致的研究，并設計了考慮服務平臺軌道和姿態控制的閉環仿真系統[6]。然而上述研究存在一定的局限性，首先是研究對象均選擇了ENVISAT衛星，ENVISAT的太陽翼面積較大且為矩形，利于減小羽流的非線性效應；其次ENVISAT的運動形式為繞次大慣性軸自旋，章動效應不明顯，太陽翼法線與當地線速度近似重合，對羽流撞擊消旋而言是很理想的工況；此外，Thomas定義瞄準策略(pointing strategy)為消旋推力器噴管相對目標應具備的位置和姿態，基于噴管軸線方向與目標的角速度矢量的幾何關系給出了一種瞄準策略，該方法精度有限且受目標自身特點影響較大，不利于推廣。

因此，本文重點研究以太陽翼為主要受力對象的羽流撞擊消旋問題，提出了基于優化方法的羽流撞擊消旋瞄準策略設計方法，重點解決一般情形下的瞄準策略。首先針對常見GEO衛星平臺和太陽翼構型對消旋問題進行了歸納并給出了典型工況，其次根據參數化羽流撞擊模型分析羽流撞擊力相對目標位置姿態的敏感性，在此基礎上提出了基于優化方法的瞄準策略規劃模型，仿真結果表明該瞄準策略不但能顯著提高消旋力矩精度，而且對不同的消旋工況具有良好的適用性。

1 數學模型

1.1 消旋動力學模型

高軌航天器中GEO通信衛星占據較大比例，同時由于軌位資源限制，對在軌服務的需求尤其旺盛。一般的GEO通信衛星平臺采用中心體加附件的基本構型，附件包括對稱或非對稱安裝的太陽翼和各類天線設備等。天線類設備體積較小且剛度較大，在本文中不單獨考慮，簡化為中心剛體的一部分。在三軸穩定衛星的動力學分析中，太陽翼一般作為柔性附件考慮，由于其展開機構處通常配置有阻尼器，因此在本文中為了降低動力學建模的復雜度，假設太陽翼為剛體，不發生彈性變形。

自由翻滾衛星通常處于整星斷電狀態，但太陽翼驅動機構具有一定的斷電保持力矩，因此假設太陽翼保持與星本體相對轉角不變，不會因羽流撞擊產生的力矩轉動。

在此基礎上,建立衛星模型,并定義坐標系如下：

{b},目標本體坐標系，原點位于星體質心，zb沿星體縱向，xb軸為太陽翼安裝方向，yb軸由右手坐標系確定。

{A},附件坐標系，初始狀態下，原點位于星體質心，xA軸為太陽翼的縱向，zA軸沿太陽翼法線方向且與星本體zb軸重合，yA軸沿太陽翼橫向。

單軸式太陽翼可繞xA軸轉動，在初始安裝位置下，xA軸與xb軸重合。目標星基本構型以及坐標系關系如圖1所示。

圖1 目標星構型Fig.1 Configuration of target

根據最大慣性軸定理，非對稱剛體的自由運動趨于演變為繞最大慣量軸的自旋運動[7]。一些針對在軌失效航天器的研究也表明，失效衛星往往形成繞最大或者次大主慣性軸的近似自旋運動[6]。基于以上分析，針對本文選取的衛星構型，假設衛星近似自旋運動，角動量矢量指向星本體xb軸附近。剛體衛星的姿態動力學方程為:

(1)

式中：ω為衛星瞬時角速度；I為目標星慣量陣；T為繞質心的力矩。假定目標衛星本體系與慣性主軸重合，則在{b}系下I為對角矩陣。

根據文獻[6]，期望消旋力矩矢量可以按照如下方法確定：依據目標星姿態運動參數，選取適當的控制器來確定期望角加速度，從而根據式(1)得到期望消旋力矩矢量。針對于處于自旋狀態的目標星，考慮比例控制則期望角加速度為：

(2)

式中：Kp為比例系數；ωref為控制目標，在消旋問題中應令ωref=0。考慮太陽翼相對初始位置的轉角a，{b}到{A}的轉換矩陣為：

(3)

顯然，在{A}中研究消旋推力器相對太陽翼的位置姿態關系更加直觀簡便，可以直接建立期望消旋力矩與消旋推力器噴管理想位置姿態之間的對應關系，因此下文中如不加特殊說明，標量運算均在{A}中進行。

為了簡化問題便于分析計算，本文忽略太陽翼的表面特征，認為太陽翼為連續的無厚度矩形平板。在消旋操作過程中，總是對某一側太陽翼施加羽流撞擊作用，因此不妨以+xA一側的太陽翼為研究對象，尺寸參數如表1所示。

表1 太陽翼幾何參數

1.2 羽流場模型

Simons模型[8]是一種參數化羽流分布模型。模型指出，羽流傳播距離遠大于噴管出口截面尺寸時，羽流場中任一點處的壓強可以用流場的幾何參數表示。基于分析建立的羽流密度分布函數為：

(4)

式中：ρ*為噴管喉部的壓強；ρ為流場中任一點處的羽流密度大小；r為該點距離噴管出口的距離標量；R*為噴管喉部半徑；f(θ)為羽流角分布函數，θ為沿噴管出口到該點視線方向與推力器軸線方向的夾角；羽流常數Ap由氣體性質所決定，

(5)

式中：θlim為推力器羽流的極限偏轉角；U*/Ulim為氣體極限膨脹狀態下的當地馬赫數。根據等熵流動特性，得到近似關系：

(6)

式中：γ為氣體比熱比。

羽流的角分布函數f(θ)為：

(7)

式中：θ∞為超音速無粘流時的θlim值。Simons考慮噴管出口處附面層效應得到羽流隨偏轉角度增大的衰減滿足指數規律[7]：

f(θ)=f(θ0)e-β(θ-θ0)

(8)

式中：θ0和β的取值取決于噴管出口處的狀態。如果流動是理想的無粘流，附面層中的流動也符合余弦規律：

(9)

式中：δ為噴管出口處的附面層厚度；Re為噴管出口處的截面半徑。當θ0接近θ∞時，有以下近似關系成立：

(10)

由此可得β：

(11)

1.3 表面作用模型

羽流充分膨脹后形成自由分子流動，在與物面接觸的過程中通過碰撞、反彈和吸附等行為將動量傳遞給目標。Feshe模型對于羽流動量傳播采用的假設為：沿噴管出口軸線方向上微元受力應該等于噴管出口處的壓力。文獻[6]對羽流模型進行了分析比較，將Simons模型的分布規律代入Feshe模型：

(12)

式中：C為常數，

(13)

羽流氣體與物體表面的相互作用是復雜的物理過程，通常用鏡面反射系數和散射系數來描述傳播的動量[9]：

(14)

式中：Uw為完全熱適應情況下分子離開壁面時的反射速度；v為入射角；n為當地法向；ev為來流方向;σn和στ分別為法向和切向的動量傳播系數，根據物體表面特性和羽流特性確定。

本文中，根據LEO軌道上試驗數據[10]給出以下假設：

1)化學推力器的羽流是高熱自由分子流動；

2)壁面溫度300 K，得到速度比為0.23；

3)氣體采用完全熱適應假設，存在97%漫反射。

針對消旋推力器參數的選取，本文以某單組元10 N推力器作為參考，給出設計參數如表2所示。

表2 消旋推力器設計參數

2 消旋策略

2.1 敏感性分析方法

除了在軌服務技術中普遍存在的技術難點[11-13]外，羽流撞擊消旋技術的復雜性主要體現在羽流模型的非線性。在消旋問題中這種非線性體現在兩個方面：第一，隨著入射流的偏轉以及目標表面的特征，羽流撞擊所產生合力的大小和方向存在非線性；第二，由于羽流的擴散效應以及目標尺寸有限，導致總有部分羽流不能被目標捕獲而造成非線性。羽流敏感性分析是研究消旋力矩相對噴管相對位置姿態的變化趨勢，對非線性進行初步分析，從而為瞄準策略設計提供基礎。

根據式(14)給出的流場中某一面積微元所受羽流撞擊力，在附件系下對太陽翼表面區域積分得到太陽翼所受的羽流撞擊力為：

(15)

表面力作用的積分屬于第一類曲面積分，可以轉化為每個分量在該區域上的二重積分，即：

(16)

類似有羽流撞擊力所產生的相對目標星質心的力矩為:

(17)

因此消旋力矩可以表示為噴管相對位置姿態的函數：

T=T(xN,yN,zN,α,β)

(18)

式中：(xN,yN,zN)為消旋推力器噴管的坐標；(α,β)為推力器噴管指向的方位角，定義方式與極坐標相同，偏航角α∈[0,2π]表示噴管指向矢量與xA軸正向夾角，俯仰角β∈[0,π]表示噴管指向矢量與zA軸正向夾角。

由式(18)可知消旋力矩矢量由5個變量確定，對應噴管相對太陽翼的5自由度運動。因此在敏感性分析中分別考慮推力器典型的平動和轉動對消旋力矩的影響。

平動情況下令噴管垂直指向太陽翼，即β=π，噴管沿3個坐標軸平動，從而分析傳播距離和捕獲面積對消旋力矩的影響。轉動情況下，選取太陽翼中心區域某一定點作為瞄準點，改變噴管軸線相對指向，從而研究在不同入射方位下的消旋力矩變化規律。

2.2 瞄準優化方法

(1)約束條件

在整個消旋任務中，服務星應該滿足安全性約束，本文的安全約束區域按照球形安全區域設計，給出極限安全距離約束條件為:

(19)

式中：RN為在附件坐標系下的位置矢量；dmin為最小安全半徑。

(2)目標函數

根據分析，瞄準策略的設計可用優化問題描述，優化目標是使消旋力矩的方向盡量接近期望方向的同時使幅值最大化。根據指數函數特點不妨令目標函數為：

(20)

式中：B為常數;‖T‖表示消旋力矩矢量的范數，

(21)

ε為實際力矩與期望力矩方向的夾角，即誤差角:

(22)

該目標函數具備以下性質：

1)消旋力矩幅值‖T‖增大以及誤差夾角ε均使得目標函數降低。

2)在‖T‖增大的初段，目標函數下降迅速，此時力矩幅值所占的權重較大；在‖T‖增長的末段，指數項和三角函數項的變化速度相近，此時精度與幅值具有相近的權重。

3)調控參數B可以改善函數特性，對優化結果有一定影響。

結合約束條件得到瞄準策略的優化模型為:

(23)

3 仿真分析

3.1 敏感性仿真

根據第2.1節中的分析，可以設計仿真如下：僅平動情況為噴管軸線始終指向-z軸方向，控制變量β=π。消旋力矩相對噴管坐標的變化規律如圖2所示。

圖2(a)為噴管坐標y=0時，不同z向距離下，Ty相對x坐標的變化規律，隨著羽流傳播距離越短，Ty極值點對應的x坐標逐漸增大。由于太陽翼關于x軸的對稱性，Tx和Tz的大小不受x改變影響。同時，Ty最大值與z的大小近似成平方反比。

圖2(b)(d)分別展示了在不同z向距離下，Tx和Tz相對y坐標的變化規律。其中Tx和Tz均隨y增大而增大，羽流撞擊的切向力所產生的力矩要大于等效力臂增大產生的力矩，從而說明針對當前構型的太陽翼，羽流撞擊消旋能夠很好地適應太陽翼傾斜的工況。

圖2(c)表示Ty相對y坐標的變化規律。Ty隨著y增大而減小，可以認為是由于太陽翼捕獲的羽流減少了。

此外，可以看出Tx量級小于Ty和Tz，這是由于當噴管指向位置沿+y偏移時，切向合力產生+x方向的力矩，法向合力產生-x方向的力矩，兩者存在一定的抵消效應，這種效應使得繞太陽翼x軸方向的轉動相比繞另外兩軸的轉動運動更難消除。另一方面也可以認為羽流撞擊的方式產生的Tx較小，能夠滿足太陽翼不發生轉動的基本假設。

其次分析消旋力矩相對入射角的變化規律。令噴管始終指向點(9.8,0,1.5)，傳播距離為8 m，在不同的噴管指向角下，根據噴管相對位置姿態關系，可以反解出噴管的坐標，從而可以分別計算噴管沿不同指向時消旋力矩各分量的變化規律。利用對稱性可以降低計算量，當β∈[90°,180°]時，表示噴管垂直向下，力矩變化規律與β∈[0,90°]時對稱；同理，在a∈[0°,180°]時的力矩分布規律與a∈[180°,270°]對稱，因此只考慮a∈[0°,180°]，β∈[120°,180°]時力矩變化規律。

圖3展示了3個方向的力矩分量相對噴管方位角的變化規律。從仿真結果可以看出，推力器越接近于垂直方向，即俯仰角β接近180°時，Ty分量到達最大，同時合力矩最大。當噴管在yoz平面附近以較大入射角施加羽流時，偏航角α在90°或270°附近，而俯仰角β較小，此時對應的羽流捕獲面積最小的情況，合力矩明顯減小，由于y方向的羽流撞擊力增大，此時Tz有最大值，意味著改變噴管傾角能夠產生相當的切向力從而產生較大的z向消旋力矩。此外，在偏轉入射時Tx的量級同樣小于Ty和Tz，說明傾斜入射工況下Tx的抵消作用仍舊存在，與之前的分析一致。

圖2 平動時消旋力矩變化規律Fig.2 Torque relates to translational motions

圖3 噴管轉動對消旋力矩影響Fig.3 Influence of thruster rotation on torques

綜合以上分析可以得到如下結論：羽流撞擊所產生的消旋力矩能夠產生足夠大的Ty分量，是消旋實現的核心要求；噴管傾斜時能夠產生同一量級的Ty和Tz，較難產生相同量級的Tx，這一特點符合羽流撞擊消旋的任務需求，同時對于任務規劃具有指導作用。在羽流撞擊消旋的過程中，應當盡量保持噴管軸線指向目標平板的橫向中心區域以減輕非線性效應。

3.2 消旋仿真

由于在實際消旋過程中，化學推力器羽流的傳播速度遠大于服務平臺與目標之間的相對運動速度，因此忽略羽流場的瞬態效應以及雙星相對運動造成的相對位置姿態改變。同時，羽流撞擊消旋總是針對一側的太陽翼開展的。隨著目標周期性轉動運動，總以逐漸接近的太陽帆板為主要受力對象，因此假設目標其他部分受到的羽流撞擊忽略不計。

根據式(23)所示的瞄準策略規劃模型，設計仿真算例如表3所示。情況1為期望消旋力矩方向Tdes在y軸附近，情況2和3為z方向期望消旋力矩分量較大；同時為考察目標函數參數不同所造成的影響，令B取值為0.01和0.5。

表3 仿真參數設置

仿真計算的結果如表4和圖4所示。對于期望力矩方向接近y軸正向的情況，該優化方法能夠使實際消旋力矩與期望方向夾角收斂到0°，同時能保持較大的力矩幅值，這對于消旋是非常有利的。對于Tz較大的情況，調整噴管相對位置姿態能夠產生足夠的z方向消旋力矩。此時可能出現問題是，力矩的幅值和指向精度需要權衡，如圖4(b)所示，追求較小的誤差角度可能導致力矩幅值減小，此時可以通過合理配置參數A進行修正。

從仿真結果中可以看出A對精度和消旋力矩幅值有明顯的調整作用，因此有必要分析參數A的影響。保持期望消旋力矩矢量為(-0.01,1,0.4)，令A在[0.01,10]的范圍內變化，得到參數A對消旋力矩各分量以及模值大小Tabc的變化規律如圖5所示。參數A的取值對消旋力矩有明顯影響，隨著A取值的減小，消旋力矩的方向誤差和消旋力矩的幅值均減小，相當于增加了角度誤差的權重，使得尋優結果傾向于角度誤差更小。當A逐步增大，當力矩幅值收斂時，角度誤差約為12°。因此，可以認為A的理想取值范圍應大約在[1,10]，在保證力矩較大的情況下有足夠精度。

表4 仿真計算結果

圖4 仿真計算結果Fig.4 Results of simulation

圖5 參數B不同取值對結果的影響Fig.5 Effect of various B to results

4 結束語

本文基于參數化羽流撞擊模型，針對中心剛體加太陽翼形式的典型目標，開展了羽流撞擊消旋的瞄準策略的優化設計，主要結論如下：

1)首先分析了消旋力矩對相對位置姿態的敏感性，可以看出羽流撞擊能夠產生同量級的法向和切向力，因此適用于受力表面存在一定傾斜的情況，針對細長形狀的太陽翼，受限于形狀很難產生較大的繞長軸方向的消旋力矩。

2)利用指數和余弦函數特性建立的目標函數可以實現瞄準策略的優化設計，給出理想的推力器位置姿態，從而為制導策略的設計提供參考，同時通過調整目標函數中的參數，能夠對精度和幅值權重合理分配，改善了優化過程和實際效果。

3)基于優化搜索方法實現的瞄準策略提高了羽流撞擊消旋的瞄準精度，并給出了一系列仿真結果作為參考，為后續在線瞄準策略設計提供了基礎，但由于優化方法需要約束條件才能保證收斂性，且計算時間較長，因此仍需對算法進一步改進，提升計算效率和易用度。