洪振強(qiáng)，俞 潔，劉 偉，楊立峰，陸國(guó)平，施晨康，趙 輝

（1.上海衛(wèi)星工程研究所，上海 201109；2.上海航天技術(shù)研究院，上海 201109）

0 引言

地球同步軌道（Geosynchronous Earth Orbit，GEO）三軸穩(wěn)定衛(wèi)星在軌穩(wěn)態(tài)輪控模式下，通過(guò)飛輪吸收環(huán)境干擾力矩，保持衛(wèi)星姿態(tài)穩(wěn)定。在環(huán)境力矩的作用下，飛輪角動(dòng)量存在隨時(shí)間累積項(xiàng)，需要定期進(jìn)行角動(dòng)量卸載，避免飛輪轉(zhuǎn)速飽和。目前，國(guó)內(nèi)外針對(duì)航天器環(huán)境干擾和動(dòng)力學(xué)特性的在軌辨識(shí)已經(jīng)開(kāi)展了廣泛研究。文獻(xiàn)［1］提出了一種利用光壓力矩輔助衛(wèi)星太陽(yáng)電池翼角度調(diào)整進(jìn)行角動(dòng)量管控的方法，但該卸載方法僅對(duì)對(duì)日定向衛(wèi)星適用。文獻(xiàn)［2］設(shè)計(jì)了一種基于以飛輪極限角動(dòng)量為參考的零運(yùn)動(dòng)力矩分配輪系角動(dòng)量管控策略，主要目的是避免系統(tǒng)角動(dòng)量未達(dá)到包絡(luò)面時(shí)部分飛輪角動(dòng)量出現(xiàn)飽和。文獻(xiàn)［3］的研究對(duì)象限定在采用慣性系為控制基準(zhǔn)的航天器角動(dòng)量管理。文獻(xiàn)［4］主要解決引力梯度力矩和氣動(dòng)力矩常值部分帶來(lái)的角動(dòng)量積累問(wèn)題。文獻(xiàn)［5］采用引力梯度力矩平衡姿態(tài)，設(shè)計(jì)了基于極點(diǎn)配置的空間站角動(dòng)量管理控制器。文獻(xiàn)［6］對(duì)環(huán)月地軌道環(huán)繞衛(wèi)星所受重力梯度力矩進(jìn)行了分析。文獻(xiàn)［7］研究了當(dāng)軌控推力器存在大力矩?cái)_動(dòng)時(shí)系統(tǒng)角動(dòng)量管控方法。文獻(xiàn)［8］研究了角動(dòng)量管理對(duì)衛(wèi)星軌道的影響，但均不涉及太陽(yáng)光壓力矩。文獻(xiàn)［9］研究電推進(jìn)衛(wèi)星通過(guò)將推力器的指向略微偏離質(zhì)心來(lái)產(chǎn)生控制力矩，完成角動(dòng)量卸載。文獻(xiàn)［10］設(shè)計(jì)了一種基于氣動(dòng)力矩和模型預(yù)測(cè)的低地軌道（Low Earth Orbit，LEO）衛(wèi)星角動(dòng)量管控算法，但不適用于GEO 衛(wèi)星。文獻(xiàn)［11］用幾何方法詳細(xì)討論了航天器各表面受曬的情況，推導(dǎo)了不同太陽(yáng)光入射條件下航天器各面所受到的太陽(yáng)光壓力矩?cái)?shù)學(xué)模型。文獻(xiàn)［12］提出一種精確求取太陽(yáng)光壓有效作用面積的方法。文獻(xiàn)［13］對(duì)太陽(yáng)光壓力矩與大氣阻力力矩的建模與仿真提出了新的方法，準(zhǔn)確地對(duì)太陽(yáng)光壓力矩與大氣阻力力矩進(jìn)行了估計(jì)。文獻(xiàn)［14］論述了大型三軸氣浮臺(tái)轉(zhuǎn)動(dòng)慣量和干擾力矩高精度聯(lián)合辨識(shí)技術(shù)，但需要進(jìn)行主動(dòng)激勵(lì)，不便于在軌操作。文獻(xiàn)［15］闡述了導(dǎo)航衛(wèi)星太陽(yáng)光壓建模方法與模型特性分析，主要研究太陽(yáng)光壓對(duì)軌道攝動(dòng)的影響，而不是姿態(tài)控制的影響。本文針對(duì)對(duì)稱布局GEO衛(wèi)星角動(dòng)量緩慢累積的特點(diǎn)，開(kāi)展在軌環(huán)境力矩辨識(shí)和角動(dòng)量管控研究，可為姿控系統(tǒng)的精細(xì)化設(shè)計(jì)提供參考依據(jù)，并以此為基礎(chǔ)提出使用3 臺(tái)飛輪接入閉環(huán)控制狀態(tài)下的角動(dòng)量管控策略，延長(zhǎng)角動(dòng)量卸載周期。

