长征运载火箭发射地火直接转移轨道研究

耿光有，王珏，侯锡云，余梦伦，王建明，张志国

（1.北京航空航天大学 宇航学院，北京100083； 2.北京宇航系统工程研究所，北京100076；3.南京大学 天文与空间科学学院，南京210023）

为充分发挥新一代长征运载火箭［1］氢氧入轨级动力系统的高比冲优势，中国首次探测火星工程（简称探火工程）［2］中，采用了长征运载火箭在近地停泊轨道滑行至预定位置后，二次点火工作将探测器直接送入地火转移轨道的发射方式，显著提高了探测器有效仪器设备的装载量，但同时也明显增加了运载火箭的发射难度，如需要大型低温长征运载火箭在发射日实现相应发射轨道的零窗口准时发射等。而由于液氢约－253℃的严酷低温工作条件限制，使得其点火工作条件保障困难，主要包括地面起飞时预冷好的点火条件，以及在停泊轨道上较长时间滑行后满足二次点火工作的条件等［3-4］。国际上，虽然美国已经实现了低温入轨级具备近地滑行一圈后再点火的能力，但即便其成熟火箭Delta IV（采用氢氧推进剂），在近年来也多有媒体报道其地面点火发射一再推迟。而为了满足低温氢氧入轨级的严酷二次点火工作要求，长征运载火箭在当前采用了较为保守的设计处理措施［3-4］，因此需要严格限制最长滑行时间不能超过设计门限。为了能够在每26个月才能有一组较合适奔火任务的发射机会中实现及时发射，国际上普遍采用了持续2～3周的发射日期窗口与每天至少30min的发射窗口，以及相应的发射轨道设计以保证深空探测任务的顺利实施［5-6］。新研制的大型长征运载火箭在滑行时间严格受限下，再兼顾测控布站限制等，完成持续2～3周发射日期窗口的探火发射轨道设计尤为艰难，本文即研究如何解决这一问题。

不考虑滑行时间严格受限和发射场航落区严格约束下，地火转移轨道的设计已有较成熟的理论［7］，通常根据Lambert求解逃逸速度并采用圆锥曲线拼接法完成初值设计，再借助B平面矢量法［8］与微分修正完成行星际转移轨道设计。B平面矢量法由于能够使计算初值快速收敛而得以普遍应用。国内外此方面的文献较多，如文献［9］结合B平面矢量法给出了地球到火星转移轨道的近地发射段双曲线轨道的几何约束分析，但没有长征运载火箭需要关注的停泊轨道最长滑行时间等具体约束限制下的深入分析研究等。文献［10］分析了2009—2024年地球到火星间往返发射窗口机会及发射能量C3与逃逸赤纬等参数的数值列表及扼要分析等，但同样不考虑运载火箭停泊轨道的滑行时间等约束的限制。文献［11］主要以运载火箭发射重量及C3为约束，分析了探测器轨道设计，指出火星探测器全部的双曲线出发轨道近地点包络为一沿地球逃逸速度矢量VE∞为轴、以r0sinβ为半径的圆轨迹（r0为近地点高度，β为近地点拱线与双曲渐近线夹角）。文献［12］根据探测器火星目标轨道高度、倾角及天顶距约束分析了地球出发轨道速度、近地点幅角与赤经等参数设计需求，然而发射场位置、射向及运载火箭滑行时间等约束对发射轨道的限制关系并未涉及。文献［13］将B平面矢量法推广至制导策略设计中。文献［14］分析了火星探测转移轨道的偏差传播，指出从近地点到近火点的位置误差放大约100万倍、速度误差放大约10万倍，实际飞行中途修正不可缺省。此外，一般的近地轨道卫星发射任务，只要求运载火箭的起飞时刻、星箭分离点的半长轴、偏心率、轨道倾角和近地点幅角等5个参数满足要求即可以完成发射轨道设计，而火星等深空探测发射任务的最终目标点是近火点目标轨道，因此除要求运载火箭以上参数严格满足要求外，还需要星箭分离点的升交点经度与真近点角等2个参数也必须满足探测器要求，而这2个参数与实际飞行路径直接关联，因此如果按照传统的轨道界面分工只能采用双方反复拼接的设计方法，计算效率与精度等显然难以满足当前设计需求。实际上，复杂、多约束的运载火箭发射轨道设计对探测任务能否顺利实施影响显著。例如，2003年欧洲航天局采用联盟号及Fregat上面级火箭发射火星快车探测器时，即遭遇到上述问题［15］：起初仅有10天的预期发射窗口及1 120 kg的运载能力，发射风险较高，后来经过联盟号火箭优化停泊轨道倾角及二次点火条件，以及优化调整了整流罩与火箭三子级落区，最终将运载能力增加到1 200 kg，发射日期窗口拓展到了30天，并且优化后简化了控制系统装订参数——采用9套飞行程序满足了持续30天的发射日期窗口，仅需要探测器在飞行第10天时深空机动10m／s，这些措施有力保障了探测工程的顺利实施。美国DeltaⅡ7925型火箭发射探测火星等深空任务时，通过固定93°和99°两个射向保障了航落区安全，简化了控制系统装订参数，优化确定了连续多日的发射日期窗口［5-6］。新一代长征运载火箭从海南发射火星探测器时，由于低温入轨级2次起动间最长滑行时间严格受限［3-4］，使入轨点位置（即轨道近地点幅角）的调节范围明显减小，为满足相应发射年份下2～3周的有效发射日期窗口，分析发现不排除结合探测器深空速度脉冲机动以完成发射轨道总体优化［15-16］，故当前文献资料中的设计算法难以直接满足长征运载火箭发射轨道深入优化的需要。依此背景，本文重点研究了长征运载火箭在当前滑行时间、发射场射向及落区等多约束限制下，发射探测器行星际转移轨道的优化设计算法。该算法直接采用探测器近火点目标轨道参数，并结合近地出发轨道诸多设计约束开展设计，在高精度力学计算模型基础上，提出双向微分修正算法，从而也一并解决了在转移轨道的深空对接点处增加速度脉冲机动，以进一步完成发射轨道优化的问题。

