2016, Vol. 3
近年来，小行星探测逐渐成为人类深空探测的又一热点[1, 2, 3]。开展对小行星的探测，可为研究太阳系的形成与演化、生命的起源与进化、防御小行星的撞击等重大基础科学问题提供原始资料[4, 5]。而小行星的自转给小行星探测、登陆小行星和防御小行星等带来了很多麻烦。为便于很好地跟踪探测小行星，学者们纷纷把注意力集中于近小行星航天器的轨道动力学问题，对不规则的小行星引力场进行建模[6, 7, 8]，分析复杂引力场中航天器的轨道特点[9, 10, 11]。
而小行星的自转给登陆和防御小行星带来的困难仍未完全解决。为此，本文提出一种利用太阳帆绳系系统逐渐减小小行星自转速率的方法，其构型如图 1所示。其中，太阳帆的初始位置A位于小行星的同步轨道上，系绳与小行星的连接点B固定在偏离平衡点正下方一段距离的位置。在系绳长度不变的情况下，利用太阳帆受到的太阳光压力使其始终保持与小行星同步，这样就避免了由于小行星自转而引起的系绳缠绕问题。通过控制太阳帆使系绳始终拉紧，系绳中的拉力便可以持续提供一个与小行星自转方向相反的力矩，从而减小其自转速率。
|
图 1 太阳帆绳系系统构型 Fig. 1 Configuration of the tethered solar sail system |
本中首先建立小行星固连坐标系下太阳帆的动力学模型。其次，选定一颗具有代表性的小行星，利用其质量、自转角速率、近似椭球体的三轴长度等参数计算得到二阶二次引力场中稳定平衡点的初始位置。最后，仿真得到使该小行星停止旋转过程中太阳帆姿态的控制律以及所需的时间，验证了方案的可行性。
1 动力学模型假设系绳与太阳帆始终在小行星的赤道平面(Oxy平面)内，且太阳帆受到的太阳光压力在自转轴方向上的拉力为零。由于系绳拉力的影响，小行星的自转角速度逐渐减小。本文对文献[6]中航天器在小行星2阶2次引力场中的动力学模型进行了改进，则太阳帆在小行星2阶2次引力场中的动力学模型为
| $ \left\{ {\begin{aligned} \ddot x = 2{ω} \dot y + {{ω} ^2}x - \frac{{{μ} x}}{{{r^3}}} + \frac{{3{μ} {C_{20}}x}}{{2{r^5}}} + \frac{{6{μ} {C_{22}}x}}{{{r^5}}} \hspace{10pt}\\ - \frac{{15{μ} {C_{22}}x({x^2} - {y^2})}}{{{r^7}}} + \frac{{{F}\cos {φ} - {T}\sin \theta }}{m}{\text{ }} \\ \ddot y = - 2{ω} \dot x + {{ω} ^2}y - \frac{{{μ} y}}{{{r^3}}} + \frac{{3{μ} {C_{20}}y}}{{2{r^5}}} - \frac{{6{μ} {C_{22}}y}}{{{r^5}}} \hspace{5pt}\\ - \frac{{15{μ} {C_{22}}y({x^2} - {y^2})}}{{{r^7}}} + \frac{{{F}\sin {φ} - {T}\cos \theta }}{m} \hspace{6pt}\\ \dot {ω} = - \frac{{ {{T}}y\sin \theta }}{J} \hspace{125pt} \end{aligned}} \right. $ | (1) |
其中：ω为小行星的自转角速度； $r = \sqrt {{x^2} + {y^2}}$ ；μ为小行星的引力常数；C20，C22分别为小行星的球谐系数；F为太阳帆受到的太阳光压力大小；φ为F与x轴的夹角；T为系绳中的拉力大小；θ为T与y轴的夹角；m 为太阳帆的质量；J为小行星的转动惯量。
| $ \ddot x\cos (\frac{{π} }{2} - \theta ) + \ddot y\cos \theta = 0 $ | (2) |
由于系绳的长度不变，所以太阳帆在x轴和y轴方向上的加速度在系绳方向上的分量和为零，即
由此可以得到系绳中拉力T的表达式。
由于太阳光压力的大小与太阳帆到太阳的距离有关，所以先将太阳帆在小行星固连坐标系中的位置坐标转换为日心惯性系中的位置坐标
| $ \left[{\begin{array}{*{20}{c}} \ {{X_A}} \\ {{Y_A}} \\ {{Z_A}} \end{array}} \right] = \left[{\begin{array}{*{20}{c}} 1& 0& 0 \\ 0& {\cos {λ} }& {\sin {λ} } \\ 0& { - \sin {λ} }& {\cos {λ} } \end{array}} \right]\\ \left[{\begin{array}{*{20}{c}} {\cos {ω} t}& { - \sin {ω} t}& 0 \\ {\sin {ω} t}& {\cos {ω} t}& 0 \\ 0& 0& 1 \end{array}} \right]\left[{\begin{array}{*{20}{c}} x \\ y \\ 0 \end{array}} \right] + \left[{\begin{array}{*{20}{c}} {{X_O}} \\ {{Y_O}} \\ {{Z_O}} \end{array}} \right] $ | (3) |
其中：[XA YA ZA]T为太阳帆在日心惯性系下的坐标；[XO YO ZO]T为小行星质心在日心惯性系下的坐标；λ为小行星固连坐标系与日心惯性系下的一个旋转角。
则太阳光压力的表达式为
| $ {F} = \eta {P_0}A\left( {\frac{{R_{se}^2}}{{R_A^2}}} \right){\cos ^2}{α} $ | (4) |
其中：η为太阳帆总推进系数；P0为在地球附近绝对黑体受到的太阳光压；Rse为日地平均距离，且
| $ {R_A} = \sqrt {X_A^2 + Y_A^2 + Z_A^2} $ | (5) |
| $ \cos {α} = \left( {\frac{{{{[{X_A}{\text{ }}{Y_A}{\text{ }}{Z_A}]}^{\rm T}}[{X_n}{\text{ }}{Y_n}{\text{ }}{Z_n}]}}{{\sqrt {X_A^2 + Y_A^2 + Z_A^2} \sqrt {X_n^2 + Y_n^2 + Z_n^2} }}} \right) $ | (6) |
| $ \left[{\begin{array}{*{20}{c}} {{X_n}} \\ {{Y_n}} \\ {{Z_n}} \end{array}} \right] = \left[{\begin{array}{*{20}{c}} 1& 0& 0 \\ 0& {\cos {λ} }& {\sin {λ} } \\ 0& { - \sin {λ} }& {\cos {λ} } \end{array}} \right]\\ \left[{\begin{array}{*{20}{c}} {\cos {ω} t}& { - \sin {ω} t}& 0 \\ {\sin {ω} t}& {\cos {ω} t}& 0 \\ 0& 0& 1 \end{array}} \right]\\ \left( {\left[{\begin{array}{*{20}{c}} {\cos {φ} } \\ {\sin {φ} } \\ 0 \end{array}} \right] + \left[{\begin{array}{*{20}{c}} x \\ y \\ 0 \end{array}} \right]} \right) + \left[{\begin{array}{*{20}{c}} {{X_O}} \\ {{Y_O}} \\ {{Z_O}} \end{array}} \right] $ | (7) |
角度φ的变化量Δφ为太阳帆绳系系统中的控制量。初始时刻，太阳帆位于小行星的稳定平衡点上，太阳光压力方向与系绳方向相同，即 ${{φ} _0} = \frac{{π} }{2} - {\theta _0}$ 。利用PID控制方法对系统进行控制，则控制律可表示为
| $ \begin{gathered} \Delta {φ} = {K_{\rm P}}e + {K_{\rm D}}\dot e \qquad \qquad \qquad \qquad \qquad \\ \ \qquad {\text{ = }}[{k_{\rm P1}}{\text{ }}{k_{\rm P2}}{\text{ }}{k_{\rm P3}}]{[x - {x^ * }{\text{ }}y - {y^ * }{\text{ }}{ω} - {{ω} ^ * }]^{\rm T}} \\ {\text{ }} + [{k_{\rm D1}}{\text{ }}{k_{\rm D2}}]{[\dot x - {{\dot x}^ * }{\text{ }}\dot y - {{\dot y}^ * }]^{\rm T}} \quad \,\\ \end{gathered} $ | (8) |
其中：x*、y*为小行星固连坐标系中平衡点的位置坐标；w*为经过一段时间后小行星的预期自转角速度； ${\dot x^ * }$ 、 ${\dot y^ * }$ 皆为零。
3 仿真与分析以阿波菲斯小行星(小行星99942)为例进行仿真验证，相关参数如表 1所示。太阳帆的相关参数如表 2所示。则太阳帆在小行星固连坐标系中的初始位置坐标为(0,754.528 6,0)。为便于计算，令系绳与小行星的边缘刚好相切，即θ0 = 18.363 2°，L = 754.824 0 m。当PID参数调整为kP1 = 0.003，kP2 = -0.09，kP3 = 0，kD1 = 0.24，kD2 = 0，仿真时长为85.9 d时，各变量仿真结果如图 2～图 5所示。
|
表 1 小行星99942相关参数 Table 1 Parameters of the asteroid 99942 |
|
图 2 x轴和y轴方向上的位置变化曲线 Fig. 2 Graph of the position change in the direction of x axis and y axis |
|
|
表 2 阳帆相关参数 Table 2 Parameters of the solar sail |
|
图 3 x轴和y轴方向上的速度变化曲线 Fig. 3 Graph of the velocity change in the direction of x axis and y axis |
|
图 4 小行星自转角速度变化曲线 Fig. 4 Graph of the change of the asteroid’s rotation rate |
|
图 5 控制量变化曲线 Fig. 5 Graph of the change of the control variable |
图 2中所示两条曲线分别为太阳帆相对初始位置的偏移量在x轴和y轴上的分量，可看出这两个分量的量级分别为10-3 m和10-4 m。
图 3中所示两条曲线分别为太阳帆相对初始时刻的速度变化量在x轴和y轴上的分量，可看出这两个分量的量级都为10-6 m/s。
图 4所示为小行星自转角速度的变化曲线，当仿真时长为85.9 d时，小行星的自转角速率能够减小至3.435 4×10-8 rad/s。
图 5所示为控制量Δφ的变化曲线。
由图 2到图 5可以看出，在适量调整太阳帆姿态的情况下，可以将太阳帆稳定在小行星的初始平衡点处，通过系绳传递的拉力有效减缓甚至消除小行星的自旋。
4 结束语本文通过建立绳系太阳帆在小行星固连坐标系中的动力学模型，控制太阳帆姿态将其稳定在小行星初始平衡点处，借由系绳传递的拉力给小行星持续提供一个与其自转方向相反的力矩，使小行星自转角速率逐渐减小。仿真结果表明，面积为106 m2的太阳帆，经过约86 d可将小行星的自转消除，验证了该方法的有效性。
当然，系绳长度、不同小行星附着点等对太阳帆控制及对小行星的消旋效果也会产生影响，我们将在下一步的工作中对这些内容进行详细研究。
| [1] | 李爽, 崔平远. 着陆小行星的滑模变结构控制[J]. 宇航学报, 2005, 26(6): 808-812.Li S, Cui P Y. Variable structure with sliding-mode control for landing on asteroids[J]. Journal of Astronautics, 2005, 26(6): 808-812.( 1) |
| [2] | 王开强, 李志海, 张柏楠. 载人小行星探测的飞行模式[J]. 载人航天, 2014, 20(1): 89-94.Wang K Q, Li Z H, Zhang B N. Mission mode of the human asteroid exploration[J]. Manned Spaceflight, 2014, 20(1): 89-94.(1) |
| [3] | 赵建文, 李人杰, 王晓慧, 等. 一种载人小行星探测目标星初选方法[J]. 载人航天, 2014, 20(6): 574-579.Zhao J R, Li R J, Wang X H, et al. A primary target selection method for human asteroid exploration[J]. Manned Spaceflight, 2014, 20(6): 574-579.(1) |
| [4] | 徐伟彪, 赵海斌. 小行星深空探测的科学意义和展望[J]. 地球科学进展, 2005, 20(11): 1183-1190.Xu W B, Zhao H B. Deep space exploration of asteroids: the science perspectives[J]. Advances in Earth Science, 2005, 20(11): 1183-1190.(1) |
| [5] | 于洋, 宝音贺西. 小天体附近的轨道动力学研究综述[J]. 深空探测学报, 2014, 1(2): 93-104.Yu Y, Baoyin H X. Review of orbital dynamics in the vicinity of solar system small celestial bodies[J]. Journal of Deep Space Exploration, 2014, 1(2): 93-104.(1) |
| [6] | Hu W D. Orbital motion in uniformly rotating second degree and order gravity fields[D]. Michigan: University of Michigan Ph.D. thesis, 2002: 13-17.(2) |
| [7] | 张振江, 崔祜涛, 任高峰. 不规则形状小行星引力环境建模及球谐系数求取方法[J]. 航天器环境工程, 2010, 27(3): 383-388.Zhang Z J, Cui H T, Ren G F. Modeling for the gravitation potential environment of an irregular-shaped asteroid and the spherical harmonic coefficient estimation[J]. Spacecraft Environment Engineering, 2010, 27(3):383-388.(1) |
| [8] | Hu W D, Scheeres D J. Numerical determination of stability regions for orbital motion in uniformily rotating second degree and order gravity fields[J]. Planetary and Space Science, 2004(52): 685-692.(1) |
| [9] | Hu W D, Scheeres D J. Spacecraft orbital motion around asteroids[J]. Journal of Astronautics. 2006, 27(5): 961-964.(1) |
| [10] | Zhang Z J, Cui H T, Cui P Y. Determining of a novel periodic orbits around equilibrium points of asteroids[J]. Transactions of the Japan Society for Aeronautical and Space Sciences, 2012, 55(1):44-50.(1) |
| [11] | 倪彦硕, 宝音贺西, 李俊峰. 考虑太阳摄动的小行星附近轨道动力学[J]. 深空探测学报, 2014, 1(1): 67-74.Ni Y S, Baoyin H X, Li J F. Orbit dynamics in the vicinity of asteroids with solar perturbation[J]. Journal of Deep Space Exploration, 2014, 1(1): 67-74.(1) |