1 慣性系下角動(dòng)量模型

靜止軌道衛(wèi)星運(yùn)行在穩(wěn)態(tài)輪控模式時(shí)，通過(guò)飛輪吸收外界干擾力矩而保持衛(wèi)星姿態(tài)三軸穩(wěn)定。由于地球靜止軌道高度達(dá)到36 000 km，地磁力矩和大氣阻尼力矩可以忽略不計(jì)，只需要考慮太陽(yáng)光壓和重力梯度力矩的影響。

1.1 太陽(yáng)光壓干擾力矩

關(guān)于太陽(yáng)光壓力矩，已經(jīng)在很多文獻(xiàn)中進(jìn)行了詳細(xì)的推導(dǎo)和描述。本文的目的是在慣性坐標(biāo)系中建立太陽(yáng)光壓干擾模型，提取特征參數(shù)，為在軌精確辨識(shí)奠定基礎(chǔ)。單邊太陽(yáng)翼與雙邊對(duì)稱太陽(yáng)翼的太陽(yáng)光壓力矩分析方法一致，其最大區(qū)別是太陽(yáng)光壓壓心和衛(wèi)星質(zhì)心距離存在較大差別。

以雙邊對(duì)稱太陽(yáng)翼為例建立太陽(yáng)光壓干擾力矩模型，如圖1 所示。由于角動(dòng)量守恒定理須在慣性系下描述，因此，首先建立慣性參考系

O

X

Y

Z

，其原點(diǎn)

O

位于衛(wèi)星質(zhì)心，

O

X

指向衛(wèi)星星下點(diǎn)時(shí)刻為12：00 對(duì)應(yīng)的飛行方向，

O

Y

沿軌道負(fù)法向，

O

Z

、

O

X

和

O

Y

滿足右手定則。設(shè)太陽(yáng)光壓壓心在

O

X

Y

Z

坐標(biāo)系中的初始位置矢量為

P

=[

P P P

]，太陽(yáng)光壓法向力為

F

=[0 0

F

]，切向力為

F

=[0

F

0]，力的大小與太陽(yáng)光入射角

α

的關(guān)系見(jiàn)文獻(xiàn)［1］。

圖1 太陽(yáng)光壓作用示意圖Fig.1 Schematic diagram of solar pressure effects

設(shè)慣性系中太陽(yáng)光壓力矩為

T

，其可根據(jù)太陽(yáng)光壓力和壓心矢量計(jì)算得到。由于壓心矢量隨衛(wèi)星的軌道運(yùn)動(dòng)相對(duì)慣性空間旋轉(zhuǎn)，所以先計(jì)算時(shí)變的壓心矢量：

式中：

A

(

θ

)為繞

y

軸的旋轉(zhuǎn)矩陣；旋轉(zhuǎn)角

θ

與軌道角速度相關(guān)。對(duì)于GEO 衛(wèi)星，可直接由軌道角速度

ω

乘以時(shí)間得到。進(jìn)一步可得慣性系下時(shí)變的太陽(yáng)光壓力矩

T

(

t

)為

1.2 重力梯度干擾力矩

對(duì)于圓軌道衛(wèi)星，在小姿態(tài)情況下，衛(wèi)星本體坐標(biāo)系下的重力梯度力矩

T

可表示為

式中：

φ

、

θ

分別為衛(wèi)星的滾動(dòng)角和俯仰角。從而可得慣性系中的重力梯度力矩

T

(

t

)為

穩(wěn)態(tài)輪控下，衛(wèi)星本體坐標(biāo)系與軌道坐標(biāo)系近似重合，即

φ

≈

θ

≈0，因此，式（4）可簡(jiǎn)化為

1.3 整星角動(dòng)量積累模型

根據(jù)角動(dòng)量守恒定理，在慣性系下對(duì)干擾力矩進(jìn)行積分即可得到慣性系中的角動(dòng)量變化曲線。設(shè)慣性系中的角動(dòng)量為

H

=[

H

H

H

]，有

進(jìn)一步，結(jié)合衛(wèi)星的軌道和姿態(tài)運(yùn)動(dòng)，可將慣性系中的角動(dòng)量映射到衛(wèi)星本體系中，得到衛(wèi)星本體系下的角動(dòng)量

H

=[

H

H

H

]，對(duì)于GEO衛(wèi)星，有

經(jīng)上述分析，可得GEO 軌道衛(wèi)星環(huán)境干擾力矩特性，見(jiàn)表1。

X

、

Z

軸合成角動(dòng)量反映了太陽(yáng)光壓干擾力矩的累積效應(yīng)，

Y

軸角動(dòng)量反映了重力梯度力矩的累積效應(yīng)。

表1 GEO 軌道環(huán)境干擾力矩分析總結(jié)表（相對(duì)軌道系零姿態(tài)時(shí)）Tab.l Analysis summary of GEO orbit environment disturbance torques（relative to the zero attitude of the orbit system）

2 環(huán)境干擾力矩辨識(shí)

2.1 辨識(shí)算法

根據(jù)慣性系下的角動(dòng)量模型式（7）和式（8），其中，式（7）反映了法向太陽(yáng)光壓引起的角動(dòng)量累積和切向太陽(yáng)光壓及重力梯度力矩引起的角動(dòng)量交變，式（8）反映了重力梯度力矩引起的角動(dòng)量累積和切向太陽(yáng)光壓引起的角動(dòng)量交變。提取特征參數(shù)，改寫(xiě)環(huán)境干擾角動(dòng)量累積模型如下：

式中：

H

、

H

分別為

X

向和

Y

向的初始角動(dòng)量；特征參數(shù)

a

、

a

、

c

、

a

、

a

、

c

分別為

根據(jù)上述模型，結(jié)合在軌遙測(cè)數(shù)據(jù)，采用最小二乘擬合即可得到各特征系數(shù)。根據(jù)特征參數(shù)，可計(jì)算得到法向和切向太陽(yáng)光壓力引起的干擾力矩辨識(shí)結(jié)果如下：

結(jié)果1，法向太陽(yáng)光壓力引起的繞

O

X

軸干擾力矩，

P

?

F

=

c

；結(jié)果2，法向太陽(yáng)光壓力引起的繞

O

Y

軸干擾力矩，

P

?

F

=-

a

?

ω

；結(jié)果3，切向太陽(yáng)光壓力引起的繞

O

Z

軸干擾力矩，

P

?

F

=-

a

?

ω

；

2.2 辨識(shí)結(jié)果

根據(jù)飛輪在軌轉(zhuǎn)速遙測(cè)數(shù)據(jù)，處理得到衛(wèi)星本體系中和慣性系中的角動(dòng)量變化曲線如圖2 所示。

圖2 衛(wèi)星在軌角動(dòng)量變化示意圖Fig.2 Angular momentum change of the on-orbit satellite

對(duì)慣性系下的角動(dòng)量曲線進(jìn)行擬合，得到擬合結(jié)果如圖3 和圖4 所示。

圖3 慣性系下整星角動(dòng)量擬合結(jié)果Fig.3 Fitting results of the angular momentum of the whole satellite in the inertial frame

圖4 慣性系下整星角動(dòng)量擬合殘差Fig.4 Fitting residual error of the whole satellite angular momentum in the inertial frame

對(duì)應(yīng)特征參數(shù)的數(shù)值辨識(shí)結(jié)果如下：

已知衛(wèi)星慣性角速度

ω

=7.272 2×10rad/s，根據(jù)上述辨識(shí)結(jié)果可解得，在星下點(diǎn)時(shí)刻為12：00 對(duì)應(yīng)的慣性參考系中：結(jié)果1，法向太陽(yáng)光壓力引起的繞

O

X

軸干擾力矩，

P

?

F

=7.91×10N ?m；結(jié)果2，法向太陽(yáng)光壓力引起的繞

O

Y

軸干擾力矩，

P

?

F

=3.54×10N ?m；結(jié)果3，切向太陽(yáng)光壓力引起的繞

O

Z

軸干擾力矩，

P

?

F

=2.40×10N ?m；結(jié)果4，切向太陽(yáng)光壓力引起的繞

O

X

軸干擾力矩，

P

?

F

=-2.24×10N ?m；結(jié)果5，繞

O

Y

軸重力梯度力矩，

T

=3

ω

?

I

=-2.20×10N ?m。

將上述辨識(shí)結(jié)果代入式（10），將仿真結(jié)果與在軌遙測(cè)結(jié)果比對(duì)，如圖5 所示，誤差在0.5 N·m·s 以內(nèi)，表明所建角動(dòng)量模型和特征參數(shù)辨識(shí)算法正確且有效。

圖5 角動(dòng)量仿真模型與在軌實(shí)際角動(dòng)量比對(duì)示意圖Fig.5 Comparison diagram of the angular momentum simulation model and the actual on-orbit angular momentum