1 长征运载火箭海南发射场约束

工程中从飞行安全角度考虑，一般允许的射向A0范围约为90°～115°。

注意到，地球自转（赤道上转速为465.1m／s）对实现长征运载火箭目标入轨速度的贡献Vlaunch为

式中：фL为发射场纬度。

依据齐奥尔科夫斯基公式［17］，动力比冲Isp，火箭起始质量m0，燃烧结束质量mk产生的速度增量ΔV为

对纬度约20°的发射场，如果射向从90°增加到115°，按照氢氧发动机比冲约4 300 m／s，结合探测器与火箭末级入轨质量的大致相对关系［18］，易知发射火星探测器等深空任务的运载能力损失约3%，故可能情况下需要优先选择东射向发射。

2 火星探测任务发射轨道设计

2.1 发射约束下初始轨道搜索

2.1.1 Lambert问题求解

忽略地球和火星的引力及实际尺寸，认为转移轨道直接从地球（在日心坐标系下的位置）出发，到达火星（在日心坐标系下的位置），如图1（a）所 示，该 问 题 为 典 型 的 日 心 Lambert问题：

式中：tM和tE分别为抵达和出发时刻；μ为日心天体引力常数；rM和rE分别为抵达和出发时刻地球与火星距日心距离；c为地球与火星的直线距离；a为转移轨道半长轴。

关于Lambert问题的求解已有非常成熟的算法，本文不再赘述，推荐采用Izzo的计算方法［19］，如果轨道近似为抛物线，推荐采用Gooding的算法以避免Izzo算法中的奇点问题［20-21］。图1（a）中，Vd和Va分别为Lambert转移轨道影响球处出发和抵达速度矢量，VE和VM分别为地球和火星的速度矢量，V0和V1分别为Lambert转移轨道相对于地球与火星影响球的速度矢量。通过V0和V1可以进一步设计地球影响球（SOI）内的双曲线飞行轨迹及火星影响球内的双曲线飞行轨迹，从而将整个转移轨道由地球影响球内的双曲线弧段＋日心Lambert弧段＋火星影响球内双曲线弧段构成（见图1（b）），此即经典圆锥曲线拼接模型。

图1 从地球转移至火星的飞行轨道示意图Fig.1 Schematic diagram of Earth-to-Mars flight trajectory

2.1.2 出入影响球的双曲线轨道初值设计

本文着重介绍地球影响球内双曲线弧段设计方法，火星影响球内双曲线弧段设计方法类似。有别于传统的借助B平面矢量法，本文算法直接根据近地端发射场的约束及火星端目标轨道的约束设计双曲线弧段。通过日心Lambert问题求解，得探测器相对于地球影响球边界速度矢量Vo（地球与火星端分别对应图1中V0和V1），地球段采用下述计算单位：

