太阳帆推进任务的快速仿真方法*

钱 航，郑建华，吴 霞，高 东，刘宇飞

(1.中国运载火箭技术研究院总体设计部，北京100076; (2.中国科学院国家空间科学中心，北京100190; (3.中国空间技术研究院钱学森空间技术实验室，北京100094)

太阳帆推进任务的快速仿真方法*

钱 航1，郑建华2，吴 霞2，高 东2，刘宇飞3

(1.中国运载火箭技术研究院总体设计部，北京100076; (2.中国科学院国家空间科学中心，北京100190; (3.中国空间技术研究院钱学森空间技术实验室，北京100094)

研究太阳帆的力学特性和轨道控制设计方法，导出太阳帆的无奇点控制律.提出通过STK中MATLAB语言编写的嵌入式脚本(Plug-in Script)来将由控制律得到的光压力加速度矢量，添加到STK轨道计算力学模型中，从而进行轨道控制的方法.仿真结果表明，对于常规方法难以进行仿真分析的航天器动力学模型(如太阳帆)，所提出的方法能快速灵活地支持其相应的任务，并增强任务场景的可视化，从而实现利用STK丰富的功能特性进行复杂航天任务的设计、分析和验证.

复杂轨道设计;卫星工具包;太阳帆;小推力控制

0 引 言

随着深空探测活动的广泛开展和深空探测技术的进步，采用普通的推进技术已经很难满足日益发展的深空科学任务要求，例如将探测器推进到太阳系深处，或者使探测器脱离黄道面到达太阳极轨进行对日观测，或者实现日心悬浮轨道等非开普勒轨道等.开发更为有效的推进技术逐渐成为人们研究的焦点，其中太阳帆就是一个重点的研究对象.太阳帆航天器是指利用太阳光压力获得推进力进行宇宙航行的新型飞行器，它被认为是未来人类探索深空的最实际而有效的途径之一.

目前航天任务中绝大部分轨道设计和控制的方法都是基于开普勒轨道理论.随着深空探测和国防建设的深入展开，传统开普勒轨道理论已经不能满足人类未来航天活动的需要.利用太阳帆可以实现非开普勒轨道，比如绕日极轨［1-2］、日心或地心悬浮轨道［3-5］，在空间科学探测和空间态势感知上具有广阔的应用前景.

开展太阳帆非开普勒轨道设计与控制方面的研究将可以为我国未来深空探测任务提供创新性的轨道设计方案，具有重要意义.由于太阳帆力学模型，所涉及空间关系，以及姿态和轨道耦合特性等比较复杂，使得开发用于太阳帆探测任务设计、分析与验证软件成本较高.若能利用现有成熟的软件或者对现有成熟软件进行再次开发，以完成上述设计与分析等任务，则将不仅能够缩短研究周期、减少开发成本，更能够拓展任务设计和分析能力.

卫星工具包STK以其在时域和空域的专业分析方面具有极高的准确性和高度集成的三维可视化等优势，成为航天工业领先的商品化仿真分析软件［6］.但STK不能直接应用于太阳帆探测任务设计、分析与仿真验证.原因是多方面的，STK只提供脉冲式和小推力式推进动力学模型，没有提供太阳帆推进动力学模型;即使应用STK自带的太阳辐射压模型，因太阳帆轨道和姿态耦合的特点，则必须编写复杂的姿态数据文件.

本文提出太阳帆轨道设计方法，并以此基于STK提供的基本功能模块，利用STK自身的灵活性与MATLAB协同工作，进行复杂太阳帆任务的仿真分析.本文对太阳帆探测任务做如下合理简化:1)仅考虑二体模型，多体效应可以从STK选项中加入;2)仅考虑理想太阳帆帆面模型，非理想帆面因素可以在MATLAB的M文件中编写加入;3)不考虑太阳帆挠性振动对姿态和轨道的影响，可以通过本文的方法加入计算分析.在此基础上，着重解决应用STK对特殊的太阳帆动力学模型进行任务设计和分析的难题.

1 太阳帆轨道控制设计方法

由于太阳帆所受的光压力很小，其轨道动力学在很多方面与小推力类型航天器相似，即连续的小推力作用在航天器上持续改变轨道，但太阳帆光压力来源于太阳，故推力方向不可以朝向太阳.仅考虑二体模型，则太阳帆只受太阳引力和光压力:

式中:r为从太阳质心指向太阳帆质心的位置矢量; μ为太阳引力常数;FSRP为太阳光压力.

定义轨道坐标系，轴向为从太阳与太阳帆连线指向太阳帆方向，横向在轨道平面内与轴向垂直并指运动方向，法向与轴向和横向构成右手坐标系.光压力FSRP在该坐标系中轴向、横向和法向的加速度分量分别为［7］

式中:β为太阳帆的光压因子;α和δ分别为太阳帆的锥角和转角，其空间位置关系如图1所示［4］;图1中坐标系上的r⌒、q⌒和 p⌒分别为轴向、横向和法向单位矢量;n为帆面法向矢量.

图1 轨道坐标系中太阳帆锥角和转角Fig.1 Solar sail of cone and clock angles

太阳帆轨道控制归根结底是姿态控制，即控制太阳帆帆面朝向来控制推力矢量，进而实现对飞行轨道的控制.联系具体任务时，可视需要改变轨道根数来设计控制律，如太阳帆绕日极轨设计就需要抬升轨道倾角i.由于太阳帆不用考虑燃料消耗问题，仅从时间最优的角度出发，控制律使轨道倾角变化率取到极大值甚至是最大值，而这可以对太阳帆轨道六根数形式的动力学方程求导得出［8］

式中:sign为符号函数;f为真近点角;ω为近心点角距.

同样，增大轨道长半轴a的无奇点控制律为

式中e为轨道偏心率.

增大轨道偏心率e的无奇点控制律为

其锥角和转角形式同公式(5)和(6).

增大升交点赤经Ω和近心点角距ω的控制律与公式(3)相同.如果需要减小轨道根数，则在相应锥角α的控制律反号即可.轨道根数耦合作用的影响在轨道设计时是需要引起注意的.

对给定的复杂太阳帆轨道任务进行分阶段轨道设计，每个轨道段只调节一种轨道根数，调节轨道根数顺序应遵循整体飞行时间最短和调整次数最少的原则进行.太阳帆采用上述相应的轨道控制律使得该轨道根数以最短时间变化到目标值，再转接到下一个轨道段，拼接起来得到完整的任务轨道，最后视情况对轨道根数进行微调.

2 STK嵌入式脚本

太阳帆探测任务轨道设计的另一关键在于如何将随时间变化的光压力也就是控制力作用于太阳帆动力学模型，本文提出通过STK的嵌入式脚本来实现.STK模块中的Astrogator轨道分析模块允许用户通过MALAB、VB和Perl语言编写，它有自己特定的语法结构.以MALAB语言编写的脚本为例，其函数体结构如下［9］:

其中，XXXX为用户自定义的函数名，PART-1为变量定义部分，PART-2为变量计算部分.STK脚本中输入输出变量的类型和定义方式是固定的，可以查阅STK脚本帮助文档.在太阳帆动力学仿真中，输入变量为轨道六根数，输出变量为计算得到的光压力加速度矢量(脚本不能输出控制力矢量).脚本文件与场景文件放在同一个文件夹下.

(1)变量定义

变量定义方式如下所示:

任意一个变量定义都是由若干个关键词(keyword)与数值对(value)构成的.如上名为“accelerationDesc”的变量是由4个关键词与数值对组成的.加速度变量的关键字'RefName'只有4个值可选，对于太阳帆来说值为'CbiLVLH'，表示选择的是相对中心天体惯性系(CBI)的轨道坐标系(LVLH).

(2)变量计算

如果需要引用定义的变量，就采用output.Acceleration的格式，Acceleration即为accelerationDesc的ArgumentName.将第1节导出的太阳帆轨道控制律公式(3)～(5)代入式(2)计算就可以得到光压力加速度矢量，输出给output.Acceleration(3×1的矢量).

3 太阳帆探测任务仿真

太阳帆在实施深空任务中，需要不断改变其轨道根数以达到目标轨道，同时也要考虑工程约束和天文现象影响.下面以太阳极轨任务(solar polar orbiter)为例，说明太阳帆探测任务的轨道任务仿真及轨道特性分析.

太阳极轨任务为ESA提出设想在2014年5月16日发射的太阳帆航天器，其特征加速度ac= 0.42mm/s2，从运载火箭送入的停泊轨道出发后经过内螺旋飞行至长半轴a=0.48 AU(1AU= 1.495 978 70×1011m，表示地球到太阳的平均距离)，目的是为了将太阳帆送到离太阳较近处，以获得更大的光压力，为下一个阶段抬升轨道倾角节约时间.考虑到材料的耐热性，太阳帆不能距离太阳太近，0.5 AU左右最合适.接着抬升轨道倾角至90°达到极轨，“居高临下”展开太阳日冕物质抛射观测［10］.

图2 太阳帆任务仿真流程Fig.2 Flow chart of solar sail mission

3.1 设置太阳帆初始状态

(1)在卫星Sail的任务控制序列MCS窗口下Initial State页面，选择坐标系Coord.System为Sun MeanEclpJ2000.

(2)在Coordinate Type中选择Keplerian.

(3)Orbit Epoch和轨道根数设置如表1.

表1 太阳帆的初始轨道根数Tab.1 Initial orbit elements of solar sail

3.2 轨道内螺旋

1.1.3 仪器:日立7020型全自动生化分析仪;Beckman-Coulter MAXM血球分析仪;电子称(常熟市双杰测试仪器厂产品);万分之一电子天平(梅特勒-托利多仪器(上海)有限公司产品);Thermo BIOFUGE PRIMO R型离心机。

该段设置太阳帆从地球的停泊轨道(即长半轴a=1 AU)出发后经过内螺旋飞行至长半轴a= 0.48 AU的近日飞行轨道，采用减小轨道长半轴控制律.

(1)在View菜单下选择Astrogator Browser，进入文件夹Calculation Objects下的子文件夹Mean Elems界面，选中Mean Semimajor Axis并点击上面的Duplicate，复制一个新的组件命名为Sun Mean Semimajor Axis.双击该组件，选择坐标系CoordSystem为Sun MeanEclpJ2000，根数类型ElementType为Kozai-Izsak Mean.按照同样的方法添加组件Sun Mean Eccentricity、Sun Mean Inclination、Sun Mean RAAN、Sun Mean Argument of Perigee和Sun Mean True Anomaly.这些组件将在MATLAB的M文件中用到.

(2)按照第2节的方法编写M文件，变量定义为上面新建的组件，根据公式(4)的减小轨道长半轴控制律编写变量计算部分，命名函数为SailSemimajorAxisEval.

(3)进入文件夹Propagator Functions下的子文件夹Plugins界面，选中Astrogator Script Driver，并点击上面的Duplicate，复制一个新的组件命名为Sail Semimajor Axis Astrogator Script Driver.双击该组件，点击项PluginConfig.UseEvaluateFile的值由false改为 true，输 入 PluginConfig.EvaluateFile 值 为SailSemimajor Axis-Eval.m.

(4)进入文件夹Propagators中，选中Helio-centric，并点击上面的Duplicate，复制一个新的组件命名为Sail Semimajor Axis Heliocentric.双击该组件，删除TwoBody Force外的其他天体引力，再插入刚才建的Sail Semimajor Axis Astrogator Script Driver，点击OK关闭.

(5)选中MCS窗口中的DecreSemimajor-Axis，设置 Propagator项为 Sail Semimajor Axis Heliocentric.单击Propagator项后面的“Advanced...”按钮，将Maximum Propagation Time项设置为不选中状态.

(6)点击Stopping Conditions旁边“Insert...”按钮，在弹出的页面中选择UserSelect，单击OK按钮插入.选中Duration，点击Remove按钮删去.点击项User Calc Object右边的“...”按钮，在弹出的页面中选择文件夹Mean Elems中的组件Sun Mean Semimajor Axis，点击OK关闭.设置Trip为0.48 AU，设置Criterion为Cross Either.

3.3 抬升轨道倾角

按照3.2节同样的过程，命名M文件为SailInclinationEval.m，终止条件为90°.

3.4 设置自由飞行段

(1)进入文件夹Propagators中，选中Heliocen-tric，并点击上面的Duplicate，复制一个新的组件命名为Heliocentric Point.双击该组件，删除TwoBody Force外的其他天体引力，点击OK关闭.

(2)选中MCS窗口中的Cruise，设置Propagator项为 Heliocentric Point.单击 Propagator项后面的“Advanced...”按钮，将Maximum Propagation Time项设置为不选中状态.

(3)设置Stopping Conditions下方的Trip为10yr.

3.5 运行和分析MCS

单击Run按钮运行整个MCS，观测三维图形窗口中太阳帆从地球附近出发，内螺旋飞往太阳一定距离时，然后每飞一圈就抬升一定轨道高度直至形成极轨的整个动态过程，最终形成图3.

图3 太阳极轨任务轨道的三维图形显示Fig.3 3D orbit plot of solar polar orbit

4 轨道特性分析

通过上述方法不仅可以实现太阳帆探测任务的快速设计，并且可以利用STK众多分析资源和数据、图形输出功能对太阳帆的轨道特性进行分析［11］.下面以跟踪性能和日凌影响分析为例说明.

4.1 跟踪性能分析

添加喀什站Kashi_China后，单击Facility Tools-Access-Graph下面的“AER...”，即得到喀什站对太阳帆的高度角、方位角和距离的变化历程.

可见在2020年6月到7月的这段时间内，喀什站每天的跟踪弧度约为10 h，喀什站跟踪太阳帆的最远距离为1.680×108km，最近距离为1.200×108km.

4.2 日凌影响分析

从地面站看，当探测器和太阳的角度位置相近时称为日凌，探测器和太阳相对于地面站的角度定义为SEP角.当SEP角度很小时，会对测控性能产生较大的影响.当SEP角小于3°时，X频段上行链路和下行链路的性能都会恶化.日凌期间测控策略就是使探测器工作在既不需要遥控也不需要遥测的模式.

(1)在Object Browser中选中喀什站.

(2)在弹出的页面中单击Create New Vector...按钮，在弹出的页面中设置Name为Sail，Type为Displacement.

(3)单击Origin Point下面的Select...按钮，在弹出的页面中选择Kashi_China下面的Center;单击Destination Point下面的Select...按钮，在弹出的页面中选择Sail下面的Center，单击OK按钮.

(4)返回到 Vector Geometry Tool页面，单击Create New Angle...按钮，在弹出的页面中设置Name为SEP，Type为Between Vector.

(5)单击From Vector下面的“Select...”按钮，在弹出的页面中选择新建的由喀什站指向太阳帆的Sail矢量，单击OK按钮;单击To Vector下面的“Select...”按钮，在弹出的页面中选择Sun矢量，单击OK按钮返回到角度创建页面，再单击OK按钮返回到Vector Geometry Tool页面，单击Close按钮关闭该页面.至此，SEP角就建立了.

图4 喀什站对太阳帆的跟踪弧段Fig.4 Track segment of Kashi facility to solar-sail

(6)在Object Browser中选中喀什站，单击Facility Tool-Graph，在出现的页面中点击Style下面的New按钮.在弹出页面中，单击Angle左边的“+”号展开下面的数据，选中SEP Angle，单击Insert Y Axis按钮，此时在Y Axis下面的方框中出现Angle-SEP Angle.

(7)单击Apply按钮后再单击OK按钮回到STK Graph Tool页面，在Change Name按钮左边的框中输入SEP，再单击Change Name按钮.

(8)单击Graph下面的“Create...”按钮，得到在该段时间内SEP角度变化历程.由图5可见，一年中SEP角度小于3°有两次，每次都不到5天，且SEP角最小为1°.

图5 喀什站的SEP角Fig.5 SEP angle of Kashi facility

5 结 论

仿真结果表明，根据本文提出的太阳帆轨道控制设计方法，利用STK/MATLAB协同工作可以将复杂太阳帆探测任务分解为轨道根数快速拼接，并以此进行轨道特性分析，增强了任务场景的可视化，提高了研究效率.该方法可以推广到任何类型复杂航天器如国际空间站挠性振动的任务分析，太阳帆与深空通信网络之间的可见性分析与链路计算，以及多任务探测器中途与行星交会的可能性分析，为复杂航天任务设计和分析提供了一条便捷高效的途径，具有重要的现实意义和较高的工程参考价值.

［1］GARBE G，MONTGOMERY E.An overview of NASA’s solar sail propulsion project［C］//The 39thAIAA/ ASME/SAE//-ASEE Joint Propulsion Conference＆Exhibit.Huntsville，USA:AIAA，2003.

［2］WIE B，THOMAS S，PALUSZEK M，et al.Propellantless AOCS design for a 160 m，450 kg sailcraft of the solar polar imager mission［C］//The 41stAIAA/ASME/ SAE//-ASEE Joint Propulsion Conference＆ Exhibit.Tucson，USA:AIAA，2005.

［3］钱航，郑建华，于锡峥，等.太阳帆航天器悬浮轨道动力学与控制［J］.空间科学学报，2013，33(4): 458-464.QIAN H，ZHENG J H，YU X Z，et al.Dynamics and control of displaced orbits for solar sail spacecraft.Chinese Journal of Space Science，2013，33(4):458-464.

［4］MATTEO C，MCINNES C R.Generation of optimal trajectories for earth hybrid pole sitters［J］.Journal of Guidance Control and Dynamics，2011，34(3):847-859.

［5］BOOKLESS J，MCINNES C R.Dynamics and control of displaced periodic orbits using solar-sail propulsion［J］.Journal of Guidance Control and Dynamics，2006，29 (3):527-537.

［6］黄文德，郗晓宁，王威.基于STK的月球探测任务分析、设计与验证方法［J］.中国空间科学技术，2010，(6):72-80.HUANG W D，XI X N，WANG W.Methods of mission analysis，design and validation for lunar exploration based on STK［J］.Chinese Space Science and Technology，2010，(6):72-80.

［7］MCINNES C R.Solar sailing:technology，dynamics and mission applications［M］.Chichester:Springer Praxis，1999:118-120.

［8］王松霞，徐世杰，陈统.太阳帆飞行器轨道动力学分析［J］.中国空间科学技术，2006，(2):30-37.WANG S X，XU S J，CHEN T.Orbit dynamics analysis of solar spacecraft［J］.Chinese Space Science and Technology，2006，(2):30-37.

［9］吴霞.小卫星编队飞行队形控制与仿真［D］.北京:中国科学院空间科学与应用研究中心，2006.WU X.Formation control and simulation of satellites formation flying［D］.Beijing:Center for Space Science and Applied Research，Chinese Academy of Sciences，2011.

［10］MACDONALD M，MCINNES C R.The Solar polar orbiter:a solar sail technology reference study［J］.Journal of Spacecraft and Rockets，2006，43(5):960-972.

［11］丁溯泉，张波，刘世勇.STK在航天任务仿真分析中的应用［M］.北京:国防工业出版社，2011:206-212.DING S Q，ZHANG B，LIU S Y.Application of simulation analysis for space mission based on STK［M］.Beijing:National Defense Industry Press，2011:206-212.

Method of Mission Analysis for Solar Sail Exploration Based on STK

QIAN Hang1，ZHENG Jianhua2，WU Xia2，GAO Dong2，LIU Yufei3
(1.Systems Engineering Division，China Academy of Launch Vehicle Technology，Beijing 100190，China;2.National Space Science Center，CAS，Beijing 100190，China;3.Qian Xuesen Laboratory of Space Technology，China Academy of Space Technology，Beijing 100094，China)

Dynamic of solar sail and method of trajectory design have been analyzed.Then the control forces calculated by the control strategy are added into the STK orbit propagator，which is realized through the STK plug-in coded by MATLAB.The simulation results show that the methods proposed can realize expeditiously and efficiently that STK/MATLAB coordinated work to neatly support control simulation of special force model such as solar sail dynamic，and the abundant capabilities of STK are suitable for designing and analysising the complex space mission.

complex trajectory design;STK;solar sail;low thrust control

V529.1

:A

:1674-1579(2016)02-0008-06

10.3969/j.issn.1674-1579.2016.02.002

钱 航(1987—)，男，博士，研究方向为太阳帆航天器动力学、控制与仿真;郑建华(1966—)，女，研究员，研究方向为飞行器动力学、深空探测与轨道设计、鲁棒控制理论及应用等;吴 霞(1981—)，女，副研究员，研究方向为飞行器动力学、深空探测轨道设计与仿真;高 东(1978—)，男，副研究员，研究方向为小卫星编队飞行姿态控制与自主导航算法设计;刘宇飞(1979—)，男，高级工程师，研究方向为太空发电站、太阳帆等新型薄膜航天器的动力学与控制.

*中国科学院战略性先导科技专项资金资助项目(XDA04060303)和太阳帆民用航天资助项目.

2015-12-17