3 角動(dòng)量管控

根據(jù)前文的角動(dòng)量積累機(jī)理分析可知，為了延長(zhǎng)GEO 衛(wèi)星在軌的角動(dòng)量卸載周期，一方面應(yīng)盡量減小衛(wèi)星質(zhì)心與太陽(yáng)光壓壓心之間的距離（尤其是在軌道法向）；另一方面應(yīng)盡量減小衛(wèi)星慣量積

I

，從而減小角動(dòng)量積累的梯度。但是，隨著季節(jié)的變化和星上剩余燃料的減少，衛(wèi)星的壓心、質(zhì)心、轉(zhuǎn)動(dòng)慣量和慣性積均不可避免發(fā)生變化。角動(dòng)量管控的目的是為了延長(zhǎng)卸載周期，根據(jù)角動(dòng)量增長(zhǎng)趨勢(shì)設(shè)置合理的角動(dòng)量初值（即對(duì)應(yīng)飛輪初始轉(zhuǎn)速）。因此，有必要對(duì)太陽(yáng)光壓力矩和重力梯度力矩進(jìn)行辨識(shí)，并根據(jù)辨識(shí)結(jié)果合理設(shè)置飛輪卸載轉(zhuǎn)速（而不是簡(jiǎn)單卸載到零轉(zhuǎn)速），從而延長(zhǎng)角動(dòng)量卸載周期，提升衛(wèi)星的在軌應(yīng)用效能。假設(shè)衛(wèi)星配置三正一斜總共4 臺(tái)飛輪，正裝

XYZ

飛輪與衛(wèi)星本體軸一致，斜裝

S

飛輪與

XYZ

夾角一致，均為54.74°。以上一節(jié)的太陽(yáng)光壓力矩和重力梯度力矩辨識(shí)結(jié)果為例，將其代入整星動(dòng)力學(xué)模型，以飛輪轉(zhuǎn)速超過(guò)±2 000 r/min（對(duì)應(yīng)角動(dòng)量約15 N·m·s）為卸載標(biāo)志，則可得不同飛輪組合接入狀態(tài)下，角動(dòng)量管控前（對(duì)應(yīng)飛輪初始轉(zhuǎn)速均為0 r/min）對(duì)應(yīng)的卸載周期和進(jìn)行角動(dòng)量管控后的卸載周期見(jiàn)表2。由表2 可見(jiàn)，角動(dòng)量管控后使得角動(dòng)量卸載周期增加1 倍，管控前后的整星角動(dòng)量累積曲線如圖6 和圖7 所示。圖7 中，從左到右4張圖分別展示了4 種接入工況下的飛輪轉(zhuǎn)速變化曲線，每張圖的最左側(cè)為設(shè)置的飛輪轉(zhuǎn)速初值，最右側(cè)為卸載時(shí)對(duì)應(yīng)的飛輪轉(zhuǎn)速，橫坐標(biāo)表示卸載周期。

圖6 角動(dòng)量管控前不同飛輪組合下的飛輪轉(zhuǎn)速變化仿真曲線Fig.6 Simulation curves of flywheel speed changes under different flywheel combinations before angular momentum management and control

圖7 角動(dòng)量管控后不同飛輪組合下的飛輪轉(zhuǎn)速變化仿真曲線Fig.7 Simulation curves of flywheel speed changes under different flywheel combinations after angular momentum management and control

表2 不同飛輪組合下角動(dòng)量管控前后卸載周期統(tǒng)計(jì)表Tab.2 Statistics of unloading periods before and after angular momentum management and control under different flywheel combinations

4 結(jié)束語(yǔ)

本文建立了慣性系下雙對(duì)稱太陽(yáng)電池陣GEO衛(wèi)星的角動(dòng)量積累模型，并投影到衛(wèi)星本體系中，將角動(dòng)量積累轉(zhuǎn)化為實(shí)際飛輪轉(zhuǎn)速的變化。結(jié)合實(shí)際在軌飛輪轉(zhuǎn)速遙測(cè)數(shù)據(jù)，精準(zhǔn)辨識(shí)GEO 衛(wèi)星在軌環(huán)境干擾力矩，獲取定量結(jié)果，并作為輸入條件，仿真了不同飛輪組合接入狀態(tài)下的飛輪轉(zhuǎn)速變化情況。以角動(dòng)量卸載周期最長(zhǎng)為原則，優(yōu)化了角動(dòng)量管控策略，將飛輪的角動(dòng)量卸載周期提升為原來(lái)的2 倍，大幅提升衛(wèi)星在軌的應(yīng)用效能。