式中：aE为地球赤道半径；mE为地球质量；［］表示括号内变量在不同区域飞行段是不同参量，此即指无量纲化，如地球段、日心飞行段、火星飞行段等。

设定近地点高度hpE，从式（6）得偏心率eo：

双曲线轨道倾角io与发射场射向AL及发射场纬度фL关系为

从式（8）中，得双曲线轨道近地点幅角ωo和升交点赤经Ωo为

式中：p＝cosiosinωosinhEo－ΓcosiocosωocoshEo；q＝ΓsinωocoshEo＋cosωosinhEo；μE为地球引力常数为速度矢量Vo的分量；Eo为地球影响球半径ro处偏近点角，即

结合式（10）得对应的真近点角fo：

式（8b）一般同时求得2个近地点幅角ωo，代入式（8a）分别得对应的升交点赤经Ωo。根据式（11）计算双曲线轨道纬度幅角uo：

为便于区别，工程上通常以uo∈（－π／2，π／2）作为升交出发，以uo∈（π／2，3π／2）作为降交出发。

由式（5）～式（11）可以给出地球至火星转移轨道在脱离地球影响球处的Kepler轨道根数ao、eo、io、Ωo、ωo、fo，而星箭分离点的地心轨道根数与该轨道只是真近点角fI不同，fI由火箭发射段确定。研究表明，不同射向及滑行时间对应的fI变化不大，可近似认为是常值。

综上，在地球影响球内，以地球和探测器组成二体问题动力学模型，根据日心Lambert弧段给出的地球影响球边界处的速度矢量Vo，可得J2000系下双曲线出发轨道星箭分离点轨道根数ao、eo、io、Ωo、ωo、fI，并由之得位置矢量rJ2000和速度矢量再通过式（12）坐标转换［22］得星箭分离时刻与地固CGCS2000系瞬时重合（之后保持惯性空间不动）的地心赤道惯性系（以下标G表示）下位置rG和速度

式中：Q（t）为岁差章动矩阵；R（t）为地球转角矩阵；W（t）为极移矩阵。由再得相应的轨道根数需求：ao、eo、ioE、ΩoE、ωoE、fI。

抵达近火点双曲线轨道计算在火星J2000平赤道坐标系（用下标ME表示）下描述，该坐标系与J2000火心坐标系（用下标MJ表示）转换关系［22］见式（13），IAU年度工作报告对火星北极指向α0、δ0的精确数据每三年互联网发布一次。

式中：矩阵M下标x，y，z分别表示绕相应坐标轴的转换矩阵［17，22］。火星影响球内的双曲线弧段设计方法与地球影响球的双曲弧段设计方法类似但不是本文关注的重点，所以火星端的双曲线轨道设计这里不再赘述。

2.1.3 长征运载火箭滑行时间及发射场设计约束

长征运载火箭发射轨道星下点轨迹几何示意见图2，采用圆停泊轨道，忽略各摄动因素，得

图2 探测火星地球出发轨道的几何图示Fig.2 Earth-to-Mars launch trajectory geometry illustration

设θ1为从长征运载火箭起飞到进入停泊轨道前的工作弧段，θ2为停泊轨道滑行结束至探测器分离的弧段（分析表明θ1、θ2可以近似为常值），其他符号含义同前，据各参量物理含义，易得所需长征运载火箭入轨级2次起动间滑行时间Th为

图2中，Υ表示春分点，αo和δo分别为地球至火星转移轨道在地球影响球处的逃逸赤经和赤纬，其他符号含义均在上下文中有明确说明。根据长征运载火箭发射探测器入轨时刻升交点赤经Ω及起飞时刻发射场恒星时角αL定义，结合图2球面三角关系易得

式中：λL为发射场经度；Ωe为探测器分离点时刻发射轨道计算的升交点经度；ωe为地球自转角速度；tf0为火箭滑行时间Th＝0 s时从起飞到探测器分离的时间。对选定的射向AL，预先选择很小的阈值εi＞0（i＝1，2），调整发射当日出发时刻tE使：

易得满足长征运载火箭最长滑行时间限制的飞行轨道参数及起飞时刻初解。

2.1.4 射向对飞行轨道及航程角的影响

长征运载火箭关注的发射日出发轨道及航程角与射向的初步关系分析如下。结合图2，得飞行航程角θ（θ＝ωoE＋fo－d）与地球影响球处逃逸赤纬δo、射向AL及发射场纬度фL关系：

设фL＝20°，可得图3。最长滑行时间限定了近地点幅角ωoE与飞行航程角θ的变化范围，每一组航程角对应一个起飞时刻。

