CN114684389A - 考虑再入约束的月地转移窗口及精确转移轨道确定方法 - Google Patents

考虑再入约束的月地转移窗口及精确转移轨道确定方法 Download PDF

Info

Publication number
CN114684389A
CN114684389A CN202210335695.3A CN202210335695A CN114684389A CN 114684389 A CN114684389 A CN 114684389A CN 202210335695 A CN202210335695 A CN 202210335695A CN 114684389 A CN114684389 A CN 114684389A
Authority
CN
China
Prior art keywords
point
reentry
moon
transfer
lunar
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202210335695.3A
Other languages
English (en)
Other versions
CN114684389B (zh
Inventor
乔栋
贾飞达
李翔宇
庞博
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Institute of Technology BIT
Original Assignee
Beijing Institute of Technology BIT
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing Institute of Technology BIT filed Critical Beijing Institute of Technology BIT
Priority to CN202210335695.3A priority Critical patent/CN114684389B/zh
Publication of CN114684389A publication Critical patent/CN114684389A/zh
Application granted granted Critical
Publication of CN114684389B publication Critical patent/CN114684389B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B64—AIRCRAFT; AVIATION; COSMONAUTICS
    • B64G—COSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
    • B64G1/00—Cosmonautic vehicles
    • B64G1/22—Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
    • B64G1/24—Guiding or controlling apparatus, e.g. for attitude control
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/04—Constraint-based CAD
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Remote Sensing (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Geometry (AREA)
  • Evolutionary Computation (AREA)
  • Chemical & Material Sciences (AREA)
  • Combustion & Propulsion (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Computer Hardware Design (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Navigation (AREA)

Abstract

本发明公开的一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，属于航空航天技术领域。本发明实现方法为：根据地面着陆场约束、再入约束确定固连系下再入点状态，提高着陆可靠性；预设月地转移时长和任务搜索日期区间，确定区间内地月引力模型下每天的最优月地转移轨道对应的近月点高度，根据预设的近月点高度偏差限确定月地转移窗口；在高精度摄动模型下，修正月球影响球边界的探测器状态以及对应的月球影响球内高精度转移轨道；在再入点前施加速度修正使得探测器的转移轨道在月球影响球边界处位置精确逆向拼接，得到月球影响球边界至再入点的高精度转移轨道。本发明具有窗口计算效率高、转移轨道精度高、再入舱着陆可靠性高的优点。

Description

考虑再入约束的月地转移窗口及精确转移轨道确定方法
技术领域
本发明涉及一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，尤其适用于考虑月球端出发倾角和高度约束、地球端大气再入角和再入航程约束以及着陆场约束的转移窗口及精确转移轨道的高效确定，属于航空航天技术领域。
背景技术
近年来，我国探测月球的步伐逐渐加快，嫦娥五号成功执行无人月球采样返回任务，我国的载人登月计划也在筹备论证中。月地转移窗口的确定是上述任务可行的前提，而精确月地转移轨道的计算能够确保任务的安全可靠执行。然而，当前发展的关于月地转移窗口及精确转移轨道的计算方法存在计算效率低、窗口精度低的问题，实现转移窗口及精确转移轨道的一体化高效计算具有很强的挑战性。
在已发展的关于月地转移窗口计算方法先技术[1](参见[1]郑爱武,周建平.直接再入大气的月地返回窗口搜索策略[J].航空学报,2014,35(8):2243-2250.)提出基于双二体模型的月地返回窗口的搜索策略，根据当天最小再入角初步确定返回窗口，并在受摄模型下求解精确转移轨道。该方法能够在大范围内快速搜索返回窗口，但双二体模型的应用使得窗口精度较低，且动力学切换处的拼接状态较难给出。
因此在上述背景下，本发明提出一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，实现了考虑月球端出发倾角和高度约束、地球端大气再入角和再入航程约束以及着陆场约束的月地转移窗口及精确转移轨道的高效确定。
发明内容
本发明的主要目的是提供一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，根据地面着陆场约束、再入约束确定固连系下再入点状态，提高再入舱的着陆可靠性；预设月地转移时长和任务搜索日期区间，确定区间内地月引力模型下每天的最优月地转移轨道对应的近月点高度，根据预设的近月点高度偏差限确定月地转移窗口；在高精度摄动模型下，对选定日期的最优月地转移轨道，分别在月球影响球边界和再入前施加两次速度修正，逆向拼接得到近月点至地球再入点的高精度转移轨道。本发明具有窗口计算效率高、转移轨道精度高、再入舱着陆可靠性高的优点。
本发明的目的是通过下述技术方案实现的：
本发明公开的一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，根据地面着陆场经纬度约束、再入点高度约束、再入轨迹航程和倾角约束、预设再入速度大小确定地心固连系下再入点处探测器的状态，提高再入舱着陆的可靠性；预设再入日期和月地转移时长，确定地月引力模型下当天近月点高度最小的最优月地转移轨道；预设任务搜索日期区间，遍历得到区间内地月引力模型下每天的最优月地转移轨道对应的近月点高度，根据预设的近月点高度偏差限确定月地转移窗口，保证窗口确定的可靠性和高效性；预设近月端出发倾角和高度约束，基于得到的月地转移窗口，修正得到高精度摄动模型下月球影响球边界的探测器状态以及对应的月球影响球内高精度转移轨道；采用微分修正，在再入点前一天处施加速度修正使得探测器的转移轨道在月球影响球边界处位置精确拼接，得到月球影响球边界至再入点的高精度转移轨道。
本发明公开的一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，包括如下步骤：
步骤一：根据地面着陆场经纬度约束、再入点高度约束、再入轨迹航程和倾角约束、预设再入速度大小确定地心固连系下再入点处探测器的状态，提高再入舱着陆的可靠性。
步骤1.1：根据地面着陆场纬度、再入轨迹倾角和航程确定再入点的纬度幅角。
预设地面着陆场纬度为
Figure BDA0003574284720000021
再入轨迹倾角为ir，则升轨再入时，根据球面三角形正弦定理，着陆场的纬度幅角UF为
Figure BDA0003574284720000022
降轨再入时，着陆场的纬度幅角UF为
Figure BDA0003574284720000023
预设再入轨迹的航程为Sr，则再入轨迹的航程幅角Ur为
Figure BDA0003574284720000024
其中，Re为地球平均半径。
从而得到再入点的纬度幅角UP为
UP＝UF-Ur (4)
步骤1.2：根据着陆点经度和再入点纬度幅角确定再入点的飞行方位角和再入点经纬度。
预设地面着陆场经度为λF，根据再入点的纬度幅角UP和再入轨迹轨道倾角ir，由球面三角形纳皮尔法则确定再入点到升交点的经度差ΔλP为
Figure BDA0003574284720000031
从而再入点的飞行方位角Q为
Figure BDA0003574284720000032
基于着陆场纬度幅角UF计算着陆场到升交点的经度差ΔλF为
Figure BDA0003574284720000033
则确定再入点的经度λP为
λP＝λF-ΔλF+ΔλP (8)
再入点的纬度
Figure BDA0003574284720000034
为
Figure BDA0003574284720000035
步骤1.3：根据再入点高度约束、再入点经纬度和飞行方位角确定再入点处探测器在地心固连系中的状态。
预设再入点高度为H，则再入点的地心距R为
R＝H+Re (10)
从而确定再入点在地心固连系下的位置矢量为
Figure BDA0003574284720000036
定义再入坐标系原点位于探测器上，Y轴指向天顶方向，X轴在再入轨迹面内垂直Y轴指向飞行方向，Z轴与X轴和Y轴组成右手坐标系。记γ为再入角，vre为预设再入速度大小，则再入坐标系中探测器在再入点处的速度矢量vre为
vre＝[vrecosγ,vresinγ,0] (12)
再入坐标系到地心固连系的旋转矩阵MEF为
Figure BDA0003574284720000037
则地心固连系下再入点处探测器的速度矢量为
ve＝(MEF·vre T)T (14)
其中，上标“T”为转置符号。
进而，地心固连系下再入点处探测器的状态为
Xe＝[re,ve] (15)
上述基于地面着陆场的再入点状态确定保证了再入舱着陆的可靠性。
步骤二：预设再入日期和月地转移时长，确定地月引力模型下当天近月点高度最小的最优月地转移轨道。
步骤2.1：选择预设再入日期一天内的选定时刻，根据预设的月地转移时长，在地月引力模型下将探测器状态从再入点逆推至近月点得到选定时刻的可行月地转移轨道。
由于月球公转速度相比于地球自转速度更慢，所以随着地球自转，每天之中存在一条能到达月球附近的最优月地转移轨道。
由于探测器在月地转移过程中受到的引力作用主要来自地球和月球，故在转移窗口确定前采用地月引力模型进行转移轨道计算，在保证可靠性的前提下大大提高计算效率。建立地月引力模型为
Figure BDA0003574284720000041
其中，
Figure BDA0003574284720000042
为探测器在地月引力模型中的加速度，r为探测器在地心惯性系下的位置矢量，rM为月球在地心惯性系下的位置矢量，μE为地球引力常数，μM为月球引力常数。
预设再入日期时段的儒略日范围为D1零时至D1+1零时，选择再入点时刻为T0∈[D1,D1+1]。月地转移时长预设为ΔT，则近月点时刻为Tm＝T0-ΔT，通过月球星历表获得Tm时刻的月球位置为rM＝[xM,yM,zM]。
再入时刻T0对应的格林尼治时间为GST，则地心固连系至地心惯性系旋转矩阵MFI为
Figure BDA0003574284720000043
则地心惯性系下再入点处探测器的状态为
Xi＝[(MFI·re T)T,(MFI·ve T)T] (18)
将探测器由再入点状态Xi逆推至近月点得到实际转移时长为
Figure BDA0003574284720000047
的最优月地转移轨道。
步骤2.2：根据实际转移时长与预设的转移时长之差，调整预设再入速度大小，得到实际转移时长为预设转移时长的可行月地转移轨道。
记再入速度大小为vre时，实际再入点至近月点转移时长为
Figure BDA0003574284720000044
若
Figure BDA0003574284720000045
说明再入速度设置偏大，探测器提前到达近月点，需减小再入速度大小vre；若
Figure BDA0003574284720000046
说明再入速度设置偏小，探测器滞后到达近月点，需增大再入速度大小vre。
通过迭代调整再入速度ve的大小vre使得近月点至再入点转移时长刚好为ΔT。确定每次迭代中再入速度大小vre的调整量(步长)Δvre是关键，步长过大或过小都会导致转移时长存在偏差。所以，提出一种根据实际转移时长与期望转移时长偏差确定调整步长的方法。
迭代中月地实际转移时长与预设转移时长的偏差为
Figure BDA0003574284720000051
则下次迭代时再入速度大小的调整步长
Δvre＝δt/100(δt单位为天) (19)
再入速度调整为
vre＝vre+Δvre (20)
当δt＜0时，Δvre＜0，即再入速度调小，当δt＞0时，Δvre＞0，即再入速度调大。据此进行再入速度大小迭代，直至满足转移时长精度要求
Figure BDA0003574284720000052
其中，εT为转移时长精度误差限，根据实际工程需要进行设置。
至此，得到转移时长固定为ΔT的再入点至近月点的可行月地转移轨道。
步骤2.3：采用黄金分割法，以近月点高度为指标，确定预设日期一天内的最优月地转移轨道。
地球一天自转一周，所以T0在[D1,D1+1]的一天内变化时，惯性系下可行月地转移的轨道方向会旋转一周，所以步骤2.2得到可行月地转移轨道的近月点高度在一天内近似周期变化，最大近月点高度甚至可达数万km。所以要寻找[D1,D1+1]一天内近月点高度最低的月地转移轨道对应的再入点时间T0。一天内转移时长固定为ΔT的月地转移轨道近月点高度变化曲线呈单峰+单谷趋势，所以考虑采用黄金分割法求得一天内最小近月点高度对应的再入点时间T0。
对于最优月地转移轨道求解问题，黄金分割法实现的具体步骤为：在时间区间[a,b](a＝D1,b＝D1+1)内取点a1,a2把区间分为三段，其中
a1＝a+0.382(b-a)，a2＝a+0.618(b-a) (22)
记F(x)为再入点时间为x时可行月地转移轨道对应的近月点高度，若F(a1)＜F(a2)，令b＝a2，a不变；若F(a1)＞F(a2)，令a＝a1，b不变。实施上述迭代过程，直至满足最优再入点时间精度
b-a＜εrp (23)
其中，εrp为最优再入时刻精度误差限。至此，转移时长为ΔT时的使得近月点高度最小的最优再入点时间T0＝a，从而得到[D1,D1+1]一天内近月点高度最小的最优月地转移轨道。
步骤三：预设任务搜索日期区间，遍历得到区间内地月引力模型下每天的最优月地转移轨道对应的近月点高度，根据预设的近月点高度偏差限确定月地转移窗口，保证窗口确定的可靠性和高效性。
预设任务搜索日期区间为[D2,D3](D3-D2＞1)，采用步骤一和步骤二遍历得到每一天τ∈[D2,D3]中的最优月地转移轨道及其对应的近月点高度H。
经过轨道修正实现月球低轨捕获所需的燃耗与修正前最优月地转移轨道的近月点高度相关。若最优月地转移轨道近月点高度过大，会导致修正燃耗过大，超出探测器能力范围或任务约束。所以，根据每天最小的近月点半径来判断窗口的存在性。
预设近月点高度偏差限为χ，若日期τ的最优月地转移轨道近月点高度H＜χ，则判定日期τ在月地转移窗口中，若H＞χ，则判定日期τ不在月地转移窗口中。从而，可靠高效地得到日期区间[D2,D3]内的月地转移窗口Γ为每一个在转移窗口中的日期τ的集合。
步骤四：预设近月端出发倾角和高度约束，基于得到的月地转移窗口，修正得到高精度摄动模型下月球影响球边界的探测器状态以及对应的月球影响球内高精度转移轨道。
步骤4.1：以步骤三得到的最优月地转移轨道的再入状态为初值，在高精度动力学模型下逆推确定近月点处参数。
在月地转移窗口中任选T0∈Γ，基于步骤一、步骤二、步骤三得到最优月地转移轨道，采用微分修正设计符合月球端出发倾角和高度约束的月球影响球内高精度转移轨道。
以得到的最优月地转移轨道的再入点处状态为初值，在高精度动力学下进行精确转移轨道设计。建立考虑日地月中心引力、日月非球形摄动、太阳光压的高精度动力学模型为
Figure BDA0003574284720000061
其中，μS为太阳引力常数，rS为根据星历文件得到的地心惯性系下太阳的位置矢量，aNE为地球非球形摄动加速度，aNM为月球非球形摄动加速度，asrp为太阳光压摄动加速度。
将最优转移轨道作为设计初值，将初值Xi在高精度动力学式(24)中逆向积分后，得到地心惯性系下的近月点状态并将其转化到月心赤道惯性坐标系中，表示为Xf＝[rf,vf]，以轨道角动量h、半长轴a、偏心率e、半通径p、真近点角f、偏近点角H为中间变量，计算距近月点时间τ，轨道高度rpm，轨道倾角im。
记标量与矢量对应关系为
a＝|a|,e＝|e|,h＝|h|,rf＝|rf|,vf＝|vf| (25)
则通过近月点状态rf和vf推导得到近月点处以下式子
Figure BDA0003574284720000071
记单位向量
n＝[0,0,1] (27)
则距近月点时间τ，轨道高度rpm，轨道倾角im可以通过式(28)获得：
Figure BDA0003574284720000072
步骤4.2：预设月地转移轨道的月球端出发近月点高度约束和倾角约束，采用微分修正确定月球影响边界处探测器的状态和月球影响球内高精度转移轨道。
预设月地转移轨道的月球端出发近月点高度为
Figure BDA0003574284720000073
倾角为
Figure BDA0003574284720000074
将距近月点时间τ、轨道高度rpm、轨道倾角im三个量与预设值之差作为目标约束变量(标称距近月点时间为0)，故月地转移轨道求解转化为如下非线性方程组求解：
[rpm-r′pm,im-i′m,τ]＝0 (29)
将初值Xi在高精度动力学式(24)中逆向积分后，得到地心惯性系下的月球影响球边缘处探测器的状态并将其转换为轨道根数形式，选取修正量为月球影响球边缘处初始状态的半长轴a、升交点赤经Ω、近心点幅角w，通过数值方法分析修正量对目标约束变量的影响。分别给每个修正量初始状态一个小的扰动δ，经过逆向积分到末端近月点后得到目标约束变量对于扰动δ的响应，即得到目标约束变量相对修正量初始状态的敏度矩阵K1为
Figure BDA0003574284720000075
则修正量的微分与目标约束变量误差之间的关系为
Figure BDA0003574284720000081
迭代数次即可得到满足非线性方程组(29)的月球影响球边缘处探测器的初始状态a、Ω、w的值。从而得到近月点倾角和高度满足预设约束的月球影响球边缘处的探测器状态为XSOI＝[rSOI,vSOI]。以此状态为初值，对月球影响球内轨道进行高精度递推，得到近月点状态XM＝[rM,vM]。
探测器从高度为rpm的月球圆轨道施加离轨脉冲进入高精度月地转移轨道时，所需施加的速度增量为
Figure BDA0003574284720000082
其中，RM为月球平均半径。
步骤五：采用微分修正，在再入点前一天处施加速度修正使得探测器的转移轨道在月球影响球边界处位置精确拼接，得到月球影响球边界至再入点的高精度转移轨道。
修正后的月球影响球边缘处探测器状态与修正前状态相比，位置和速度均有小幅度变化，所以需对再入点-月球影响球边缘段轨道进行修正以对准月球影响球边缘处探测器状态。
采用再入点前一天施加修正脉冲的方式使得月球影响球边缘处探测器位置对准，即整条月地转移轨道分为(1)近月点-月球影响球边缘段、(2)月球影响球边缘-再入前一天段、(3)再入前一天-再入点段这三段轨道，施加的脉冲共有三次，即除了每两段之间的连接脉冲之外，还包括式(32)的近月点处离轨脉冲。
第三段轨道递推以再入点初值Xi在高精度动力学模型(24)下逆推得到再入前一天处探测器的状态为Xrr＝[rrr,vrr]，其中rrr＝[xrr,yrr,zrr]，vrr＝[vxrr,vyrr,vzrr]。选择终端约束为步骤四得到的第一段轨道月球影响球边缘处探测器的位置rSOI＝[xSOI,ySOI,zSOI]。分别给再入点前一天处探测器速度vrr的每个分量一个小的扰动δ，经过逆向积分后得到月球影响球边缘处探测器位置受扰动的响应，即可得到月球影响球边缘处探测器位置相对再入点前一天处速度扰动的敏度矩阵K2为
Figure BDA0003574284720000091
从而得到速度修正量的微分和月球影响球边缘处探测器位置误差的关系为
Figure BDA0003574284720000092
迭代数次得到使得月球影响球处探测器位置精确对准的再入点前一天处探测器速度增量Δvrr＝[Δvxrr,Δvyrr,Δvzrr]的值，从而，第二段轨道中再入点前一天处初始速度为v′rr＝vrr+Δvrr。相应的第二段转移轨道中月球影响球边缘处探测器的末端速度v′SOI，则为了与第三段转移轨道中探测器的初始状态对应，月球影响球边缘处探测器应施加的速度增量为ΔvSOI＝vSOI-v′SOI。
至此，得到月地转移窗口Γ内再入点时刻为T0时的高精度月地转移轨道，探测器所需施加的总速度增量为
|Δvall|＝|ΔvM|+|ΔvSOI|+|Δvrr| (35)
还包括步骤六：结合步骤四得到的月球影响球内高精度转移轨道和步骤五得到的月球影响球边界至再入点的高精度转移轨道，得到考虑月球端出发倾角和高度约束、地球端大气再入角和再入航程约束以及着陆场约束的高精度月地转移轨道，根据得到的高精度月地转移轨道对探测器施加相应机动，实现高效率、高精度、高可靠性的月地转移。
有益效果
(1)本发明的考虑再入约束的月地转移窗口及精确转移轨道确定方法，由于是从地面着陆场考虑再入轨迹航程和倾角约束确定的再入点状态，所以该计算方法确保再入舱精确着陆的同时使得再入轨迹精确满足再入角和再入航程约束，再入点状态精度高，对再入约束的适应性强。
(2)本发明的考虑再入约束的月地转移窗口及精确转移轨道确定方法，在简化动力学中进行大范围转移轨道搜索得到月地返回窗口，在高精度动力学模型下确定单条精确月地转移轨道，实现月地转移窗口和高精度月地转移轨道的高效一体化计算，窗口计算效率高，精确轨道精度高。
(3)本发明的考虑再入约束的月地转移窗口及精确转移轨道确定方法，在月球影响球处施加机动能够满足任意近月点高度和倾角出发的高精度月地转移，出发状态变化时对再入段轨迹影响较小，所以工程适应性强。
(4)本发明的考虑再入约束的月地转移窗口及精确转移轨道确定方法，相较于传统的无机动自然转移窗口，得到的月地转移窗口范围更大，能够给月球探测任务提供更多的可行时间选择。
附图说明
图1为本发明公开的考虑再入约束的月地转移窗口及精确转移轨道确定方法流程图。
图2为本发明公开的考虑再入约束的月地转移窗口及精确转移轨道确定方法求解得到的不考虑近月端约束的最优月地转移轨道；
图3为本发明公开的考虑再入约束的月地转移窗口及精确转移轨道确定方法求解得到的2030年每天最优转移轨道近月点高度与转移窗口确定图；
图4为本发明公开的考虑再入约束的月地转移窗口及精确转移轨道确定方法求解得到的满足所有设定约束的高精度月地转移轨道。
具体实施方式
为了更好地说明本发明的目的和优点，下面结合具体实施示例对本发明做出详细解释。
实施例1：
如图1所示，本实施例提供一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，具体实现步骤如下：
步骤一：根据地面着陆场经纬度约束、再入点高度约束、再入轨迹航程和倾角约束、预设再入速度大小确定地心固连系下再入点处探测器的状态，保证再入舱着陆的可靠性。
步骤1.1：根据地面着陆场纬度、再入轨迹倾角和航程确定再入点的纬度幅角。
预设地面着陆场纬度为
Figure BDA0003574284720000101
再入轨迹倾角为ir＝45°，则升轨再入时，根据球面三角形正弦定理，着陆场的纬度幅角UF为
Figure BDA0003574284720000102
降轨再入时，着陆场的纬度幅角UF为
Figure BDA0003574284720000103
预设再入轨迹的航程为Sr，则再入轨迹的航程幅角Ur为
Figure BDA0003574284720000111
其中，Re为地球平均半径。
从而得到再入点的纬度幅角UP为
UP＝UF-Ur (39)
步骤1.2：根据着陆点经度和再入点纬度幅角确定再入点的飞行方位角和再入点经纬度。
预设地面着陆场经度为λF＝101.5°E，根据再入点的纬度幅角UP和再入轨迹轨道倾角ir，由球面三角形纳皮尔法则确定再入点到升交点的经度差ΔλP为
Figure BDA0003574284720000112
从而再入点的飞行方位角Q为
Figure BDA0003574284720000113
基于着陆场纬度幅角UF计算着陆场到升交点的经度差ΔλF为
Figure BDA0003574284720000114
则确定再入点的经度λP为
λP＝λF-ΔλF+ΔλP (43)
再入点的纬度
Figure BDA0003574284720000115
为
Figure BDA0003574284720000116
步骤1.3：根据再入点高度约束、再入点经纬度和飞行方位角确定再入点处探测器在地心固连系中的状态。
预设再入点高度为H＝120km，则再入点的地心距R为
R＝H+Re (45)
从而确定再入点在地心固连系下的位置矢量为
Figure BDA0003574284720000117
定义再入坐标系原点位于探测器上，Y轴指向天顶方向，X轴在再入轨迹面内垂直Y轴指向飞行方向，Z轴与X轴和Y轴组成右手坐标系。记γ为再入角，vre为预设再入速度大小，则再入坐标系中探测器在再入点处的速度矢量vre为
vre＝[vrecosγ,vresinγ,0] (47)
再入坐标系到地心固连系的旋转矩阵MEF为
Figure BDA0003574284720000121
则地心固连系下再入点处探测器的速度矢量为
ve＝(MEF·vre T)T (49)
其中，上标“T”为转置符号。
进而，地心固连系下再入点处探测器的状态为
Xe＝[re,ve]＝[[4314.9,4783.6,851.5]km,[-7.033,3.535,7.248]km/s] (50)
上述基于地面着陆场的再入点状态确定保证了再入舱着陆的可靠性。
步骤二：预设再入日期和月地转移时长，确定地月引力模型下当天近月点高度最小的最优月地转移轨道。
步骤2.1：选择预设再入日期一天内的选定时刻，根据预设的月地转移时长，在地月引力模型下将探测器状态从再入点逆推至近月点得到选定时刻的可行月地转移轨道。
由于月球公转速度相比于地球自转速度更慢，所以随着地球自转，每天之中存在一条能到达月球附近的最优月地转移轨道。所以对每天的月地转移轨道分别进行计算。
由于探测器在月地转移过程中受到的引力作用主要来自地球和月球，故在转移窗口确定前采用地月引力模型进行可行转移轨道计算，在保证可靠性的前提下大大提高计算效率。建立地月引力模型为
Figure BDA0003574284720000122
其中，
Figure BDA0003574284720000123
为探测器在地月引力模型中的加速度，r为探测器在地心惯性系下的位置矢量，rM为月球在地心惯性系下的位置矢量，μE为地球引力常数，μM为月球引力常数。
预设再入日期时段的儒略日范围为[D1,D1+1]＝[2462777.5,2462778.5](2030/10/0300:00:00–2030/10/04 00:00:00)，选择再入点时刻为T0∈[D1,D1+1]。月地转移时长预设为ΔT＝3.0天，则近月点时刻为Tm＝T0-ΔT，通过月球星历表获得Tm时刻的月球位置为rM＝[xM,yM,zM]。
再入时刻T0对应的格林尼治时间为GST，则地心固连系至地心惯性系旋转矩阵MFI为
Figure BDA0003574284720000131
则地心惯性系下再入点处探测器的状态为
Xi＝[(MFI·re T)T,(MFI·ve T)T] (53)
将探测器由再入点状态Xi逆推至近月点得到实际转移时长为
Figure BDA0003574284720000132
的可行月地转移轨道。
步骤2.2：根据实际转移时长与预设的转移时长之差，调整预设再入速度大小，得到实际转移时长为预设转移时长的可行月地转移轨道。
记再入速度大小为vre时，实际再入点至近月点转移时长为
Figure BDA0003574284720000133
若
Figure BDA0003574284720000134
说明再入速度设置偏大，探测器提前到达近月点，需减小再入速度大小vre；若
Figure BDA0003574284720000135
说明再入速度设置偏小，探测器滞后到达近月点，需增大再入速度大小vre。
通过迭代调整再入速度ve的大小vre使得近月点至再入点转移时长刚好为ΔT。确定每次迭代中再入速度大小vre的调整量(步长)Δvre是关键，步长过大或过小都会导致转移时长存在偏差。所以，提出一种根据实际转移时长与期望转移时长偏差确定调整步长的方法。
假设迭代中月地实际转移时长与预设转移时长的偏差为
Figure BDA0003574284720000136
则下次迭代时再入速度大小的调整步长
Δvre＝δt/100(δt单位为天) (54)
再入速度调整为
vre＝vre+Δvre (55)
当δt＜0时，Δvre＜0，即再入速度调小，当δt＞0时，Δvre＞0，即再入速度调大。据此进行再入速度大小迭代，直至满足转移时长精度要求
Figure BDA0003574284720000137
其中，εT为转移时长精度误差限，根据实际工程需要进行设置。
至此，得到转移时长固定为ΔT＝3.0天的再入点至近月点的可行月地转移轨道。
步骤2.3：采用黄金分割法，以近月点高度为指标，确定特定日期一天内的最优月地转移轨道。
地球一天自转一周，所以T0在[D1,D1+1]的一天内变化时，惯性系下可行月地转移的轨道方向会旋转一周，所以步骤2.2得到可行月地转移轨道的近月点高度在一天内近似周期变化，最大近月点高度甚至可达数万km。所以要寻找[D1,D1+1]一天内近月点高度最低的月地转移轨道对应的再入点时间T0。一天内转移时长固定为ΔT的月地转移轨道近月点高度变化曲线呈单峰+单谷趋势，所以考虑采用黄金分割法求得一天内最小近月点高度对应的再入点时间T0。
对于最优月地转移轨道求解问题，黄金分割法实现的具体步骤为：在时间区间[a,b](a＝D1,b＝D1+1)内取点a1,a2把区间分为三段，其中
a1＝a+0.382(b-a)，a2＝a+0.618(b-a) (57)
记F(x)为再入点时间为x时可行月地转移轨道对应的近月点高度，若F(a1)＜F(a2)，令b＝a2，a不变；若F(a1)＞F(a2)，令a＝a1，b不变。实施上述迭代过程，直至满足最优再入点时间精度
b-a＜εrp (58)
其中，εrp为最优再入时刻精度误差限。至此，转移时长为ΔT时的使得近月点高度最小的最优再入点时间T0＝a，从而得到[D1,D1+1]一天内近月点高度最小的最优月地转移轨道。图2给出了2030年10月3日一天中的最优月地转移轨道，对应的再入点时间为2030/10/0322:26:01，再入速度大小vre＝10.6541km/s。
步骤三：预设任务搜索日期区间，遍历得到区间内地月引力模型下每天的最优月地转移轨道对应的近月点高度，根据设定的近月点高度偏差限确定月地转移窗口。
预设任务搜索日期区间为[D2,D3]＝[2462502.5,2462867.5](2030/01/01 00:00:00–2031/01/01 00:00:00)，采用步骤一和步骤二遍历得到每一天τ∈[D2,D3]中的最优月地转移轨道及其对应的近月点高度H。
经过轨道修正实现月球低轨捕获所需的燃耗与修正前最优月地转移轨道的近月点高度相关。若最优月地转移轨道近月点高度过大，会导致修正燃耗过大，超出探测器能力范围或任务约束。所以，根据每天最小的近月点半径来判断窗口的存在性。
预设近月点高度偏差限为χ＝50000km，若日期τ的最优月地转移轨道近月点高度H＜χ，则判定日期τ在月地转移窗口中，若H＞χ，则判定日期τ不在月地转移窗口中。从而，可靠高效地得到日期区间[D2,D3]内的月地转移窗口Γ为每一个在转移窗口中的日期τ的集合。如图3所示，2030年一年中月地转移窗口Γ共有74天。
步骤四：预设近月端出发倾角和高度约束，基于得到的月地转移窗口，修正得到高精度摄动模型下月球影响球边界的探测器状态以及对应的月球影响球内高精度转移轨道。
步骤4.1：以步骤三得到的最优月地转移轨道的再入状态为初值，在高精度动力学模型下逆推确定近月点处参数。
在月地转移窗口中选取T0＝2462778.43474(2030/10/03 22:26:01)，基于步骤一、步骤二、步骤三得到最优月地转移轨道，采用微分修正设计符合月球端出发倾角和高度约束的月球影响球内高精度转移轨道。
以得到的最优月地转移轨道的再入点处状态为初值，在高精度动力学下进行精确转移轨道设计。建立考虑日地月中心引力、日月非球形摄动、太阳光压的高精度动力学模型为
Figure BDA0003574284720000151
其中，μS为太阳引力常数，rS为根据星历文件得到的地心惯性系下太阳的位置矢量，aNE为地球非球形摄动加速度，aNM为月球非球形摄动加速度，asrp为太阳光压摄动加速度。
将最优转移轨道作为设计初值，将初值Xi在高精度动力学式(59)中逆向积分后，得到地心惯性系下的近月点状态并将其转化到月心赤道惯性坐标系中，表示为Xf＝[rf,vf]，以轨道角动量h、半长轴a、偏心率e、半通径p、真近点角f、偏近点角H为中间变量，计算距近月点时间τ，轨道高度rpm，轨道倾角im。
记标量与矢量对应关系为
a＝|a|,e＝|e|,h＝|h|,rf＝|rf|,vf＝|vf| (60)
则通过近月点状态rf和vf推导得到近月点处以下式子
Figure BDA0003574284720000152
记单位向量
n＝[0,0,1] (62)
则距近月点时间τ，轨道高度rpm，轨道倾角im可以通过式(63)获得：
Figure BDA0003574284720000153
步骤4.2：预设月地转移轨道的月球端出发近月点高度约束和倾角约束，采用微分修正确定月球影响边界处探测器的状态和月球影响球内高精度转移轨道。
预设月地转移轨道的月球端出发近月点高度为
Figure BDA0003574284720000161
倾角为
Figure BDA0003574284720000162
将距近月点时间τ、轨道高度rpm、轨道倾角im这三个量与预设值之差作为目标约束变量(标称距近月点时间为0)，故月地转移轨道求解可转化为如下非线性方程组求解：
[rpm-r′pm,im-i′m,τ]＝0 (64)
将初值Xi在高精度动力学式(59)中逆向积分后，得到地心惯性系下的月球影响球边缘处探测器的状态并将其转换为轨道根数形式，选取修正量为月球影响球边缘处初始状态的半长轴a、升交点赤经Ω、近心点幅角w，通过数值方法分析修正量对目标约束变量的影响。分别给每个修正量初始状态一个小的扰动δ，经过逆向积分到末端近月点后得到目标约束变量对于扰动δ的响应，即可得到目标约束变量相对修正量初始状态的敏度矩阵K1为
Figure BDA0003574284720000163
则修正量的微分与目标约束变量误差之间的关系为
Figure BDA0003574284720000164
迭代数次即可得到满足非线性方程组(64)的月球影响球边缘处探测器的初始状态a、Ω、w的值。从而，得到使得近月点倾角和高度满足预设约束的月球影响球边缘处的探测器状态为XSOI＝[rSOI,vSOI]。以此状态为初值，对月球影响球内轨道进行高精度递推，得到近月点状态XM＝[rM,vM]。
探测器从高度为rpm的月球圆轨道施加离轨脉冲进入高精度月地转移轨道时，所需施加的速度增量为
Figure BDA0003574284720000165
其中，RM为月球平均半径。
步骤五：采用微分修正，在再入点前一天处施加速度修正使得探测器的转移轨道在月球影响球边界处位置精确拼接，得到月球影响球边界至再入点的高精度转移轨道。
修正后的月球影响球边缘处探测器状态与修正前状态相比，位置和速度均有小幅度变化，所以需对再入点-月球影响球边缘段轨道进行修正以对准月球影响球边缘处探测器状态。
采用再入点前一天施加修正脉冲的方式使得月球影响球边缘处探测器位置对准，即整条月地转移轨道分为(1)近月点-月球影响球边缘段、(2)月球影响球边缘-再入前一天段、(3)再入前一天-再入点段这三段轨道，施加的脉冲共有三次，即除了每两段之间的连接脉冲之外，还包括式(67)的近月点处离轨脉冲。
第三段轨道递推以再入点初值Xi在高精度动力学模型(59)下逆推得到再入前一天处探测器的状态为Xrr＝[rrr,vrr]，其中rrr＝[xrr,yrr,zrr]，vrr＝[vxrr,vyrr,vzrr]。选择终端约束为步骤四得到的第一段轨道月球影响球边缘处探测器的位置rSOI＝[xSOI,ySOI,zSOI]。分别给再入点前一天处探测器速度vrr的每个分量一个小的扰动δ，经过逆向积分后得到月球影响球边缘处探测器位置受扰动的响应，即可得到月球影响球边缘处探测器位置相对再入点前一天处速度扰动的敏度矩阵K2为
Figure BDA0003574284720000171
从而得到速度修正量的微分和月球影响球边缘处探测器位置误差的关系为
Figure BDA0003574284720000172
迭代数次得到使得月球影响球处探测器位置精确对准的再入点前一天处探测器速度增量Δvrr＝[Δvxrr,Δvyrr,Δvzrr]＝[-5.2,-1.2,5.8]m/s的值，从而，第二段轨道中再入点前一天处初始速度为v′rr＝vrr+Δvrr。相应的第二段转移轨道中月球影响球边缘处探测器的末端速度v′SOI，则为了与第三段转移轨道中探测器的初始状态对应，月球影响球边缘处探测器应施加的速度增量为ΔvSOI＝vSOI-v′SOI＝[4.1,1.6,-2.7]m/s。
至此，得到月地转移窗口Γ内再入点时刻为T0时的高精度月地转移轨道，探测器所需施加的总速度增量为
|Δvall|＝|ΔvM|+|ΔvSOI|+|Δvrr|＝888.2m/s (70)
还包括步骤六：结合步骤四得到的月球影响球内高精度转移轨道和步骤五得到的月球影响球边界至再入点的高精度转移轨道，得到考虑月球端出发倾角和高度约束、地球端大气再入角和再入航程约束以及着陆场约束的高精度月地转移轨道，根据得到的高精度月地转移轨道对探测器施加相应机动，实现高效率、高精度、高可靠性的月地转移。
以上所述的具体描述，对发明的目的、技术方案和优点进行了进一步详细说明。需要理解的是，以上所述仅为本发明的具体实施示例，用于解释本发明，并不用于限定本发明的保护范围。凡在本发明的精神和原则之内所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (7)

1.一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，其特征在于：包括如下步骤，
步骤一：根据地面着陆场经纬度约束、再入点高度约束、再入轨迹航程和倾角约束、预设再入速度大小确定地心固连系下再入点处探测器的状态，提高再入舱着陆的可靠性；
步骤二：预设再入日期和月地转移时长，确定地月引力模型下当天近月点高度最小的最优月地转移轨道；
步骤三：预设任务搜索日期区间，遍历得到区间内地月引力模型下每天的最优月地转移轨道对应的近月点高度，根据预设的近月点高度偏差限确定月地转移窗口，保证窗口确定的可靠性和高效性；
步骤四：预设近月端出发倾角和高度约束，基于得到的月地转移窗口，修正得到高精度摄动模型下月球影响球边界的探测器状态以及对应的月球影响球内高精度转移轨道；
步骤五：采用微分修正，在再入点前施加速度修正使得探测器的转移轨道在月球影响球边界处位置精确拼接，得到月球影响球边界至再入点的高精度转移轨道。
2.如权利要求1所述的一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，其特征在于：还包括步骤六，结合步骤四得到的月球影响球内高精度转移轨道和步骤五得到的月球影响球边界至再入点的高精度转移轨道，得到考虑月球端出发倾角和高度约束、地球端大气再入角和再入航程约束以及着陆场约束的高精度月地转移轨道，根据得到的高精度月地转移轨道对探测器施加相应机动，实现高效率、高精度、高可靠性的月地转移。
3.如权利要求1或2所述的一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，其特征在于：步骤1实现方法为，
步骤1.1：根据地面着陆场纬度、再入轨迹倾角和航程确定再入点的纬度幅角；
预设地面着陆场纬度为
Figure FDA0003574284710000013
再入轨迹倾角为ir，则升轨再入时，根据球面三角形正弦定理，着陆场的纬度幅角UF为
Figure FDA0003574284710000011
降轨再入时，着陆场的纬度幅角UF为
Figure FDA0003574284710000012
预设再入轨迹的航程为Sr，则再入轨迹的航程幅角Ur为
Figure FDA0003574284710000021
其中，Re为地球平均半径；
从而得到再入点的纬度幅角UP为
UP＝UF-Ur (4)
步骤1.2：根据着陆点经度和再入点纬度幅角确定再入点的飞行方位角和再入点经纬度；
预设地面着陆场经度为λF，根据再入点的纬度幅角UP和再入轨迹轨道倾角ir，由球面三角形纳皮尔法则确定再入点到升交点的经度差ΔλP为
Figure FDA0003574284710000022
从而再入点的飞行方位角Q为
Figure FDA0003574284710000023
基于着陆场纬度幅角UF计算着陆场到升交点的经度差ΔλF为
Figure FDA0003574284710000024
则确定再入点的经度λP为
λP＝λF-ΔλF+ΔλP (8)
再入点的纬度
Figure FDA0003574284710000025
为
Figure FDA0003574284710000026
步骤1.3：根据再入点高度约束、再入点经纬度和飞行方位角确定再入点处探测器在地心固连系中的状态；
预设再入点高度为H，则再入点的地心距R为
R＝H+Re (10)
从而确定再入点在地心固连系下的位置矢量为
Figure FDA0003574284710000027
定义再入坐标系原点位于探测器上，Y轴指向天顶方向，X轴在再入轨迹面内垂直Y轴指向飞行方向，Z轴与X轴和Y轴组成右手坐标系；记γ为再入角，vre为预设再入速度大小，则再入坐标系中探测器在再入点处的速度矢量vre为
vre＝[vrecosγ,vresinγ,0] (12)
再入坐标系到地心固连系的旋转矩阵MEF为
Figure FDA0003574284710000031
则地心固连系下再入点处探测器的速度矢量为
ve＝(MEF·vre T)T (14)
其中，上标“T”为转置符号；
进而，地心固连系下再入点处探测器的状态为
Xe＝[re,ve] (15)
上述基于地面着陆场的再入点状态确定保证再入舱着陆的可靠性。
4.如权利要求3所述的一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，其特征在于：步骤2实现方法为，
步骤2.1：选择预设再入日期一天内的选定时刻，根据预设的月地转移时长，在地月引力模型下将探测器状态从再入点逆推至近月点得到选定时刻的可行月地转移轨道；
由于月球公转速度相比于地球自转速度更慢，所以随着地球自转，每天之中存在一条能到达月球附近的最优月地转移轨道；
由于探测器在月地转移过程中受到的引力作用主要来自地球和月球，故在转移窗口确定前采用地月引力模型进行转移轨道计算，在保证可靠性的前提下大大提高计算效率；建立地月引力模型为
Figure FDA0003574284710000032
其中，
Figure FDA0003574284710000033
为探测器在地月引力模型中的加速度，r为探测器在地心惯性系下的位置矢量，rM为月球在地心惯性系下的位置矢量，μE为地球引力常数，μM为月球引力常数；
预设再入日期时段的儒略日范围为D1零时至D1+1零时，选择再入点时刻为T0∈[D1,D1+1]；月地转移时长预设为ΔT，则近月点时刻为Tm＝T0-ΔT，通过月球星历表获得Tm时刻的月球位置为rM＝[xM,yM,zM]；
再入时刻T0对应的格林尼治时间为GST，则地心固连系至地心惯性系旋转矩阵MFI为
Figure FDA0003574284710000041
则地心惯性系下再入点处探测器的状态为
Xi＝[(MFI·re T)T,(MFI·ve T)T] (18)
将探测器由再入点状态Xi逆推至近月点得到实际转移时长为
Figure FDA0003574284710000042
的最优月地转移轨道；
步骤2.2：根据实际转移时长与预设的转移时长之差，调整预设再入速度大小，得到实际转移时长为预设转移时长的可行月地转移轨道；
记再入速度大小为vre时，实际再入点至近月点转移时长为
Figure FDA0003574284710000043
若
Figure FDA0003574284710000044
说明再入速度设置偏大，探测器提前到达近月点，需减小再入速度大小vre；若
Figure FDA0003574284710000045
说明再入速度设置偏小，探测器滞后到达近月点，需增大再入速度大小vre；
通过迭代调整再入速度ve的大小vre使得近月点至再入点转移时长刚好为ΔT；确定每次迭代中再入速度大小vre的调整量Δvre是关键，步长过大或过小都会导致转移时长存在偏差；所以，提出一种根据实际转移时长与期望转移时长偏差确定调整步长的方法；
迭代中月地实际转移时长与预设转移时长的偏差为
Figure FDA0003574284710000046
则下次迭代时再入速度大小的调整步长
Δvre＝δt/100(δt单位为天) (19)
再入速度调整为
vre＝vre+Δvre (20)
当δt＜0时，Δvre＜0，即再入速度调小，当δt＞0时，Δvre＞0，即再入速度调大；据此进行再入速度大小迭代，直至满足转移时长精度要求
Figure FDA0003574284710000047
其中，εT为转移时长精度误差限，根据实际工程需要进行设置；
至此，得到转移时长固定为ΔT的再入点至近月点的可行月地转移轨道；
步骤2.3：采用黄金分割法，以近月点高度为指标，确定预设日期一天内的最优月地转移轨道；
地球一天自转一周，所以T0在[D1,D1+1]的一天内变化时，惯性系下可行月地转移的轨道方向会旋转一周，所以步骤2.2得到可行月地转移轨道的近月点高度在一天内近似周期变化，最大近月点高度甚至可达数万km；所以要寻找[D1,D1+1]一天内近月点高度最低的月地转移轨道对应的再入点时间T0；一天内转移时长固定为ΔT的月地转移轨道近月点高度变化曲线呈单峰+单谷趋势，所以考虑采用黄金分割法求得一天内最小近月点高度对应的再入点时间T0；
对于最优月地转移轨道求解问题，黄金分割法实现的具体步骤为：在时间区间[a,b](a＝D1,b＝D1+1)内取点a1,a2把区间分为三段，其中
a1＝a+0.382(b-a)，a2＝a+0.618(b-a) (22)
记F(x)为再入点时间为x时可行月地转移轨道对应的近月点高度，若F(a1)＜F(a2)，令b＝a2，a不变；若F(a1)＞F(a2)，令a＝a1，b不变；实施上述迭代过程，直至满足最优再入点时间精度
b-a＜εrp (23)
其中，εrp为最优再入时刻精度误差限；至此，转移时长为ΔT时的使得近月点高度最小的最优再入点时间T0＝a，从而得到[D1,D1+1]一天内近月点高度最小的最优月地转移轨道。
5.如权利要求4所述的一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，其特征在于：步骤3实现方法为，
预设任务搜索日期区间为[D2,D3](D3-D2＞1)，采用步骤一和步骤二遍历得到每一天τ∈[D2,D3]中的最优月地转移轨道及其对应的近月点高度H；
经过轨道修正实现月球低轨捕获所需的燃耗与修正前最优月地转移轨道的近月点高度相关；若最优月地转移轨道近月点高度过大，会导致修正燃耗过大，超出探测器能力范围或任务约束；所以，根据每天最小的近月点半径来判断窗口的存在性；
预设近月点高度偏差限为χ，若日期τ的最优月地转移轨道近月点高度H＜χ，则判定日期τ在月地转移窗口中，若H＞χ，则判定日期τ不在月地转移窗口中；从而，可靠高效地得到日期区间[D2,D3]内的月地转移窗口Γ为每一个在转移窗口中的日期τ的集合。
6.如权利要求5所述的一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，其特征在于：步骤4.1：以步骤三得到的最优月地转移轨道的再入状态为初值，在高精度动力学模型下逆推确定近月点处参数；
在月地转移窗口中任选T0∈Γ，基于步骤一、步骤二、步骤三得到最优月地转移轨道，采用微分修正设计符合月球端出发倾角和高度约束的月球影响球内高精度转移轨道；
以得到的最优月地转移轨道的再入点处状态为初值，在高精度动力学下进行精确转移轨道设计；建立考虑日地月中心引力、日月非球形摄动、太阳光压的高精度动力学模型为
Figure FDA0003574284710000061
其中，μS为太阳引力常数，rS为根据星历文件得到的地心惯性系下太阳的位置矢量，aNE为地球非球形摄动加速度，aNM为月球非球形摄动加速度，asrp为太阳光压摄动加速度；
将最优转移轨道作为设计初值，将初值Xi在高精度动力学式(24)中逆向积分后，得到地心惯性系下的近月点状态并将其转化到月心赤道惯性坐标系中，表示为Xf＝[rf,vf]，以轨道角动量h、半长轴a、偏心率e、半通径p、真近点角f、偏近点角H为中间变量，计算距近月点时间τ，轨道高度rpm，轨道倾角im；
记标量与矢量对应关系为
a＝|a|,e＝|e|,h＝|h|,rf＝|rf|,vf＝|vf| (25)
则通过近月点状态rf和vf推导得到近月点处以下式子
Figure FDA0003574284710000062
记单位向量
n＝[0,0,1] (27)
则距近月点时间τ，轨道高度rpm，轨道倾角im可以通过式(28)获得：
Figure FDA0003574284710000063
步骤4.2：预设月地转移轨道的月球端出发近月点高度约束和倾角约束，采用微分修正确定月球影响边界处探测器的状态和月球影响球内高精度转移轨道；
预设月地转移轨道的月球端出发近月点高度为
Figure FDA0003574284710000064
倾角为
Figure FDA0003574284710000065
将距近月点时间τ、轨道高度rpm、轨道倾角im三个量与预设值之差作为目标约束变量，故月地转移轨道求解转化为如下非线性方程组求解：
[rpm-r′pm,im-i′m,τ]＝0 (29)
将初值Xi在高精度动力学式(24)中逆向积分后，得到地心惯性系下的月球影响球边缘处探测器的状态并将其转换为轨道根数形式，选取修正量为月球影响球边缘处初始状态的半长轴a、升交点赤经Ω、近心点幅角w，通过数值方法分析修正量对目标约束变量的影响；分别给每个修正量初始状态一个小的扰动δ，经过逆向积分到末端近月点后得到目标约束变量对于扰动δ的响应，即得到目标约束变量相对修正量初始状态的敏度矩阵K1为
Figure FDA0003574284710000071
则修正量的微分与目标约束变量误差之间的关系为
Figure FDA0003574284710000072
迭代数次即可得到满足非线性方程组(29)的月球影响球边缘处探测器的初始状态a、Ω、w的值；从而得到近月点倾角和高度满足预设约束的月球影响球边缘处的探测器状态为XSOI＝[rSOI,vSOI]；以此状态为初值，对月球影响球内轨道进行高精度递推，得到近月点状态XM＝[rM,vM]；
探测器从高度为rpm的月球圆轨道施加离轨脉冲进入高精度月地转移轨道时，所需施加的速度增量为
Figure FDA0003574284710000073
其中，RM为月球平均半径。
7.如权利要求6所述的一种考虑再入约束的月地转移窗口及精确转移轨道确定方法，其特征在于：步骤5实现方法为，
修正后的月球影响球边缘处探测器状态与修正前状态相比，位置和速度均有小幅度变化，所以需对再入点-月球影响球边缘段轨道进行修正以对准月球影响球边缘处探测器状态；
采用再入点前一天施加修正脉冲的方式使得月球影响球边缘处探测器位置对准，即整条月地转移轨道分为(1)近月点-月球影响球边缘段、(2)月球影响球边缘-再入前一天段、(3)再入前一天-再入点段这三段轨道，施加的脉冲共有三次，即除了每两段之间的连接脉冲之外，还包括式(32)的近月点处离轨脉冲；
第三段轨道递推以再入点初值Xi在高精度动力学模型(24)下逆推得到再入前一天处探测器的状态为Xrr＝[rrr,vrr]，其中rrr＝[xrr,yrr,zrr]，vrr＝[vxrr,vyrr,vzrr]；选择终端约束为步骤四得到的第一段轨道月球影响球边缘处探测器的位置rSOI＝[xSOI,ySOI,zSOI]；分别给再入点前一天处探测器速度vrr的每个分量一个小的扰动δ，经过逆向积分后得到月球影响球边缘处探测器位置受扰动的响应，即可得到月球影响球边缘处探测器位置相对再入点前一天处速度扰动的敏度矩阵K2为
Figure FDA0003574284710000081
从而得到速度修正量的微分和月球影响球边缘处探测器位置误差的关系为
Figure FDA0003574284710000082
迭代数次得到使得月球影响球处探测器位置精确对准的再入点前一天处探测器速度增量Δvrr＝[Δvxrr,Δvyrr,Δvzrr]的值，从而，第二段轨道中再入点前一天处初始速度为v′rr＝vrr+Δvrr；相应的第二段转移轨道中月球影响球边缘处探测器的末端速度v′SOI，则为了与第三段转移轨道中探测器的初始状态对应，月球影响球边缘处探测器应施加的速度增量为ΔvSOI＝vSOI-v′SOI；
至此，得到月地转移窗口Γ内再入点时刻为T0时的高精度月地转移轨道，探测器所需施加的总速度增量为
|Δvall|＝|ΔvM|+|ΔvSOI|+|Δvrr| (35)。
CN202210335695.3A 2022-03-31 2022-03-31 考虑再入约束的月地转移窗口及精确转移轨道确定方法 Active CN114684389B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210335695.3A CN114684389B (zh) 2022-03-31 2022-03-31 考虑再入约束的月地转移窗口及精确转移轨道确定方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210335695.3A CN114684389B (zh) 2022-03-31 2022-03-31 考虑再入约束的月地转移窗口及精确转移轨道确定方法

Publications (2)

Publication Number Publication Date
CN114684389A true CN114684389A (zh) 2022-07-01
CN114684389B CN114684389B (zh) 2024-07-12

Family

ID=82140294

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210335695.3A Active CN114684389B (zh) 2022-03-31 2022-03-31 考虑再入约束的月地转移窗口及精确转移轨道确定方法

Country Status (1)

Country Link
CN (1) CN114684389B (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117228007A (zh) * 2023-11-13 2023-12-15 北京控制工程研究所 Sso轨道航天器主动离轨策略的快速确定方法及装置

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109344449A (zh) * 2018-09-07 2019-02-15 北京空间技术研制试验中心 航天器月地转移轨道逆向设计方法
US10696423B1 (en) * 2015-12-31 2020-06-30 United States Of America As Represented By The Administrator Of Nasa Method for transferring a spacecraft from geosynchronous transfer orbit to lunar orbit
CN113310496A (zh) * 2021-05-08 2021-08-27 北京航天飞行控制中心 一种确定月地转移轨道的方法及装置

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10696423B1 (en) * 2015-12-31 2020-06-30 United States Of America As Represented By The Administrator Of Nasa Method for transferring a spacecraft from geosynchronous transfer orbit to lunar orbit
CN109344449A (zh) * 2018-09-07 2019-02-15 北京空间技术研制试验中心 航天器月地转移轨道逆向设计方法
CN113310496A (zh) * 2021-05-08 2021-08-27 北京航天飞行控制中心 一种确定月地转移轨道的方法及装置

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
周亮;胡军;: "基于双曲线B平面参数的快速微分修正月地返回轨道精确设计", 航天控制, no. 06, 15 December 2012 (2012-12-15) *
孟占峰;高珊;汪中生;周文艳;: "月地高速再入返回任务轨道设计与飞行评价", 中国科学:技术科学, no. 03, 20 March 2015 (2015-03-20) *
本立言;严玲玲;谢祥华;张锐;王国际;: "直接再入大气的月地转移轨道设计", 北京航空航天大学学报, no. 02, 29 February 2020 (2020-02-29) *
胡军;周亮;: "基于显式制导的月地返回轨道中途修正研究", 中国空间科学技术, no. 05, 25 October 2013 (2013-10-25) *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117228007A (zh) * 2023-11-13 2023-12-15 北京控制工程研究所 Sso轨道航天器主动离轨策略的快速确定方法及装置
CN117228007B (zh) * 2023-11-13 2024-01-19 北京控制工程研究所 Sso轨道航天器主动离轨策略的快速确定方法及装置

Also Published As

Publication number Publication date
CN114684389B (zh) 2024-07-12

Similar Documents

Publication Publication Date Title
CN102878995B (zh) 一种静止轨道卫星自主导航方法
US10046869B2 (en) Inertial sensing augmentation for navigation of spacecraft
CN107589756B (zh) 一种奔月卫星编队初始化方法
CN108562295B (zh) 一种基于同步卫星二体模型的三站时差定轨方法
CN102495950B (zh) 一种适用于太阳同步轨道的倾角偏置量获取方法
CN111268177B (zh) 一种静止轨道卫星分散式闭环自主位置保持控制方法
CN113740887A (zh) 一种卫星注入轨道外推及卫星理论轨道确定方法
CN112713922A (zh) 一种多波束通讯卫星的可见性快速预报算法
CN114889849A (zh) 一种航天器多约束飞越最优速度的估计方法
CN114684389A (zh) 考虑再入约束的月地转移窗口及精确转移轨道确定方法
CN109606739B (zh) 一种探测器地月转移轨道修正方法及装置
CN113447043B (zh) 一种基于gnss的卫星天文导航***误差自主标定方法及***
CN108082538B (zh) 一种考虑始末约束的多体***低能量捕获轨道方法
Zhu et al. A high-accuracy SINS/CNS integrated navigation scheme based on overall optimal correction
CN114063645B (zh) 基于偏心率倾角矢量的倾斜绕飞保持控制效果评估方法
CN115933743A (zh) 一种针对目标窗***会的轨迹规划方法
CN114802818A (zh) 晨昏轨道卫星及其对日姿态计算方法、导引方法
CN111272336B (zh) 一种基于gnss观测实现大型低轨航天器质心位移估算的方法
CN111090830A (zh) 一种高轨非合作目标在轨光压辨识方法
Xie et al. Autonomous guidance, navigation, and control of spacecraft
CN116300417B (zh) 一种大型分布式空间望远镜主次镜编队控制方法
CN116817932B (zh) 一种轨道维持、轨道确定和重力场实时测绘的一体化方法
Ibrahim Attitude and orbit control of small satellites for autonomous terrestrial target tracking
Shishov Determination of spacecraft and Phobos parameters of motion in the Phobos-Grunt project
CN115113638A (zh) 一种燃料最优主动漂移三维成像轨道控制方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant