CN111931131B - 一种行星探测器动力着陆段在线轨迹规划方法及系统 - Google Patents

一种行星探测器动力着陆段在线轨迹规划方法及系统 Download PDF

Info

Publication number
CN111931131B
CN111931131B CN202010805032.4A CN202010805032A CN111931131B CN 111931131 B CN111931131 B CN 111931131B CN 202010805032 A CN202010805032 A CN 202010805032A CN 111931131 B CN111931131 B CN 111931131B
Authority
CN
China
Prior art keywords
landing
trajectory
track
detector
elements
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.)
Active
Application number
CN202010805032.4A
Other languages
English (en)
Other versions
CN111931131A (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.)
Beihang University
Original Assignee
Beihang University
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 Beihang University filed Critical Beihang University
Priority to CN202010805032.4A priority Critical patent/CN111931131B/zh
Publication of CN111931131A publication Critical patent/CN111931131A/zh
Application granted granted Critical
Publication of CN111931131B publication Critical patent/CN111931131B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • 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
    • B64G1/242—Orbits and trajectories

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • General Physics & Mathematics (AREA)
  • Remote Sensing (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Probability & Statistics with Applications (AREA)
  • General Engineering & Computer Science (AREA)
  • Evolutionary Biology (AREA)
  • Algebra (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Operations Research (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Chemical & Material Sciences (AREA)
  • Combustion & Propulsion (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Feedback Control In General (AREA)

Abstract

本发明公开了一种行星探测器动力着陆段在线轨迹规划方法及系统。该方法利用切比雪夫伪谱法进行离散，化最优控制问题为非线性规划问题，利用逐次凸优化方法对非线性规划的求解化为逐次求解的二阶锥规划问题，进而利用原始对偶内点法作为数值求解算法迭代求解此二阶锥规划即可得到着陆过程的轨迹，同时确定终端时间。本发明的方法及系统，能够增强轨迹规划的速度和精度，适用于进行在线轨迹规划。

Description

一种行星探测器动力着陆段在线轨迹规划方法及系统
技术领域
本发明涉及轨迹规划技术领域，特别是涉及一种行星探测器动力着陆段在线轨迹规划方法及系统。
背景技术
在一般的行星探测器着陆过程中，通常会以燃料最优作为指标设计着陆轨迹引导探测器安全着陆。由于伞降过程的不确定性，使得动力着陆过程的初始条件无法事先确定。因此根据着陆阶段的实际情况在线设计轨迹是一种理想方式。在线设计着陆轨迹需要轨迹规划算法具有足够快的规划速度以及规划精度。但是着陆过程中存在着诸多实际的物理约束的存在，如能够提供的推力有限，对地的传感器的观测角度不确定等，使得优化问题的求解变得复杂。另外着陆过程所需的时长无法事先确定也会影响求解最优轨迹的速度。
现在很多数值求解最优轨迹的方法都采用直接配点法中的等间距离散法进行离散，利用SQP(逐次二次规划)进行求解所得到的NLP(非线性规划)。然而，SQP算法不具备多项式时间特性，并不适合进行在线应用。除了数值解法对求解速度的影响，为了保证离散精度，需要大量的离散节点。而大量的离散节点意味着一个大规模的优化问题，也会影响求解速度。此外，传统的确定终端时间的方法是先固定终端时间，利用数值优化方法求解此优化问题，再通过线搜索等技巧确定终端时间。因此需要求解十几个相似的优化问题才可以得到最优的终端时间，从而确定最优着陆轨迹。
因此，传统的轨迹规划方法不适合进行在线的着陆轨迹规划。如何选择少节点高精度的离散方法，利用具有多项式时间特性的数值优化方法，在求解最优终端时间的同时得到最优轨迹是提高动力着陆段轨迹优化速度的关键问题。
发明内容
本发明的目的是提供一种行星探测器动力着陆段在线轨迹规划方法及系统，将不确定的终端时间视为优化变量之一，并将其应用于着陆过程的轨迹规划中，能够增强轨迹规划的速度和精度。
为实现上述目的，本发明提供了如下方案：
一种在线轨迹规划方法，包括：
获取探测器的初始状态；所述初始状态包括着陆过程开始时探测器的质量、位置和速度；
根据所述探测器的初始状态建立着陆轨迹优化模型；所述着陆轨迹优化模型以最小化燃料消耗量为目标，并根据行星探测器着陆过程的限制因素确定约束条件；
采用切比雪夫伪谱法将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型；所述着陆轨迹非线性规划模型中的变量包括状态变量、控制变量和终端时间，所述状态变量与探测器的质量、位置和速度有关，所述控制变量与探测器发动机的推力有关，所述终端时间为着陆过程结束的时刻；
获取参考轨迹，根据所述参考轨迹采用凸优化方法将所述着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型；所述参考轨迹包括参考状态值、参考控制值和参考终端时间；
采用原始对偶内点法对所述着陆轨迹二阶锥规划模型进行求解，得到期望轨迹；所述期望轨迹包括期望状态值、期望控制值和期望终端时间；
判断着陆轨迹与所述期望轨迹的差值的切比雪夫范数是否小于预设阈值；若小于所述预设阈值，则根据所述期望轨迹采用插值法确定着陆过程的轨迹，然后输出所述着陆过程的轨迹以及所述期望终端时间；否则，根据所述期望轨迹对所述参考轨迹进行更新，并返回步骤“获取参考轨迹，根据所述参考轨迹采用凸优化方法将所述着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型”。
可选的，所述根据所述探测器的初始状态建立着陆轨迹优化模型，具体包括：
根据如下公式建立着陆轨迹优化模型的目标函数：
min J＝-m(tf)
根据如下公式建立着陆轨迹优化模型的约束条件：
Figure BDA0002628797740000021
Figure BDA0002628797740000022
fmin≤||T(t)||2≤fmax
Figure BDA0002628797740000031
m(0)＝m0,r(0)＝r0,v(0)＝v0
r(tf)＝rf,v(tf)＝vf
式中，J为着陆轨迹优化模型的目标函数，m(tf)为探测器在tf时刻的质量，
Figure BDA0002628797740000037
为r(t)的导数，r(t)为探测器在t时刻的位置，
Figure BDA0002628797740000033
为v(t)的导数，v(t)为探测器在t时刻的速度，T(t)为制动发动机在t时刻的推力，
Figure BDA0002628797740000034
为m(t)的导数，m(t)为探测器在t时刻的质量，g为重力加速度，α为发动机推力与燃料消耗量之间的系数，fmin为推力下界值，fmax为推力上界值，e1为第1个分量为1其余分量为0的单位列向量，e2为第2个分量为1其余分量为0的单位列向量，e3为第3个分量为1其余分量为0的单位列向量，γ为预设滑翔倾角，m0为着陆过程开始时探测器的质量，r0为着陆过程开始时探测器的位置，v0为着陆过程开始时探测器的速度，rf为探测器着陆点的位置，rf＝0，vf为探测器安全着陆所需的速度，vf＝0，tf为终端时间。
可选的，所述采用切比雪夫伪谱法将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型，具体包括：
根据探测器着陆过程的时域范围确定着陆轨迹中的离散节点；所述时域范围为探测器着陆过程的始端时间与终端时间形成的闭区间；
采用直接配点法确定在所述离散节点处的状态变量和控制变量；
根据所述状态变量和所述控制变量，将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型。
可选的，所述根据探测器着陆过程的时域范围确定着陆轨迹中的离散节点，具体包括：
将所述时域范围由[0,tf]转化到[-1,1]，得到
Figure BDA0002628797740000035
其中，t为实际时间，τ为伪谱时间；
获取着陆轨迹中的离散节点个数，并根据所述离散节点个数根据如下公式确定着陆轨迹中的离散节点：
Figure BDA0002628797740000036
式中，τk为第k个离散节点，N+1为离散节点总数。
可选的，所述根据所述状态变量和所述控制变量，将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型，具体包括：
根据如下公式确定所述着陆轨迹非线性规划模型的目标函数：
min J'＝-cTX
根据如下公式确定所述着陆轨迹非线性规划模型的约束条件：
Figure BDA0002628797740000041
Figure BDA0002628797740000042
||QiU||2≤mi TU
Figure BDA0002628797740000043
C1X＝X0,C2X＝Xf
其中，
X＝(R1 T R2 T R3 T V1 T V2 T V3 T ZT)T
U＝(U1 T U2 T U3 T ΣT)T
Figure BDA0002628797740000044
式中，J'为着陆轨迹非线性规划模型的目标函数，X为状态变量矩阵，U为控制变量矩阵，Rj为探测器位置向量，j＝1,2,3，R1为探测器水平横向位置向量，R2为探测器水平纵向位置向量，R3为探测器竖直方向位置向量，Vj为探测器速度向量，V1为探测器水平横向速度向量，V2为探测器水平纵向速度向量，V3为探测器竖直方向速度向量，Z为探测器质量对数向量，Uj为探测器第一控制向量，U1为探测器水平横向第一控制向量，U2为探测器水平纵向第一控制向量，U3为探测器竖直方向第一控制向量，Σ为探测器第二控制向量，rj(τN)为编号为N的离散节点处探测器的位置变量，vj(τN)为编号为N的离散节点处探测器的速度变量，z(τN)为编号为N的离散节点处探测器质量对数变量，uj(τN)为编号为N的离散节点处探测器第一控制变量，
Figure BDA0002628797740000045
u(t)为t时刻探测器第一控制变量，σ(τN)为编号为N的离散节点处探测器第二控制变量，
Figure BDA0002628797740000051
σ(t)为t时刻探测器第二控制变量，Γ(t)为松弛变量，||T(t)||≤Γ(t)，c为最后一个元素为1其余元素为0的7(N+1)维单位列向量，D为切比雪夫伪谱法对应的微分矩阵，I为7阶单位矩阵，A,B,G均为系数矩阵，
Figure BDA0002628797740000052
1(N+1)×1表示所有元素均为1的列向量，mi为除第3N+4+i个元素为1外，其余元素全为零的4(N+1)维列向量，ci为除第6N+7+i个元素为1外，其余元素全为零的7(N+1)维列向量，Tmin为最小飞行时间，Tmax为最大飞行时间，Qi为除元素(i,i)，(N+2+i,N+2+i)，(2N+3+i,2N+3+i)为1外，其余元素全为零的3(N+1)×4(N+1)阶矩阵，Pi为除元素(i,i)，(N+2+i,N+2+i)为1外，其余元素全为零的2(N+1)×7(N+1)阶矩阵，qi为除第2N+3+i个元素为1外，其余元素全为零的7(N+1)维列向量，C1为除元素(1,)1，(N+2,N+2)，(2N+3,2N+3)，(3N+4,3N+4)，(4N+5,4N+5)，(5N+6,5N+6)，(6N+7,6N+7)为1外，其余元素全为零的7×7(N+1)阶矩阵，C2为除元素(N+1,N+1)，(2(N+1),2(N+1))，(3(N+1),3(N+1))，(4(N+1),4(N+1))，(5(N+1),5(N+1))，(6(N+1),6(N+1))为1外，其余元素全为零的6×7(N+1)阶矩阵，X0为状态变量矩阵，X0＝(r0 T v0 T z0)T，z0为探测器的初始质量对数，Xf为终端状态变量矩阵，Xf＝(rf T vf T)T。
可选的，所述根据所述参考轨迹采用凸优化方法将所述着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型，具体包括：
根据如下公式确定着陆轨迹二阶锥规划模型的目标函数：
min J″＝-cTX
根据如下公式确定着陆轨迹二阶锥规划模型的约束条件：
Figure BDA0002628797740000053
Figure BDA0002628797740000054
Figure BDA0002628797740000055
||QiU||2≤mi TU
Figure BDA0002628797740000056
C1X＝X0,C2X＝Xf
式中，J″为着陆轨迹二阶锥规划模型的目标函数，{X(k),U(k),tf (k)}为参考轨迹，X(k)为第k次迭代后的状态变量矩阵，U(k)为第k次迭代后的控制变量矩阵，tf (k)为第k次迭代后的终端时间，{d0,d1,d2}为指数函数最小二乘逼近的系数。
本发明还提供一种在线轨迹规划系统，包括：
探测器的初始状态获取模块，用于获取探测器的初始状态；所述初始状态包括着陆过程开始时探测器的质量、位置和速度；
着陆轨迹优化模型确定模块，用于根据所述探测器的初始状态建立着陆轨迹优化模型；所述着陆轨迹优化模型以最小化燃料消耗量为目标，并根据行星探测器着陆过程的限制因素确定约束条件；
着陆轨迹非线性规划模型确定模块，用于采用切比雪夫伪谱法将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型；所述着陆轨迹非线性规划模型中的变量包括状态变量、控制变量和终端时间，所述状态变量与探测器的质量、位置和速度有关，所述控制变量与探测器发动机的推力有关，所述终端时间为着陆过程结束的时刻；
着陆轨迹二阶锥规划模型确定模块，用于获取参考轨迹，根据所述参考轨迹采用凸优化方法将所述着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型；所述参考轨迹包括参考状态值、参考控制值和参考终端时间；
轨迹计算模块，用于采用原始对偶内点法对所述着陆轨迹二阶锥规划模型进行求解，得到期望轨迹；所述期望轨迹包括期望状态值、期望控制值和期望终端时间；
判断模块，用于判断着陆轨迹与所述期望轨迹的差值的切比雪夫范数是否小于预设阈值；若小于所述预设阈值，则执行插值模块；否则，执行更新模块；
更新模块，用于根据所述期望轨迹对所述参考轨迹进行更新，然后执行所述着陆轨迹二阶锥规划模型确定模块；
插值模块，用于根据所述期望轨迹采用插值法确定着陆过程的轨迹，然后执行输出模块；
输出模块，用于输出所述着陆过程的轨迹以及所述期望终端时间。
可选的，所述着陆轨迹非线性规划模型确定模块，具体包括：
离散节点确定单元，用于根据探测器着陆过程的时域范围确定着陆轨迹中的离散节点；所述时域范围为探测器着陆过程的始端时间与终端时间形成的闭区间；
变量确定单元，用于采用直接配点法确定在所述离散节点处的状态变量和控制变量；
着陆轨迹非线性规划模型确定单元，用于根据所述状态变量和所述控制变量，将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型。
与现有技术相比，本发明的有益效果是：
本发明提出了一种在线轨迹规划方法及系统，利用切比雪夫伪谱法进行离散，化最优控制问题为NLP(非线性规划)问题，利用逐次凸优化方法对NLP的求解化为逐次求解的SOCP(二阶锥规划)问题，进而利用原始对偶内点法作为数值求解算法迭代求解此SOCP即可得到着陆过程的轨迹，同时确定终端时间，能够增强轨迹规划的速度和精度，适用于进行在线轨迹规划。
此外，本发明采用伪谱法的非等间距离散以及切比雪夫多项式的最佳逼近性质，使得使用切比雪夫伪谱法可以以更少的节点获得更高的离散精度；将不确定的终端时间视为优化变量可以在求解轨迹的同时求解终端时间，避免了大量求解过程；利用SOCP建模优化问题并利用原始对偶内点法作为数值求解算法可以显著提高求解精度和速度；作为数值求解精确性体现在使用伪谱法能够以更少的节点获得更高的离散精度；凸优化过程中根据变量的范围选择最小二乘逼近可以有效减少因逼近过程造成的误差，从而进一步提高精度。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。
图1为本发明实施例中在线轨迹规划方法流程图；
图2为本发明实施例中在线轨迹规划系统结构图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
本发明的目的是提供一种行星探测器动力着陆段在线轨迹规划方法及系统，将不确定的终端时间视为优化变量之一，并将其应用于着陆过程的轨迹规划中，能够增强轨迹规划的速度和精度。
为使本发明的上述目的、特征和优点能够更加明显易懂，下面结合附图和具体实施方式对本发明作进一步详细的说明。
实施例
图1为本发明实施例中在线轨迹规划方法流程图，如图1所示，一种在线轨迹规划方法，包括：
步骤101：获取探测器的初始状态；初始状态包括着陆过程开始时探测器的质量、位置和速度。
步骤102：根据探测器的初始状态建立着陆轨迹优化模型。着陆轨迹优化模型以最小化燃料消耗量为目标，并根据行星探测器着陆过程的限制因素确定约束条件。限制因素包括制动发动机产生的推力有限，以及防止飞行过程中与地面发生碰撞。
步骤102，具体包括：
根据如下公式建立着陆轨迹优化模型的目标函数：
min J＝-m(tf)
根据如下公式建立着陆轨迹优化模型的约束条件：
Figure BDA0002628797740000081
Figure BDA0002628797740000082
fmin≤||T(t)||2≤fmax
Figure BDA0002628797740000083
m(0)＝m0,r(0)＝r0,v(0)＝v0
r(tf)＝rf,v(tf)＝vf
式中，J为着陆轨迹优化模型的目标函数，m(tf)为探测器在tf时刻的质量，
Figure BDA0002628797740000091
为r(t)的导数，r(t)为探测器在t时刻的位置，
Figure BDA0002628797740000092
为v(t)的导数，v(t)为探测器在t时刻的速度，T(t)为制动发动机在t时刻的推力，
Figure BDA0002628797740000093
为m(t)的导数，m(t)为探测器在t时刻的质量，g为重力加速度，α为发动机推力与燃料消耗量之间的系数，fmin为推力下界值，fmax为推力上界值，e1为第1个分量为1其余分量为0的单位列向量，e2为第2个分量为1其余分量为0的单位列向量，e3为第3个分量为1其余分量为0的单位列向量，γ为预设滑翔倾角，可为4°，m0为着陆过程开始时探测器的质量，r0为着陆过程开始时探测器的位置，v0为着陆过程开始时探测器的速度，rf为探测器着陆点的位置，rf＝0，vf为探测器安全着陆所需的速度，vf＝0，tf为终端时间。
步骤103：采用切比雪夫伪谱法将着陆轨迹优化模型转换为着陆轨迹非线性规划模型。着陆轨迹非线性规划模型中的变量包括状态变量、控制变量和终端时间，状态变量与探测器的质量、位置和速度有关，控制变量与探测器发动机的推力有关，终端时间为着陆过程结束的时刻。
步骤103，具体包括：
根据探测器着陆过程的时域范围确定着陆轨迹中的离散节点；时域范围为探测器着陆过程的始端时间与终端时间形成的闭区间；具体包括：
将时域范围由[0,tf]转化到[-1,1]，得到
Figure BDA0002628797740000094
其中，t为实际时间，τ为伪谱时间；
获取着陆轨迹中的离散节点个数，并根据离散节点个数根据如下公式确定着陆轨迹中的离散节点：
Figure BDA0002628797740000095
式中，τk为第k个离散节点，N+1为离散节点总数，可选的，N＝30。
采用直接配点法确定在离散节点处的状态变量和控制变量；
根据状态变量和控制变量，将着陆轨迹优化模型转换为着陆轨迹非线性规划模型，具体包括：
根据如下公式确定着陆轨迹非线性规划模型的目标函数：
minJ'＝-cTX
根据如下公式确定着陆轨迹非线性规划模型的约束条件：
Figure BDA0002628797740000102
Figure BDA0002628797740000103
||QiU||2≤mi TU
Figure BDA0002628797740000104
C1X＝X0,C2X＝Xf
其中，
X＝(R1 T R2 T R3 T V1 T V2 T V3 T ZT)T
U＝(U1 T U2 T U3 T ΣT)T
Figure BDA0002628797740000101
式中，J'为着陆轨迹非线性规划模型的目标函数，X为状态变量矩阵，U为控制变量矩阵，Rj为探测器位置向量，j＝1,2,3，R1为探测器水平横向位置向量，R2为探测器水平纵向位置向量，R3为探测器竖直方向位置向量，Vj为探测器速度向量，V1为探测器水平横向速度向量，V2为探测器水平纵向速度向量，V3为探测器竖直方向速度向量，Z为探测器质量对数向量，Uj为探测器第一控制向量，U1为探测器水平横向第一控制向量，U2为探测器水平纵向第一控制向量，U3为探测器竖直方向第一控制向量，Σ为探测器第二控制向量，rj(τN)为编号为N的离散节点处探测器的位置变量，vj(τN)为编号为N的离散节点处探测器的速度变量，z(τN)为编号为N的离散节点处探测器质量对数变量，uj(τN)为编号为N的离散节点处探测器第一控制变量，
Figure BDA0002628797740000105
u(t)为t时刻探测器第一控制变量，σ(τN)为编号为N的离散节点处探测器第二控制变量，
Figure BDA0002628797740000106
σ(t)为t时刻探测器第二控制变量，Γ(t)为松弛变量，||T(t)||≤Γ(t)，c为最后一个元素为1其余元素为0的7(N+1)维单位列向量，D为切比雪夫伪谱法对应的微分矩阵，I为7阶单位矩阵，A,B,G均为系数矩阵，
Figure BDA0002628797740000111
1(N+1)×1表示所有元素均为1的列向量，mi为除第3N+4+i个元素为1外，其余元素全为零的4(N+1)维列向量，ci为除第6N+7+i个元素为1外，其余元素全为零的7(N+1)维列向量，Tmin为最小飞行时间，Tmax为最大飞行时间，Qi为除元素(i,i)，(N+2+i,N+2+i)，(2N+3+i,2N+3+i)为1外，其余元素全为零的3(N+1)×4(N+1)阶矩阵，Pi为除元素(i,i)，(N+2+i,N+2+i)为1外，其余元素全为零的2(N+1)×7(N+1)阶矩阵，qi为除第2N+3+i个元素为1外，其余元素全为零的7(N+1)维列向量，C1为除元素(1,)1，(N+2,N+2)，(2N+3,2N+3)，(3N+4,3N+4)，(4N+5,4N+5)，(5N+6,5N+6)，(6N+7,6N+7)为1外，其余元素全为零的7×7(N+1)阶矩阵，C2为除元素(N+1,N+1)，(2(N+1),2(N+1))，(3(N+1),3(N+1))，(4(N+1),4(N+1))，(5(N+1),5(N+1))，(6(N+1),6(N+1))为1外，其余元素全为零的6×7(N+1)阶矩阵，X0为状态变量矩阵，X0＝(r0 T v0 T z0)T，z0为探测器的初始质量对数，Xf为终端状态变量矩阵，Xf＝(rf T vf T)T。
步骤104：获取参考轨迹，根据参考轨迹采用凸优化方法将着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型；参考轨迹包括参考状态值、参考控制值和参考终端时间。
步骤104，具体包括：
根据如下公式确定着陆轨迹二阶锥规划模型的目标函数：
minJ″＝-cTX
根据如下公式确定着陆轨迹二阶锥规划模型的约束条件：
Figure BDA0002628797740000112
Figure BDA0002628797740000113
Figure BDA0002628797740000114
||QiU||2≤mi TU
Figure BDA0002628797740000115
C1X＝X0,C2X＝Xf
式中，J″为着陆轨迹二阶锥规划模型的目标函数，{X(k),U(k),tf (k)}为参考轨迹，X(k)为第k次迭代后的状态变量矩阵，U(k)为第k次迭代后的控制变量矩阵，tf (k)为第k次迭代后的终端时间，{d0,d1,d2}为指数函数最小二乘逼近的系数。
对于第一次迭代时的参考轨迹{X(0),U(0),tf (0)}可以选择为：X(0)由初始状态和终端状态线性插值得到；U(0)由最小推力加速度和最大推力加速度线性插值得到；tf 0为最小飞行时间Tmin和最大飞行时间Tmax之间某一个值即可，如
Figure BDA0002628797740000121
步骤105：采用原始对偶内点法(primal-dual interior point method)对着陆轨迹二阶锥规划模型进行求解，得到期望轨迹。期望轨迹包括期望状态值、期望控制值和期望终端时间。
步骤106：判断着陆轨迹与期望轨迹的差值的切比雪夫范数是否小于预设阈值。若小于预设阈值，则执行步骤108；否则，执行步骤107。
步骤107：根据期望轨迹对参考轨迹进行更新，然后返回步骤104。
步骤108：根据期望轨迹采用插值法确定着陆过程的轨迹。
步骤109：输出着陆过程的轨迹以及期望终端时间。
图2为本发明实施例中在线轨迹规划系统结构图。如图2所示，一种在线轨迹规划系统，包括：
探测器的初始状态获取模块201，用于获取探测器的初始状态；初始状态包括着陆过程开始时探测器的质量、位置和速度；
着陆轨迹优化模型确定模块202，用于根据探测器的初始状态建立着陆轨迹优化模型；着陆轨迹优化模型以最小化燃料消耗量为目标，并根据行星探测器着陆过程的限制因素确定约束条件；
着陆轨迹非线性规划模型确定模块203，用于采用切比雪夫伪谱法将着陆轨迹优化模型转换为着陆轨迹非线性规划模型；着陆轨迹非线性规划模型中的变量包括状态变量、控制变量和终端时间，状态变量与探测器的质量、位置和速度有关，控制变量与探测器发动机的推力有关，终端时间为着陆过程结束的时刻；
着陆轨迹非线性规划模型确定模块203，具体包括：
离散节点确定单元，用于根据探测器着陆过程的时域范围确定着陆轨迹中的离散节点；时域范围为探测器着陆过程的始端时间与终端时间形成的闭区间；
变量确定单元，用于采用直接配点法确定在离散节点处的状态变量和控制变量；
着陆轨迹非线性规划模型确定单元，用于根据状态变量和控制变量，将着陆轨迹优化模型转换为着陆轨迹非线性规划模型。
着陆轨迹二阶锥规划模型确定模块204，用于获取参考轨迹，根据参考轨迹采用凸优化方法将着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型；参考轨迹包括参考状态值、参考控制值和参考终端时间；
轨迹计算模块205，用于采用原始对偶内点法对着陆轨迹二阶锥规划模型进行求解，得到期望轨迹；期望轨迹包括期望状态值、期望控制值和期望终端时间；
判断模块206，用于判断着陆轨迹与期望轨迹的差值的切比雪夫范数是否小于预设阈值；若小于预设阈值，则执行插值模块；否则，执行更新模块；
更新模块207，用于根据期望轨迹对参考轨迹进行更新，然后执行着陆轨迹二阶锥规划模型确定模块；
插值模块208，用于根据期望轨迹采用插值法确定着陆过程的轨迹，然后执行输出模块；
输出模块209，用于输出着陆过程的轨迹以及期望终端时间。
对于实施例公开的系统而言，由于其与实施例公开的方法相对应，所以描述的比较简单，相关之处参见方法部分说明即可。
本文中应用了具体个例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处。综上，本说明书内容不应理解为对本发明的限制。

Claims (7)

1.一种在线轨迹规划方法，其特征在于，包括：
获取探测器的初始状态；所述初始状态包括着陆过程开始时探测器的质量、位置和速度；
根据所述探测器的初始状态建立着陆轨迹优化模型；所述着陆轨迹优化模型以最小化燃料消耗量为目标，并根据行星探测器着陆过程的限制因素确定约束条件；
采用切比雪夫伪谱法将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型；所述着陆轨迹非线性规划模型中的变量包括状态变量、控制变量和终端时间，所述状态变量与探测器的质量、位置和速度有关，所述控制变量与探测器发动机的推力有关，所述终端时间为着陆过程结束的时刻；
获取参考轨迹，根据所述参考轨迹采用凸优化方法将所述着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型；所述参考轨迹包括参考状态值、参考控制值和参考终端时间；
所述根据所述参考轨迹采用凸优化方法将所述着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型，具体包括：
根据如下公式确定着陆轨迹二阶锥规划模型的目标函数：
min J”＝-cTX
根据如下公式确定着陆轨迹二阶锥规划模型的约束条件：
Figure FDA0003597996870000011
Figure FDA0003597996870000012
Figure FDA0003597996870000013
||QiU||2≤mi TU
Figure FDA0003597996870000014
C1X＝X0,C2X＝Xf
式中，J”为着陆轨迹二阶锥规划模型的目标函数，c为最后一个元素为1其余元素为0的7(N+1)维单位列向量，X为状态变量矩阵，D为切比雪夫伪谱法对应的微分矩阵，U为控制变量矩阵，I为7阶单位矩阵，{X(k),U(k),tf (k)}为参考轨迹，X(k)为第k次迭代后的状态变量矩阵，U(k)为第k次迭代后的控制变量矩阵，tf (k)为第k次迭代后的终端时间，A,B,G均为系数矩阵，
Figure FDA0003597996870000021
{d0,d1,d2}为指数函数最小二乘逼近的系数；mi为除第3N+4+i个元素为1外，其余元素全为零的4(N+1)维列向量；ci为除第6N+7+i个元素为1外，其余元素全为零的7(N+1)维列向量；Qi为除元素(i,i)，(N+2+i,N+2+i)，(2N+3+i,2N+3+i)为1外，其余元素全为零的3(N+1)×4(N+1)阶矩阵；Pi为除元素(i,i)，(N+2+i,N+2+i)为1外，其余元素全为零的2(N+1)×7(N+1)阶矩阵；qi为除第2N+3+i个元素为1外，其余元素全为零的7(N+1)维列向量；γ为预设滑翔倾角，C1为除元素(1,1)，(N+2,N+2)，(2N+3,2N+3)，(3N+4,3N+4)，(4N+5,4N+5)，(5N+6,5N+6)，(6N+7,6N+7)为1外，其余元素全为零的7×7(N+1)阶矩阵；C2为除元素(N+1,N+1)，(2(N+1),2(N+1))，(3(N+1),3(N+1))，(4(N+1),4(N+1))，(5(N+1),5(N+1))，(6(N+1),6(N+1))为1外，其余元素全为零的6×7(N+1)阶矩阵；X0为状态变量矩阵，X0＝(r0 T v0 T z0)T，z0为探测器的初始质量对数，r0为着陆过程开始时探测器的位置，v0为着陆过程开始时探测器的速度，Xf为终端状态变量矩阵，Xf＝(rf T vf T)T，rf为探测器着陆点的位置，vf为探测器安全着陆所需的速度；
采用原始对偶内点法对所述着陆轨迹二阶锥规划模型进行求解，得到期望轨迹；所述期望轨迹包括期望状态值、期望控制值和期望终端时间；
判断着陆轨迹与所述期望轨迹的差值的切比雪夫范数是否小于预设阈值；若小于所述预设阈值，则根据所述期望轨迹采用插值法确定着陆过程的轨迹，然后输出所述着陆过程的轨迹以及所述期望终端时间；否则，根据所述期望轨迹对所述参考轨迹进行更新，并返回步骤“获取参考轨迹，根据所述参考轨迹采用凸优化方法将所述着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型”。
2.根据权利要求1所述的在线轨迹规划方法，其特征在于，所述根据所述探测器的初始状态建立着陆轨迹优化模型，具体包括：
根据如下公式建立着陆轨迹优化模型的目标函数：
min J＝-m(tf)
根据如下公式建立着陆轨迹优化模型的约束条件：
Figure FDA0003597996870000031
Figure FDA0003597996870000032
fmin≤||T(t)||2≤fmax
Figure FDA0003597996870000033
m(0)＝m0,r(0)＝r0,v(0)＝v0
r(tf)＝rf,v(tf)＝vf
式中，J为着陆轨迹非线性规划模型的目标函数，m(tf)为探测器在tf时刻的质量，
Figure FDA0003597996870000034
为r(t)的导数，r(t)为探测器在t时刻的位置，
Figure FDA0003597996870000035
为v(t)的导数，v(t)为探测器在t时刻的速度，T(t)为制动发动机在t时刻的推力，
Figure FDA0003597996870000036
为m(t)的导数，m(t)为探测器在t时刻的质量，g为重力加速度，α为发动机推力与燃料消耗量之间的系数，fmin为推力下界值，fmax为推力上界值，e1为第1个分量为1其余分量为0的单位列向量，e2为第2个分量为1其余分量为0的单位列向量，e3为第3个分量为1其余分量为0的单位列向量，m0为着陆过程开始时探测器的质量，rf＝0，vf＝0，tf为终端时间。
3.根据权利要求2所述的在线轨迹规划方法，其特征在于，所述采用切比雪夫伪谱法将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型，具体包括：
根据探测器着陆过程的时域范围确定着陆轨迹中的离散节点；所述时域范围为探测器着陆过程的始端时间与终端时间形成的闭区间；
采用直接配点法确定在所述离散节点处的状态变量和控制变量；
根据所述状态变量和所述控制变量，将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型。
4.根据权利要求3所述的在线轨迹规划方法，其特征在于，所述根据探测器着陆过程的时域范围确定着陆轨迹中的离散节点，具体包括：
将所述时域范围由[0,tf]转化到[-1,1]，得到
Figure FDA0003597996870000037
其中，t为实际时间，τ为伪谱时间；
获取着陆轨迹中的离散节点个数，并根据所述离散节点个数根据如下公式确定着陆轨迹中的离散节点：
Figure FDA0003597996870000041
式中，τk为第k个离散节点，N+1为离散节点总数。
5.根据权利要求4所述的在线轨迹规划方法，其特征在于，所述根据所述状态变量和所述控制变量，将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型，具体包括：
根据如下公式确定所述着陆轨迹非线性规划模型的目标函数：
min J'＝-cTX
根据如下公式确定所述着陆轨迹非线性规划模型的约束条件：
Figure FDA0003597996870000042
Figure FDA0003597996870000043
||QiU||2≤mi TU
Figure FDA0003597996870000044
C1X＝X0,C2X＝Xf
其中，
X＝(R1 T R2 T R3 T V1 T V2 T V3 T ZT)T
U＝(U1 T U2 T U3 T ΣT)T
Figure FDA0003597996870000045
式中，J'为着陆轨迹非线性规划模型的目标函数，Rj为探测器位置向量，j＝1,2,3，R1为探测器水平横向位置向量，R2为探测器水平纵向位置向量，R3为探测器竖直方向位置向量，Vj为探测器速度向量，V1为探测器水平横向速度向量，V2为探测器水平纵向速度向量，V3为探测器竖直方向速度向量，Z为探测器质量对数向量，Uj为探测器第一控制向量，U1为探测器水平横向第一控制向量，U2为探测器水平纵向第一控制向量，U3为探测器竖直方向第一控制向量，Σ为探测器第二控制向量，rj(τN)为编号为N的离散节点处探测器的位置变量，vj(τN)为编号为N的离散节点处探测器的速度变量，z(τN)为编号为N的离散节点处探测器质量对数变量，uj(τN)为编号为N的离散节点处探测器第一控制变量，
Figure FDA0003597996870000051
u(t)为t时刻探测器第一控制变量，σ(τN)为编号为N的离散节点处探测器第二控制变量，
Figure FDA0003597996870000052
σ(t)为t时刻探测器第二控制变量，Γ(t)为松弛变量，||T(t)||≤Γ(t)，c为最后一个元素为1其余元素为0的7(N+1)维单位列向量，1(N+1)×1表示所有元素均为1的列向量；Tmin为最小飞行时间，Tmax为最大飞行时间。
6.一种在线轨迹规划系统，其特征在于，包括：
探测器的初始状态获取模块，用于获取探测器的初始状态；所述初始状态包括着陆过程开始时探测器的质量、位置和速度；
着陆轨迹优化模型确定模块，用于根据所述探测器的初始状态建立着陆轨迹优化模型；所述着陆轨迹优化模型以最小化燃料消耗量为目标，并根据行星探测器着陆过程的限制因素确定约束条件；
着陆轨迹非线性规划模型确定模块，用于采用切比雪夫伪谱法将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型；所述着陆轨迹非线性规划模型中的变量包括状态变量、控制变量和终端时间，所述状态变量与探测器的质量、位置和速度有关，所述控制变量与探测器发动机的推力有关，所述终端时间为着陆过程结束的时刻；
着陆轨迹二阶锥规划模型确定模块，用于获取参考轨迹，根据所述参考轨迹采用凸优化方法将所述着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型；所述参考轨迹包括参考状态值、参考控制值和参考终端时间；
所述根据所述参考轨迹采用凸优化方法将所述着陆轨迹非线性规划模型转换为着陆轨迹二阶锥规划模型，具体包括：
根据如下公式确定着陆轨迹二阶锥规划模型的目标函数：
min J”＝-cTX
根据如下公式确定着陆轨迹二阶锥规划模型的约束条件：
Figure FDA0003597996870000053
Figure FDA0003597996870000054
Figure FDA0003597996870000061
||QiU||2≤mi TU
Figure FDA0003597996870000062
C1X＝X0,C2X＝Xf
式中，J”为着陆轨迹二阶锥规划模型的目标函数，c为最后一个元素为1其余元素为0的7(N+1)维单位列向量，X为状态变量矩阵，D为切比雪夫伪谱法对应的微分矩阵，U为控制变量矩阵，I为7阶单位矩阵，{X(k),U(k),tf (k)}为参考轨迹，X(k)为第k次迭代后的状态变量矩阵，U(k)为第k次迭代后的控制变量矩阵，tf (k)为第k次迭代后的终端时间，A,B,G均为系数矩阵，
Figure FDA0003597996870000063
{d0,d1,d2}为指数函数最小二乘逼近的系数；mi为除第3N+4+i个元素为1外，其余元素全为零的4(N+1)维列向量；ci为除第6N+7+i个元素为1外，其余元素全为零的7(N+1)维列向量；Qi为除元素(i,i)，(N+2+i,N+2+i)，(2N+3+i,2N+3+i)为1外，其余元素全为零的3(N+1)×4(N+1)阶矩阵；Pi为除元素(i,i)，(N+2+i,N+2+i)为1外，其余元素全为零的2(N+1)×7(N+1)阶矩阵；qi为除第2N+3+i个元素为1外，其余元素全为零的7(N+1)维列向量；γ为预设滑翔倾角，C1为除元素(1,1)，(N+2,N+2)，(2N+3,2N+3)，(3N+4,3N+4)，(4N+5,4N+5)，(5N+6,5N+6)，(6N+7,6N+7)为1外，其余元素全为零的7×7(N+1)阶矩阵；C2为除元素(N+1,N+1)，(2(N+1),2(N+1))，(3(N+1),3(N+1))，(4(N+1),4(N+1))，(5(N+1),5(N+1))，(6(N+1),6(N+1))为1外，其余元素全为零的6×7(N+1)阶矩阵；X0为状态变量矩阵，X0＝(r0 T v0 T z0)T，z0为探测器的初始质量对数，r0为着陆过程开始时探测器的位置，v0为着陆过程开始时探测器的速度，Xf为终端状态变量矩阵，Xf＝(rf T vf T)T，rf为探测器着陆点的位置，vf为探测器安全着陆所需的速度；
轨迹计算模块，用于采用原始对偶内点法对所述着陆轨迹二阶锥规划模型进行求解，得到期望轨迹；所述期望轨迹包括期望状态值、期望控制值和期望终端时间；
判断模块，用于判断所述着陆轨迹与所述期望轨迹的差值的切比雪夫范数是否小于预设阈值；若小于所述预设阈值，则执行插值模块；否则，执行更新模块；
更新模块，用于根据所述期望轨迹对所述参考轨迹进行更新，然后执行所述着陆轨迹二阶锥规划模型确定模块；
插值模块，用于根据所述期望轨迹采用插值法确定着陆过程的轨迹，然后执行输出模块；
输出模块，用于输出所述着陆过程的轨迹以及所述期望终端时间。
7.根据权利要求6所述的在线轨迹规划系统，其特征在于，所述着陆轨迹非线性规划模型确定模块，具体包括：
离散节点确定单元，用于根据探测器着陆过程的时域范围确定着陆轨迹中的离散节点；所述时域范围为探测器着陆过程的始端时间与终端时间形成的闭区间；
变量确定单元，用于采用直接配点法确定在所述离散节点处的状态变量和控制变量；
着陆轨迹非线性规划模型确定单元，用于根据所述状态变量和所述控制变量，将所述着陆轨迹优化模型转换为着陆轨迹非线性规划模型。
CN202010805032.4A 2020-08-12 2020-08-12 一种行星探测器动力着陆段在线轨迹规划方法及系统 Active CN111931131B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010805032.4A CN111931131B (zh) 2020-08-12 2020-08-12 一种行星探测器动力着陆段在线轨迹规划方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010805032.4A CN111931131B (zh) 2020-08-12 2020-08-12 一种行星探测器动力着陆段在线轨迹规划方法及系统

Publications (2)

Publication Number Publication Date
CN111931131A CN111931131A (zh) 2020-11-13
CN111931131B true CN111931131B (zh) 2022-07-15

Family

ID=73311265

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010805032.4A Active CN111931131B (zh) 2020-08-12 2020-08-12 一种行星探测器动力着陆段在线轨迹规划方法及系统

Country Status (1)

Country Link
CN (1) CN111931131B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112693631B (zh) * 2020-11-27 2022-07-29 中国人民解放军国防科技大学 飞行器在线序列凸优化中的初始轨迹生成方法及系统
CN114872934B (zh) * 2022-05-24 2023-10-13 哈尔滨工业大学 基于无损凸化和Chebyshev正交配点法的月面飞行器轨迹优化方法

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109211246A (zh) * 2018-11-07 2019-01-15 北京理工大学 不确定环境下行星着陆轨迹规划方法
CN110466804A (zh) * 2019-08-30 2019-11-19 北京理工大学 火箭动力下降着陆过程快速轨迹优化方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106647282B (zh) * 2017-01-19 2020-01-03 北京工业大学 一种考虑末端运动误差的六自由度机器人轨迹规划方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109211246A (zh) * 2018-11-07 2019-01-15 北京理工大学 不确定环境下行星着陆轨迹规划方法
CN110466804A (zh) * 2019-08-30 2019-11-19 北京理工大学 火箭动力下降着陆过程快速轨迹优化方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于凸优化和切比雪夫伪谱法的再入轨迹优化;陈嘉澍;《电子技术与软件工程》;20180831(第138期);第71-74页 *

Also Published As

Publication number Publication date
CN111931131A (zh) 2020-11-13

Similar Documents

Publication Publication Date Title
CN107102644B (zh) 基于深度强化学习的水下机器人轨迹控制方法及控制系统
CN111931131B (zh) 一种行星探测器动力着陆段在线轨迹规划方法及系统
CN111798491B (zh) 一种基于Elman神经网络的机动目标跟踪方法
CN107729706B (zh) 一种非线性机械系统的动力学模型构建方法
CN105678015B (zh) 一种高超声速三维机翼的非概率可靠性气动结构耦合优化设计方法
CN111027732A (zh) 一种多风电场出力场景的生成方法及系统
CN108168577A (zh) 基于bp神经网络的mems陀螺随机误差补偿方法
CN106441300A (zh) 一种具有自适应的协同导航滤波方法
CN111191368B (zh) 一种连续小推力行星际转移轨道优化方法和装置
CN110706265A (zh) 一种改进srckf强跟踪滤波的机动目标跟踪方法
CN106383443B (zh) 抗干扰控制方法及系统
CN104614985A (zh) 一种基于非线性规划的高阶系统最优降阶方法
CN110532621A (zh) 一种飞行器气动参数在线辨识方法
CN105675017A (zh) 一种应用于光电平台的光纤陀螺随机漂移补偿方法
CN118113064A (zh) 基于在线修正智能辨识的变体飞行器智能复合控制方法
Jameel et al. Approximate solution of high order fuzzy initial value problems
CN115338610B (zh) 双轴孔装配方法、装置、电子设备和存储介质
Iqbal et al. Comparison of various wiener model identification approach in modelling nonlinear process
Tararykin et al. Synthesizing parametrically robust control systems with state controllers and observers based on gramian method
CN104200002A (zh) 一种粘性阻尼振动信号中的模态参数提取方法
CN114548400A (zh) 一种快速灵活全纯嵌入式神经网络广域寻优训练方法
CN104038132B (zh) 具有时变测量延迟输出和噪声的伺服电机的状态观测方法
CN114137971A (zh) 一种转向系统延迟的离线辨识方法
CN117826616B (zh) 一种基于序列凸优化的飞行器快速轨迹规划方法及装置
Sun et al. A modified marginalized Kalman filter for maneuvering target tracking

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
GR01 Patent grant