结合图3，在射向范围0°～180°情况下，结合发射轨道基础理论，分析易知：

2）当δo＝фL，仅有一组射向随航程角θ减小而增大的解；当δo＝－фL，除了一组射向随航程角θ减小而增大的解外，还有航程角恰好为180°、射向任意的另外一组解。

图3 射向A L、飞行航程角θ与δo关系Fig.3 Contourmap ofδo with respect to launch azimuth A L and flight range angleθ

2.1.5 结合Pork-chop图完成初始发射轨道优选

采用上述加入滑行时间限制再精确双向微分修正的设计算法可以较快速完成初步发射轨道搜索，通过对近地点与近火点速度增量辅以Porkchop图［7］分析，即可以初步完成多设计约束下，长征运载火箭发射火星探测轨道的初步运载能力确定、连续发射日期窗口初选和初始发射轨道设计。限于篇幅，此部分内容从略。

2.2 精确设计发射轨道

2.2.1 精确动力学方程

在起飞时刻地心惯性系下，长征运载火箭飞行的质心动力学矢量方程［17，23］为

式中：m为火箭质量；rG0为惯性系原点到火箭质心矢径；P为火箭推力；R f为作用在火箭上的气动力；g为作用在火箭上的地球引力加速度。通过式（20b）转换得入轨时刻地心惯性系下位置和速度：

进而得相应的发射轨道根数：ao、eo、ioE、Ωe、ωoE、fI。

从探测器分离开始，设定整个地火转移轨道都是在日心J2000坐标系下完成积分，动力学方程为

式中：μS为太阳引力系数；j∈⊙表示太阳系八大行星及月球等；μj为太阳系下扰动天体引力系数，对于地火转移轨道设计，需要考虑八大行星及月球的影响分别为地球和火星的非球形摄动力［24-25］；δE、δM为开关函数，当探测器在地球（火星）影响球界内取1，其余取0；Fothers指地球辐射压、高层空间大气、火星稀薄大气、小行星引力等；FSRP为太阳辐射压，设定探测器面质比为S／m，S为表面积，ms为质量，忽略行星阴影区影响，则

其中：κ为探测器反射系数；c0为光速；ρ⊙为距离r⊙处的太阳辐射流。

采用RKF7（8）阶龙格库塔算法，分析计算表明，探测器在靠近地球或火星飞行段均采用10×10阶非球形引力项，比两者均采用J2项时设计出的转移轨道在近地端的速度差异小于0.1m／s、位置差异小于30 m，在近火端的速度差异约0.3m／s、位置差异约140m；光压对探测器的摄动影响量一般可能超过近地端高阶摄动影响而小于火星端高阶项影响；注意到光压影响与探测器具体形态关系较密切等因素，故综合后下面以考虑大行星及月球等中心天体引力、地球与火星飞行段均采用10×10阶次非球形引力项给出分析算例。

2.2.2 双向微分修正算法实现精确对接

记星箭分离点时刻（记为t0）实际飞行轨道与参考轨道的状态量（即位置和速度）偏差为ΔX（t0），小偏差条件下，t时刻状态量偏差ΔX（t）近似满足：

式中：Φ（t，t0）为参考轨道t0时刻至t时刻的状态转移矩阵［7］。由于从近地点到近火点转移轨道的误差传递特性急剧放大，因此直接采用近火星目标轨道根数求解近地点出发轨道极其困难，一般需要转化处理。

若分别从近地点、近火点向地球至火星转移轨道中间的某个时间点tm积分，记地球出发端的状态量偏差为ΔX（t0），火星到达端的状态量偏差为ΔX（tf），类似式（21），有

式中：ΔXE（tm）为从地球端正向积分至tm时刻轨道的状态量偏差；ΔXM（tm）为从火星端逆向积分至tm时刻轨道的状态量偏差。对同一条转移轨道，要求这2个偏差相同，即

再将式（21）中的t取为tf，并结合式（22）易得）

即时刻t0～tf的状态转移矩阵Φ（tf，t0）可分解为2个分弧段的状态转移矩阵ΦE（tm，t0）和ΦM（tm，tf）；由于将状态转移矩阵Φ（t，t0）一分为二后，两积分段的误差传递特性急剧降低，即可以较顺利地实现问题求解。实际上B平面矢量法的本质原理亦是如此，只是固定为目标天体影响球边界处的B平面参数。分析表明，由于日心转移弧段的良好解析性，使中间积分结合点可以自由移动而不影响计算结果，所以就解决了后续在结合点处增加深空机动以进一步优化长征火箭发射轨道的问题。下面对算法进一步介绍。

以rJ2000和.rJ2000代表在地心J2000坐标系下地球出发时探测器与长征运载火箭分离点的位置和速度，以rMJ和.rMJ代表探测器火心J2000坐标系下抵达近火点的位置和速度，组成一个12维的状态矢量Y，即

式中：符号R表示相对日心系的位置，下标d表示地球影响球界出发点，a表示火星影响球界抵达点。记下标E表示探测器分离点时刻地球的日心系位置和速度，下标M表示抵达近火点时刻火星的日心位置和速度，则

在地球出发点，存在3个约束：

在近火双曲线轨道段，也存在3个类似约束：

式中：上标“＋”表示正向积分轨道，上标“－”表示反向积分轨道。基于长征运载火箭发射轨道优化中，后续关注的深空速度脉冲机动优化问题，分析表明，Tm可以从地球影响球边界至火星影响球边界之间自由选择，相应根据深空对接点处速度脉冲增量为g10～12赋值［16］；对于无需深空脉冲速度机动的情况，一般可取Tm＝（T1＋T2）／2。

将式（27）～式（29）统一记为G（Y）＝（g1，g2，…，g12）＝0此时迭代方程为

式中：

可将雅可比矩阵分解为如下形式：

式中：

12项约束对应12个变量，通过迭代求解即完成火星目标轨道的地球出发轨道设计。

2.3 计算流程

考虑发射场约束，分析与精确设计长征运载火箭发射火星探测器转移轨道的设计流程见图4。

图4 长征运载火箭设计探测火星轨道示意图Fig.4 Schematic diagram of Long March launch vehicle for Earth-to-Mars trajectory design

3 分析及算例

计算采用的火星探测器目标轨道约束为：近火点高度500 km，倾角93°；结合初步选定的射向，近地约束为：射向AL＝107°，出发近地点高度200 km，以UTC时间2020-06-19出发，转移时间198天，2021-01-02T23：59升交抵达火星为算例，计算得地球影响球边界处逃逸赤纬约24°，由图3知有滑行时间Th相近的2组初解，采用双向微分修正得表1、表2精确解，表中EJ2000表示地心J2000赤道系，ME表示J2000火心IAU平赤道系［22］，SJ2000代表距近火点时刻tm＝90天时日心J2000系参数，ΔV指对应地球或火星端各自圆停泊轨道的速度增量。

分析表明：①tm可以大幅移动而不影响优化设计结果；②在结合点处速度增量可以为0 m／s或增加到超过500m／s，都能找到良好求解。

若限定最长滑行时间不能超过2 000 s［3-4］，则只有表1中解满足要求；而如果限定最长滑行时间不超过1 100 s或更短，则需要采用上述加入滑行时间限制再精确双向微分修正的设计算法，通过加大射向或结合探测器深空机动等完成发射轨道优化；具体的分析是繁琐复杂的工程问题，已不在本文讨论范围。

将表1、表2中计算结果在同等条件下，与STK仿真分析软件积分对算表明，二者计算精度与所需时长相当，算例中从近地点到深空对接点、从深空对接点到近火点，位置差异约为1 km。

表1 发射轨道T h ＝1 849 s，探测器05∶25释放，地球升交出发Table 1 T rajectory T h ＝1 849 s，p robe released on 05∶25，ascending from the Earth

表2 发射轨道T h ＝2 312 s，探测器03∶15释放，地球降交出发Table 2 T rajectory T h ＝2 312 s，p robe released on 03∶15，descending from the Ear th

4 结 论

本文给出了长征运载火箭允许滑行时间受限和发射场可用射向范围约束下，长征运载火箭发射火星探测任务从地面起飞至近火点轨道设计的完整模型算法。

1）解决了适用于工程运算所需精度的分析研究模型问题。

2）解决了火星探测任务的发射轨道优化设计问题，确定了运载能力、优选了发射日期窗口，完成了航落区优化选择等。

双向微分修正算法具有结合点可以大幅移动，以便于后续分析中探测器采用深空机动后，深入开展长征运载火箭发射轨道再优化分析的优势。工程应用表明此方法稳定、可靠好用。除了可以用于地火转移发射轨道设计外，本文算法还可用于长征运载火箭发射其他深空探测器的发射轨道设计。

致谢衷心感谢国家航天局、国家留学基金管理委员会对本课题的资助与支持，以及中国科学院大学李明涛教授在STK软件验证中给予的帮助。