CN114357853A - 基于遗传算法修正的多圈J2-Lambert转移轨道求解方法 - Google Patents

基于遗传算法修正的多圈J2-Lambert转移轨道求解方法 Download PDF

Info

Publication number
CN114357853A
CN114357853A CN202111395792.3A CN202111395792A CN114357853A CN 114357853 A CN114357853 A CN 114357853A CN 202111395792 A CN202111395792 A CN 202111395792A CN 114357853 A CN114357853 A CN 114357853A
Authority
CN
China
Prior art keywords
initial
lambert
vector
position vector
turn
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
CN202111395792.3A
Other languages
English (en)
Other versions
CN114357853B (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.)
Shanghai Aerospace Control Technology Institute
Original Assignee
Shanghai Aerospace Control Technology Institute
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 Shanghai Aerospace Control Technology Institute filed Critical Shanghai Aerospace Control Technology Institute
Priority to CN202111395792.3A priority Critical patent/CN114357853B/zh
Publication of CN114357853A publication Critical patent/CN114357853A/zh
Application granted granted Critical
Publication of CN114357853B publication Critical patent/CN114357853B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Control Of Position, Course, Altitude, Or Attitude Of Moving Bodies (AREA)

Abstract

本发明公开了一种基于遗传算法修正的多圈J2‑Lambert转移轨道求解方法，包括：获取初始参数，初始位置矢量，目标位置矢量，和航天器从初始位置矢量到达目标位置矢量时的转移时间；根据初始参数，求解无J2摄动时，航天器从初始位置矢量经过转移时间，到达目标位置矢量的多圈Lambert转移轨道；建立J2摄动下的动力学方程，根据动力学方程结合多圈Lambert转移轨道得到因J2摄动对转移轨道的最终位置产生的偏差；利用遗传算法修正多圈Lambert转移轨道，直至使得偏差达到最小值，输出最终的多圈Lambert转移轨道。本发明经过遗传算法修正的多圈J2‑Lambert转移轨道不仅适用于长周期变轨任务，且满足落点精度要求。

Description

基于遗传算法修正的多圈J2-Lambert转移轨道求解方法
技术领域
本发明涉及航天器轨道设计的技术领域，特别涉及一种基于遗传算法修正的多圈J2-Lambert转移轨道求解方法。
背景技术
空间交会对接是航天任务中一项非常重要的步骤，在载人登月，空间站建设，星际航行等任务中均有应用，空间交会对接中涉及到的一个关键技术就是Lambert变轨技术。具体来说，Lambert问题描述了“交会”过程，即航天器如何在预定的时间内，从空间中一点到达另外一个目标点。在确定航天器的引力场之后，只需知道航天器的初始速度即可推算出预定时间飞行后的到达点。该问题在航天动力学领域应用广泛，除了空间交会对接以及星际飞行之外，还有卫星拦截，弹道导弹拦截等等。
由于Lambert问题不存在解析解，因而Lambert问题的求解大多采用数值优化算法。当Lambert变轨时间较短时，Lambert转移轨道有且只有一条。但是当两点之间的变轨飞行时间很长时，满足约束条件的Lambert转移轨道的解并不唯一。若此时仍然采用单圈轨道设计，则转移轨道将呈现较大的偏心率，大大增加了Lambert交会任务所需的燃料消耗。为了避免上述情况，可采用多圈lambert转移轨道，将过长的转移时间消耗在多余的Lambert转移轨道圈数上，不以燃料消耗为代价。然而目前已有文献的多圈Lambert转移轨道设计没有考虑地球非球形摄动的影响，会增加航天器在预定飞行时间之后的落点偏差。
发明内容
本发明的目的在于提供一种基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，以弥补以往只设计单圈J2-Lambert转移轨道修正以及多圈无摄动Lambert转移轨道设计工作的不足，实现消除地球扁率对长周期变轨的影响，有利于空间交会对接任务的顺利实现的目的。
为了实现以上目的，本发明通过以下技术方案实现：
一种基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，包括：
步骤S1、获取初始参数，所述初始参数包括：初始位置矢量Ri，目标位置矢量Rf，和航天器从所述初始位置矢量Ri到达所述目标位置矢量Rf时的转移时间T。
步骤S2、根据所述初始参数，求解无J2摄动时，航天器从所述初始位置矢量Ri经过所述转移时间T，到达所述目标位置矢量Rf的多圈Lambert转移轨道。
步骤S3、建立J2摄动下的动力学方程，根据所述动力学方程结合所述多圈Lambert转移轨道得到因J2摄动对转移轨道的最终位置产生的偏差。
步骤S4、利用遗传算法修正所述多圈Lambert转移轨道，直至使得所述偏差达到最小值，输出最终的所述多圈Lambert转移轨道。
可选地，多圈Lambert转移轨道包括：初始速度矢量Vi与末端速度矢量Vf。
可选地，所述步骤S2包括：步骤S2.1、采用如下公式计算航天器的初
始地心距离ri与目标地心距离rf：
ri＝||Ri||，rf＝||Rf||
采用如下公式计算初始位置与目标位置的中心角Δθ：
Figure BDA0003370291260000021
步骤S2.2、根据所述初始地心距离ri，所述目标地心距离rf和所述中心角Δθ计算所述多圈Lambert转移轨道的六要素参数。
步骤S2.3、根据六要素参数与多圈Lambert转移轨道位置和速度状态量之间的转换关系，转化得到初始速度矢量Vi与末端速度矢量Vf。
可选地，所述步骤S2.2包括：所述多圈Lambert转移轨道的六要素参数包括：半长轴a，偏心率e，轨道倾角i，升交点赤经Ω，初始位置的纬度幅角ui，目标位置的纬度幅角uf，初始位置的真近点角θi和目标位置的真近点角θf。
其中，所述半长轴a采用如下方式进行计算，将a0＝(ri+rf)/2作为初值代入牛顿迭代法中计算所述多圈Lambert转移轨道的半长轴a。
所述偏心率e采用如下公式计算:
Figure BDA0003370291260000031
Figure BDA0003370291260000032
Figure BDA0003370291260000033
α＝π+sign(Δtm-T)(α*-π)
β＝sign(sinΔθ)β*
Figure BDA0003370291260000034
Figure BDA0003370291260000035
Figure BDA0003370291260000036
Figure BDA0003370291260000037
式中，μ＝GM＝3.986*1014[m3/s2]，其中G是万有引力常数G＝6.670*10-11[m3/(kg.s2)]；M是地球质量M＝5.976*1024[kg]。
采用如下公式计算所述真近点角θi和目标位置的真近点角θf：
Figure BDA0003370291260000038
Figure BDA0003370291260000041
Figure BDA0003370291260000042
Figure BDA0003370291260000043
采用如下公式计算所述轨道倾角i：
Figure BDA0003370291260000044
Figure BDA0003370291260000045
q＝sign(π-Δθ)
式中，qxi表示轨道法线矢量在地心惯性坐标系中的x轴分量，qyi表示轨道法线矢量在地心惯性坐标系中的y轴分量；qzi表示轨道法线矢量在地心惯性坐标系中的z轴分量，q表示轨道的法线矢量，q表示运动方向，若Δθ<π，q＝+1，否则q＝-1。
采用如下公式计算所述升交点赤经Ω：
Figure BDA0003370291260000046
采用如下公式计算所述初始位置的纬度幅角ui和目标位置的纬度幅角uf
Figure BDA0003370291260000047
Figure BDA0003370291260000048
Figure BDA0003370291260000051
Figure BDA0003370291260000052
{Rf}n＝[{Rf}nx {Rf}ny {Rf}nz]T＝LniRf
转移轨道的近地点辐角ω：ω＝ui-θi＝uf-θf。
式中，{Ri}n表示初始位置矢量在地心节点坐标系中的表达式，{Ri}nx表示初始位置矢量在地心节点坐标系中的x轴方向分量；{Ri}ny表示初始位置矢量在地心节点坐标系中的y轴方向分量；{Ri}nz表示初始位置矢量在地心节点坐标系中的z轴方向分量；{Rf}n表示目标位置矢量在地心节点坐标系中的表达式，{Rf}nx表示目标位置矢量在地心节点坐标系中的x轴方向分量；{Rf}ny表示目标位置矢量在地心节点坐标系中的y轴方向分量；{Rf}nz表示目标位置矢量在地心节点坐标系中的z轴方向分量。
可选地，所述步骤S2.3包括：初始速度矢量Vi与末端速度矢量Vf采用如下公式计算：
[Ri Vi]＝Ψ(a,e,i,ui,Ω,θi)
[Rf Vf]＝Ψ(a,e,i,uf,Ω,θf)
其中，Vi为初始位置的速度矢量，Vf为目标位置的速度矢量，Ψ是轨道六要素与轨道位置和速度状态量之间的转换关系。
可选地，所述步骤S3包括：从初始状态[Ri，Vi]开始，数值积分时间T，得到实际末端位置矢量Rf'(Vi)，计算因J2摄动对多圈Lambert转移轨道的最终位置产生的偏差：所述J2摄动下的动力学方程
Figure BDA0003370291260000053
为：
Figure BDA0003370291260000061
其中，R为位置矢量，R＝[x，y，z]T，x表示位置矢量在地心惯性坐标系中的x轴方向分量；y表示位置矢量在地心惯性坐标系中的y轴方向分量；z表示位置矢量在地心惯性坐标系中的z轴方向分量；R随时间变化，当飞行时间t＝0时，R＝Ri，表示初始位置矢量，当飞行时间等于转移时间t＝T时，R＝Rf’(Vi)表示最终位置矢量；
Figure BDA0003370291260000062
是位置矢量关于飞行时间的二阶导数；r为该位置到地心的距离，r＝||R||，μ＝GM＝3.986*1014[m3/s2]，其中G是万有引力常数G＝6.670*10-11[m3/(kg.s2)]；M是地球质量M＝5.976*1024[kg]，J2是地球扁率摄动系数，Re为地球半径。
因J2摄动对转移轨道最终位置产生的偏差ΔRf(Vi)为
ΔRf(Vi)＝||Rf'(Vi)-Rf||
其中，Rf'(Vi)为最终位置矢量。
可选地，所述步骤S4包括：步骤S4.1、遗传算法初始化：随机生成初始种群，所述初始种群为随机生成一组初始速度Vi可行解的一个集合
[Vi1,Vi2,Vi3…Vi(n-1),Vin]，其中n是所述初始种群可行解的个数。
步骤S4.2、设置目标函数为ΔRf(Vi)＝||Rf'(Vi)-Rf||，取适应度函数
Figure BDA0003370291260000063
步骤S4.3、按照所述初始种群中各个个体适应度函数值的大小，从初始种群中选出适应度函数值前若干的个体构成交配池，并将这些个体复制到下一代种群。
步骤S4.4、由交叉和变异这两个遗传算子对交配池中的若干个体进行操作，形成新一代的种群中的m个个体，其中m是新一代种群中新个体的数目；
步骤S4.5、将所述步骤S4.3复制的交配池中(n-m)个体和步骤S4.4产生的m个个体组合成新一代种群。
步骤S4.6、反复执行步骤S4.2～步骤S4.5，直至最优个体的适应度不再上升时，或者迭代次数达到预设的代数时，遗传算法终止，输出最优初始速度。
另一方面，本发明提供一种电子设备，包括处理器和存储器，所述存储器上存储有计算机程序，所述计算机程序被所述处理器执行时，实现如上文所述的方法。
再一方面，本发明还提供一种可读存储介质，所述可读存储介质内存储有计算机程序，所述计算机程序被处理器执行时，实现如上文所述的方法。
本发明至少具有以下优点之一：
本发明所提供的一种考虑J2摄动下多圈Lambert转移轨道的求解与修正方法，适用于长周期变轨需求的交会对接任务。
(1)本发明涉及的多圈J2-Lambert转移轨道为Lambert问题提供了新的轨道，可解决当变轨时间变长，单次Lambert变轨所导致的大偏心率转移轨道燃料消耗大的问题，为寻找Lambert燃料最优转移轨道提供选择。
(2)本发明利用遗传算法修正初始位置的速度生成多圈J2-Lambert转移轨道，转移轨道的动力学模型更加精确，可有效减少交会任务中因J2摄动引起的落点误差。
(3)遗传算法作为一种全局优化概率算法，可快速搜索解空间中的可行解，且利用遗传算法内在并行性的特点，可加快求解速度，节省计算时间。
附图说明
图1是本发明一实施例提供的一种基于遗传算法修正的多圈J2-Lambert转移轨道求解方法的流程图；
图2是本发明一实施例提供的Lambert转移轨道所存在的问题的示意图；
图3是本发明一实施例提供的单圈J2-Lambert转移轨道的示意图；
图4是本发明一实施例提供的一种基于遗传算法修正的多圈J2-Lambert转移轨道求解方法中的遗传算法求解步骤流程图；
图5是本发明一实施例提供的多圈J2-Lambertt转移轨道的示意图。
具体实施方式
以下结合附图和具体实施方式对本发明提出的一种基于遗传算法修正的多圈J2-Lambert转移轨道求解方法作进一步详细说明。根据下面说明，本发明的优点和特征将更清楚。需要说明的是，附图采用非常简化的形式且均使用非精准的比例，仅用以方便、明晰地辅助说明本发明实施方式的目的。为了使本发明的目的、特征和优点能够更加明显易懂，请参阅附图。须知，本说明书所附图式所绘示的结构、比例、大小等，均仅用以配合说明书所揭示的内容，以供熟悉此技术的人士了解与阅读，并非用以限定本发明实施的限定条件，故不具技术上的实质意义，任何结构的修饰、比例关系的改变或大小的调整，在不影响本发明所能产生的功效及所能达成的目的下，均应仍落在本发明所揭示的技术内容能涵盖的范围内。
如图2所示，J2-Lambert变轨问题可以描述为：考虑J2摄动，在给定航天器初始时刻的位置矢量Ri，变轨之后的目标位置矢量Rf，以及航天器从初始位置矢量Ri到目标位置矢量Rf的飞行时间间隔(变轨时间或移动时间)T，确定航天器在初始位置的速度矢量Vi与目标位置的速度矢量Vf，从而确定在两个位置处点火的速度增量。
如图1所示，本实施例提供的一种基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，包括：
步骤S1、获取初始参数，所述初始参数包括：初始位置矢量Ri，目标位置矢量Rf，和航天器从所述初始位置矢量Ri到达所述目标位置矢量Rf时的转移时间T。
步骤S2、根据所述初始参数，求解无J2摄动时，航天器从所述初始位置矢量Ri经过所述转移时间T，到达所述目标位置矢量Rf的多圈Lambert转移轨道；即求解航天器的初始速度矢量Vi与末端速度矢量Vf。
具体的，所述步骤S2包括：
步骤S2.1、采用如下公式计算航天器的初始地心距离ri与目标地心距离rf：
ri＝||Ri||，rf＝||Rf|| (1)
采用如下公式计算初始位置与目标位置的中心角Δθ：
Figure BDA0003370291260000091
步骤S2.2、根据所述初始地心距离ri，所述目标地心距离rf和所述中心角Δθ计算所述多圈Lambert转移轨道的六要素参数。
在本实施例中，所述步骤S2.2包括：所述多圈Lambert转移轨道的六要素参数包括：半长轴a，偏心率e，轨道倾角i，升交点赤经Ω，初始位置的纬度幅角ui，目标位置的纬度幅角uf，初始位置的真近点角θi和目标位置的真近点角θf。
其中，所述半长轴a采用如下方式进行计算，将a0＝(ri+rf)/2作为初值代入牛顿迭代法中计算所述多圈Lambert转移轨道的半长轴a；
牛顿迭代法的公式为：
Figure BDA0003370291260000092
式中，ak表示半长轴的第k次迭代值，ak+1表示半长轴的第k+1次迭代值；k＝0,1,2…F(ak)表示Lambert方程；F'(ak)表示Lambert方程关于半长轴的一阶导数。
F(a)计算公式如下：
Figure BDA0003370291260000093
式中，N表示转移圈数，T表示转移时间；μ＝GM＝3.986*1014[m3/s2]，其中G是万有引力常数G＝6.670*10-11[m3/(kg.s2)]；M是地球质量M＝5.976*1024[kg]；a表示半长轴；α和均为中间变量，无具体定义。
所述偏心率e采用如下公式计算:
Figure BDA0003370291260000101
Figure BDA0003370291260000102
Figure BDA0003370291260000103
α＝π+sign(Δtm-T)(α*-π) (8)
β＝sign(sinΔθ)β* (9)
Figure BDA0003370291260000104
Figure BDA0003370291260000105
Figure BDA0003370291260000106
Figure BDA0003370291260000107
sign是符号函数，
Figure BDA0003370291260000108
式中，μ为地球引力常数；p表示半通径；s表示初始位置与目标位置之间的直线距离；α*和β*均表示中间变量无具体的定义；Δtm表示沿最小能量椭圆弧的飞行时间；δm表示中间变量无具体的定义。
采用如下公式计算所述真近点角θi和目标位置的真近点角θf：
Figure BDA0003370291260000111
Figure BDA0003370291260000112
Figure BDA0003370291260000113
Figure BDA0003370291260000114
采用如下公式计算所述轨道倾角i：
Figure BDA0003370291260000115
Figure BDA0003370291260000116
q＝sign(π-Δθ)
式中，qxi表示轨道法线矢量在地心惯性坐标系中的x轴分量，qyi表示轨道法线矢量在地心惯性坐标系中的y轴分量；qzi表示轨道法线矢量在地心惯性坐标系中的z轴分量，q表示轨道的法线矢量，q表示运动方向，若Δθ<π，q＝+1，否则q＝-1。
采用如下公式计算所述升交点赤经Ω：
Figure BDA0003370291260000117
采用如下公式计算所述初始位置的纬度幅角ui和目标位置的纬度幅角uf，以及转移轨道的近地点辐角ω：：
ω＝ui-θi (22)
或者，ω＝uf-θf (23)
Figure BDA0003370291260000121
{Ri}n＝[{Ri}nx {Ri}ny {Ri}nz]T＝LniRi (25)
Figure BDA0003370291260000122
Figure BDA0003370291260000123
{Rf}n＝[{Rf}nx {Rf}ny {Rf}nz]T＝LniRf (28)
式中，{Ri}n表示初始位置矢量在地心节点坐标系中的表达式，{Ri}nx表示初始位置矢量在地心节点坐标系中的x轴方向分量；{Ri}ny表示初始位置矢量在地心节点坐标系中的y轴方向分量；{Ri}nz表示初始位置矢量在地心节点坐标系中的z轴方向分量；{Rf}n表示目标位置矢量在地心节点坐标系中的表达式，{Rf}nx表示目标位置矢量在地心节点坐标系中的x轴方向分量；{Rf}ny表示目标位置矢量在地心节点坐标系中的y轴方向分量；{Rf}nz表示目标位置矢量在地心节点坐标系中的z轴方向分量。矩阵Lni表示地心惯性坐标系到地心节点坐标系的转换矩阵，地心惯性坐标系：原点在地心，x轴在赤道平面内，指向春分点，z轴垂直于赤道平面，指向北极，y轴按右手法则确定。地心节点坐标系：原点在地心，x轴指向升交点，Z轴指向动量矩方向，y轴按右手法则确定。
步骤S2.3、根据六要素参数与多圈Lambert转移轨道位置和速度状态量之间的转换关系，转化得到初始速度矢量Vi与末端速度矢量Vf。
在本实施例中，所述步骤S2.3包括：初始速度矢量Vi与末端速度矢量Vf采用如下公式计算：
[Ri Vi]＝Ψ(a, e, i, ui, Ω, θi) (29)
[Rf Vf]＝Ψ(a, e, i, uf, Ω, θf) (30)
其中，Vi为初始位置的速度矢量，Vf为目标位置的速度矢量，Ψ是轨道六要素与轨道位置和速度状态量之间的转换关系。
步骤S3、建立J2摄动下的动力学方程，根据所述动力学方程结合所述多圈Lambert转移轨道得到因J2摄动对转移轨道的最终位置产生的偏差。
具体的，所述步骤S3包括：从初始状态[Ri，Vi]开始，数值积分时间T，得到实际末端位置矢量Rf'(Vi)，计算因J2摄动对多圈Lambert转移轨道的最终位置产生的偏差：
所述J2摄动下的动力学方程
Figure BDA0003370291260000131
为：
Figure BDA0003370291260000132
其中，R为位置矢量，R＝[x，y，z]T，x表示位置矢量在地心惯性坐标系中的x轴方向分量；y表示位置矢量在地心惯性坐标系中的y轴方向分量；z表示位置矢量在地心惯性坐标系中的z轴方向分量；R随时间变化，当飞行时间t＝0时，R＝Ri，表示初始位置矢量，当飞行时间等于转移时间t＝T时，R＝Rf’(Vi)表示最终位置矢量；
Figure BDA0003370291260000133
是位置矢量关于飞行时间的二阶导数；r为该位置到地心的距离，r＝||R||，μ＝GM＝3.986*1014[m3/s2]，其中G是万有引力常数G＝6.670*10-11[m3/(kg.s2)]；M是地球质量M＝5.976*1024[kg]；J2是地球扁率摄动系数，Re为地球半径。
因J2摄动对转移轨道最终位置产生的偏差ΔRf(Vi)为
ΔRf(Vi)＝||Rf'(Vi)-Rf|| (32)
其中，Rf'(Vi)该变量的定义为最终位置矢量，计算过程是根据初始位置矢量Ri，飞行时间T，以及J2摄动下的动力学方程，通过MATLAB软件运用四阶龙格库塔法积分得到。
步骤S4、利用遗传算法修正所述多圈Lambert转移轨道，直至使得所述偏差达到最小值，输出最终的所述多圈Lambert转移轨道。
如图4所示，所述步骤S4包括：步骤S4.1、遗传算法初始化：在一定的编码方案下，随机生成初始种群，所述初始种群为随机生成一组初始速度Vi可行解的一个集合[Vi1,Vi2,Vi3…Vi(n-1),Vin]，其中n是所述初始种群可行解的个数。
步骤S4.2、设置目标函数为ΔRf(Vi)＝||Rf'(Vi)-Rf||，由于目标函数为最小化问题，所以取适应度函数G(Vi)：
Figure BDA0003370291260000141
步骤S4.3、按照所述初始种群中各个个体适应度函数值的大小，从初始种群中选出适应度函数值前若干的个体构成交配池，并将这些个体复制到下一代种群。
步骤S4.4、由交叉和变异这两个遗传算子对交配池中的若干个体进行操作，形成新一代的种群中的m个个体，其中m是新一代种群中新个体的数目；
步骤S4.5、将所述步骤S4.3复制的交配池中(n-m)个体和步骤S4.4产生的m个个体组合成新一代种群。
步骤S4.6、反复执行步骤S4.2～步骤S4.5，直至最优个体的适应度不再上升时，或者迭代次数达到预设的代数时，遗传算法终止，输出最优初始速度。
为了更清楚的了解上述实施例，下面举一个具体的示例进行说明。步骤S1中，初始位置矢量Ri＝[934220.899；1646136.726；-3702.993]km，末端位置矢量Rf＝[-493688.816；445535.111；-5397.137]km，变轨时间T＝14.9935天。
特别地，单圈Lambert转移轨道是本发明的一个特例，针对本算例的单圈J2-Lambert转移轨道如图3所示，从初始点转移到目标点只需要大概半圈轨道，轨道偏心率较大(e＝0.917)。下面展示同输入参数下的多圈J2-Lambert转移轨道。
步骤S2.1、采用上述公式(1)计算航天器的初始地心距离ri与目标地心距离rf，以及采用上述公式(2)计算初始位置与末端位置的中心角Δθ，结果如下：
ri＝||Ri||＝1.8928·106km，rf＝||Rf||＝6.6503·105km
Δθ＝1.3528rad。
步骤S2.2中，采用上述公式(3)-(4)计算所述半长轴a，即将a0＝(ri+rf)/2＝1278900km作为初值代入牛顿迭代法计算半长轴a＝1235100km。
采用上述公式(5)～(28)进行计算得到e＝0.822，θi＝2.8546rad和θf＝4.2074rad，i＝0.0081rad，Ω＝3.9528rad，ui＝3.38497rad，uf＝4.737786rad和ω＝0.5304rad。
采用上述公式(29)和公式(30)计算结果为初始位置的速度矢量Vi＝[-1.2304；5.4568；-0.0377]km/s，目标位置的速度矢量(末端速度矢量)Vf＝[2.3213；-16.5235；0.1060]km/s。
所述步骤S3中采用上述公式(31)～(32)计算因J2摄动对转移轨道最终位置产生的偏差为ΔRf(Vi)＝||Rf'(Vi)-Rf||＝3309km。
步骤S4中，经遗传算法优化后，最终修正的初始速度Vi＝[-1.2178；5.4584；-0.0377]km/s，ΔRf(Vi)＝||Rf'(Vi)-Rf||＝0.47449m，多圈J2-Lambert转移轨道图5所示，从初始点转移到目标点需要一圈半的转移轨道，轨道偏心率变小(e＝0.822)。
另一方面，本发明还提供一种电子设备，包括处理器和存储器，所述存储器上存储有计算机程序，所述计算机程序被所述处理器执行时，实现如上文所述的方法。
再一方面，本发明还提供一种可读存储介质，其特征在于，所述可读存储介质内存储有计算机程序，所述计算机程序被处理器执行时，实现如上文所述的方法。
综上所述，本实施例公开了一种基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，该方法是在考虑J2摄动的动力学方程下，利用遗传算法修正航天器的初始速度，求解出航天器从初始位置矢量Ri，经过变轨时间T，到达目标末端位置矢量Rf的转移轨道。首先根据已知的输入参数，求解出无摄动情况下从Ri位置到Rf位置的多圈Lambert转移轨道初始速度矢量Vi；将初始状态[Ri，Vi]代入J2-Lambert动力学方程，数值积分时间T，得到实际末端位置矢量Rf’，并计算其与目标末端位置矢量Rf的偏差；以Vi为修正量，利用遗传算法，通过交叉变异等优化操作求解出满足末端位置偏差精度的初始速度矢量。经过遗传算法修正的多圈J2-Lambert转移轨道不仅适用于长周期变轨任务，且满足落点精度要求。
需要说明的是，在本文中，诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。
应当注意的是，在本文的实施方式中所揭露的装置和方法，也可以通过其他的方式实现。以上所描述的装置实施方式仅仅是示意性的，例如，附图中的流程图和框图显示了根据本文的多个实施方式的装置、方法和计算机程序产品的可能实现的体系架构、功能和操作。在这点上，流程图或框图中的每个方框可以代表一个模块、程序或代码的一部分，所述模块、程序段或代码的一部分包含一个或多个用于实现规定的逻辑功能的可执行指令，所述模块、程序段或代码的一部分包含一个或多个用于实现规定的逻辑功能的可执行指令。也应当注意，在有些作为替换的实现方式中，方框中所标注的功能也可以以不同于附图中所标注的顺序发生。例如，两个连续的方框实际上可以基本并行地执行，它们有时也可以按相反的顺序执行，这依所涉及的功能而定。也要注意的是，框图和/或流程图中的每个方框、以及框图和/或流程图中的方框的组合，可以用于执行规定的功能或动作的专用的基于硬件的系统来实现，或者可以用专用硬件与计算机指令的组合来实现。
另外，在本文各个实施方式中的各功能模块可以集成在一起形成一个独立的部分，也可以是各个模块单独存在，也可以两个或两个以上模块集成形成一个独立的部分。
尽管本发明的内容已经通过上述优选实施例作了详细介绍，但应当认识到上述的描述不应被认为是对本发明的限制。在本领域技术人员阅读了上述内容后，对于本发明的多种修改和替代都将是显而易见的。因此，本发明的保护范围应由所附的权利要求来限定。

Claims (9)

1.一种基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，其特征在于，包括：
步骤S1、获取初始参数，所述初始参数包括：初始位置矢量Ri，目标位置矢量Rf，和航天器从所述初始位置矢量Ri到达所述目标位置矢量Rf时的转移时间T；
步骤S2、根据所述初始参数，求解无J2摄动时，航天器从所述初始位置矢量Ri经过所述转移时间T，到达所述目标位置矢量Rf的多圈Lambert转移轨道；
步骤S3、建立J2摄动下的动力学方程，根据所述动力学方程结合所述多圈Lambert转移轨道得到因J2摄动对转移轨道的最终位置产生的偏差；
步骤S4、利用遗传算法修正所述多圈Lambert转移轨道，直至使得所述偏差达到最小值，输出最终的所述多圈Lambert转移轨道。
2.如权利要求1所述的基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，其特征在于，多圈Lambert转移轨道包括：初始速度矢量Vi与末端速度矢量Vf。
3.如权利要求2所述的基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，其特征在于，所述步骤S2包括：
步骤S2.1、采用如下公式计算航天器的初始地心距离ri与目标地心距离rf：
ri＝||Ri||，rf＝||Rf||
采用如下公式计算初始位置与目标位置的中心角Δθ：
Figure FDA0003370291250000021
步骤S2.2、根据所述初始地心距离ri，所述目标地心距离rf和所述中心角Δθ计算所述多圈Lambert转移轨道的六要素参数；
步骤S2.3、根据六要素参数与多圈Lambert转移轨道位置和速度状态量之间的转换关系，转化得到初始速度矢量Vi与末端速度矢量Vf。
4.如权利要求3所述的基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，其特征在于，所述步骤S2.2包括：所述多圈Lambert转移轨道的六要素参数包括：半长轴a，偏心率e，轨道倾角i，升交点赤经Ω，初始位置的纬度幅角ui，目标位置的纬度幅角uf，初始位置的真近点角θi和目标位置的真近点角θf：
其中，所述半长轴a采用如下方式进行计算，将a0＝(ri+rf)/2作为初值代入牛顿迭代法中计算所述多圈Lambert转移轨道的半长轴a；
所述偏心率e采用如下公式计算:
Figure FDA0003370291250000022
Figure FDA0003370291250000023
Figure FDA0003370291250000024
α＝π+sign(Δtm-T)(α*-π)
β＝sign(sinΔθ)β*
Figure FDA0003370291250000025
Figure FDA0003370291250000026
Figure FDA0003370291250000031
Figure FDA0003370291250000032
式中，μ＝GM＝3.986*1014[m3/s2]，其中G是万有引力常数G＝6.670*10-11[m3/(kg.s2)]；M是地球质量M＝5.976*1024[kg]；
采用如下公式计算所述真近点角θi和目标位置的真近点角θf：
Figure FDA0003370291250000033
Figure FDA0003370291250000034
Figure FDA0003370291250000035
Figure FDA0003370291250000036
采用如下公式计算所述轨道倾角i：
Figure FDA0003370291250000037
Figure FDA0003370291250000038
q＝sign(π-Δθ)
式中，qxi表示轨道法线矢量在地心惯性坐标系中的x轴分量，qyi表示轨道法线矢量在地心惯性坐标系中的y轴分量；qzi表示轨道法线矢量在地心惯性坐标系中的z轴分量，q表示轨道的法线矢量，q表示运动方向，若Δθ<π，q＝+1，否则q＝-1；
采用如下公式计算所述升交点赤经Ω：
Figure FDA0003370291250000041
采用如下公式计算所述初始位置的纬度幅角ui和目标位置的纬度幅角uf
Figure FDA0003370291250000042
{Ri}n＝[{Ri}nx {Ri}ny {Ri}nz]T＝LniRi；
Figure FDA0003370291250000043
Figure FDA0003370291250000044
{Rf}n＝[{Rf}nx {Rf}ny {Rf}nz]T＝LniRf
转移轨道的近地点辐角ω：
ω＝ui-θi＝uf-θf
式中，{Ri}n表示初始位置矢量在地心节点坐标系中的表达式，{Ri}nx表示初始位置矢量在地心节点坐标系中的x轴方向分量；{Ri}ny表示初始位置矢量在地心节点坐标系中的y轴方向分量；{Ri}nz表示初始位置矢量在地心节点坐标系中的z轴方向分量；{Rf}n表示目标位置矢量在地心节点坐标系中的表达式，{Rf}nx表示目标位置矢量在地心节点坐标系中的x轴方向分量；{Rf}ny表示目标位置矢量在地心节点坐标系中的y轴方向分量；{Rf}nz表示目标位置矢量在地心节点坐标系中的z轴方向分量。
5.如权利要求4所述的基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，其特征在于，
所述步骤S2.3包括：初始速度矢量Vi与末端速度矢量Vf采用如下公式计算：
[Ri Vi]＝Ψ(a,e,i,ui,Ω,θi)
[Rf Vf]＝Ψ(a,e,i,uf,Ω,θf)
其中，Vi为初始位置的速度矢量，Vf为目标位置的速度矢量，Ψ是轨道六要素与轨道位置和速度状态量之间的转换关系。
6.如权利要求5所述的基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，其特征在于，
所述步骤S3包括：从初始状态[Ri，Vi]开始，数值积分时间T，得到实际末端位置矢量Rf'(Vi)，计算因J2摄动对多圈Lambert转移轨道的最终位置产生的偏差：
所述J2摄动下的动力学方程
Figure FDA0003370291250000051
为：
Figure FDA0003370291250000052
其中，R为位置矢量，R＝[x，y，z]T，x表示位置矢量在地心惯性坐标系中的x轴方向分量；y表示位置矢量在地心惯性坐标系中的y轴方向分量；z表示位置矢量在地心惯性坐标系中的z轴方向分量；R随时间变化，当飞行时间t＝0时，R＝Ri，表示初始位置矢量，当飞行时间等于转移时间t＝T时，R＝Rf’(Vi)表示最终位置矢量；
Figure FDA0003370291250000053
是位置矢量关于飞行时间的二阶导数；r为该位置到地心的距离，r＝||R||，μ＝GM＝3.986*1014[m3/s2]，其中G是万有引力常数G＝6.670*10-11[m3/(kg.s2)]；M是地球质量M＝5.976*1024[kg]，J2是地球扁率摄动系数，Re为地球半径；
因J2摄动对转移轨道最终位置产生的偏差ΔRf(Vi)为
ΔRf(Vi)＝||Rf'(Vi)-Rf||
其中，Rf'(Vi)为最终位置矢量。
7.如权利要求6所述的基于遗传算法修正的多圈J2-Lambert转移轨道求解方法，其特征在于，所述步骤S4包括：
步骤S4.1、遗传算法初始化：随机生成初始种群，所述初始种群为随机生成一组初始速度Vi可行解的一个集合[Vi1,Vi2,Vi3…Vi(n-1),Vin]，其中n是所述初始种群可行解的个数；
步骤S4.2、设置目标函数为ΔRf(Vi)＝||Rf'(Vi)-Rf||，取适应度函数
Figure FDA0003370291250000061
步骤S4.3、按照所述初始种群中各个个体适应度函数值的大小，从初始种群中选出适应度函数值前若干的个体构成交配池，并将这些个体复制到下一代种群；
步骤S4.4、由交叉和变异这两个遗传算子对交配池中的若干个体进行操作，形成新一代的种群中的m个个体，其中m是新一代种群中新个体的数目；
步骤S4.5、将所述步骤S4.3复制的交配池中(n-m)个体和步骤S4.4产生的m个个体组合成新一代种群；
步骤S4.6、反复执行步骤S4.2～步骤S4.5，直至最优个体的适应度不再上升时，或者迭代次数达到预设的代数时，遗传算法终止，输出最优初始速度。
8.一种电子设备，其特征在于，包括处理器和存储器，所述存储器上存储有计算机程序，所述计算机程序被所述处理器执行时，实现权利要求1至7中任一项所述的方法。
9.一种可读存储介质，其特征在于，所述可读存储介质内存储有计算机程序，所述计算机程序被处理器执行时，实现权利要求1至7中任一项所述的方法。
CN202111395792.3A 2021-11-23 2021-11-23 基于遗传算法修正的多圈J2-Lambert转移轨道求解方法 Active CN114357853B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202111395792.3A CN114357853B (zh) 2021-11-23 2021-11-23 基于遗传算法修正的多圈J2-Lambert转移轨道求解方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202111395792.3A CN114357853B (zh) 2021-11-23 2021-11-23 基于遗传算法修正的多圈J2-Lambert转移轨道求解方法

Publications (2)

Publication Number Publication Date
CN114357853A true CN114357853A (zh) 2022-04-15
CN114357853B CN114357853B (zh) 2024-10-15

Family

ID=81096049

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202111395792.3A Active CN114357853B (zh) 2021-11-23 2021-11-23 基于遗传算法修正的多圈J2-Lambert转移轨道求解方法

Country Status (1)

Country Link
CN (1) CN114357853B (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115113638A (zh) * 2022-07-20 2022-09-27 上海交通大学 一种燃料最优主动漂移三维成像轨道控制方法
CN115571379A (zh) * 2022-09-16 2023-01-06 中电普信(北京)科技发展有限公司 一种不调整轨道面的异面交会变轨方法及系统
CN116295452A (zh) * 2023-03-23 2023-06-23 哈尔滨工业大学 一种太阳光约束下单对多星连续掠飞轨迹优化方法及装置

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
RU2013101320A (ru) * 2013-01-10 2014-07-20 Открытое акционерное общество "Информационные спутниковые системы" имени академика М.Ф. Решетнева" Способ удержания космического аппарата на геосинхронной 24-часовой орбите
CN110399652A (zh) * 2019-07-05 2019-11-01 中国人民解放军国防科技大学 基于运动伪装的航天器隐蔽机动轨迹设计方法
CN113343442A (zh) * 2021-05-21 2021-09-03 中国科学院国家空间科学中心 一种求解固定时间有限燃料多脉冲转移轨道的方法及系统

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
RU2013101320A (ru) * 2013-01-10 2014-07-20 Открытое акционерное общество "Информационные спутниковые системы" имени академика М.Ф. Решетнева" Способ удержания космического аппарата на геосинхронной 24-часовой орбите
CN110399652A (zh) * 2019-07-05 2019-11-01 中国人民解放军国防科技大学 基于运动伪装的航天器隐蔽机动轨迹设计方法
CN113343442A (zh) * 2021-05-21 2021-09-03 中国科学院国家空间科学中心 一种求解固定时间有限燃料多脉冲转移轨道的方法及系统

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
WEI QIAN 等: "A Modified Algorithm on the Midcourse Guidance Based on BP Neural Network", JOURNAL OF XI\'AN JIAOTONG UNIVERSITY, 13 April 2017 (2017-04-13) *
余培军;王维;李建成;: "对地观测卫星有限推力快速轨道机动优化方法", 北京航空航天大学学报, no. 07, 15 July 2011 (2011-07-15) *
王有亮 等: "J2摄动下卫星编队重构的多脉冲轨迹优化", 国防科技大学学报, 28 December 2018 (2018-12-28) *

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115113638A (zh) * 2022-07-20 2022-09-27 上海交通大学 一种燃料最优主动漂移三维成像轨道控制方法
CN115571379A (zh) * 2022-09-16 2023-01-06 中电普信(北京)科技发展有限公司 一种不调整轨道面的异面交会变轨方法及系统
CN116295452A (zh) * 2023-03-23 2023-06-23 哈尔滨工业大学 一种太阳光约束下单对多星连续掠飞轨迹优化方法及装置
CN116295452B (zh) * 2023-03-23 2023-11-14 哈尔滨工业大学 一种太阳光约束下单对多星连续掠飞轨迹优化方法及装置

Also Published As

Publication number Publication date
CN114357853B (zh) 2024-10-15

Similar Documents

Publication Publication Date Title
CN114357853A (zh) 基于遗传算法修正的多圈J2-Lambert转移轨道求解方法
Ovchinnikov et al. Geomagnetic field models for satellite angular motion studies
CN109592079A (zh) 一种限定时间的航天器共面交会变轨策略确定方法
CN104527994A (zh) 异面交叉快变轨道固定时间稳定姿态指向跟踪控制方法
CN102819266B (zh) 一种拟周期j2不变相对轨道编队飞行控制方法
CN111444603B (zh) 一种返回式航天器时间最短离轨轨迹快速规划方法
CN110044210B (zh) 考虑任意阶地球非球型引力摄动的闭路制导在线补偿方法
CN109606739B (zh) 一种探测器地月转移轨道修正方法及装置
CN107506505B (zh) 高精度地月自由返回轨道设计方法
CN104252548A (zh) 一种燃料最优的火星探测器入轨目标点设计方法
Wu et al. Robust sliding mode control for spacecraft global fast-tracking manoeuvre
Sokolov et al. Peculiarities of the motion of asteroid 99942 Apophis
Petukhov Quasioptimal control with feedback for multiorbit low-thrust transfer between noncoplanar elliptic and circular orbits
Djojodihardjo Influence of the earth’s dominant oblateness parameter on the low formation orbits of micro-satellites
Ding et al. Research on minimum time interception problem with a tangent impulse under relative motion models
Zhang et al. Velocity-to-be-gained deorbit guidance law using state space perturbation method
Sun et al. Piecewise attitude tracking control of a gravity gradient microsatellite for coplanar orbital transfer
Ulybyshev Discrete pseudocontrol sets for optimal control problems
Silva A formulation of the Clohessy-Wiltshire equations to include dynamic atmospheric drag
Zhang et al. Minimizing the effects of impulsive orbit maneuver uncertainty
Giordano et al. Analysis and Optimization of Robust Trajectories in Cislunar Environment with Application to the LUMIO Cubesat
Pontani et al. Variable-time-domain neighboring optimal guidance and attitude control for low-thrust orbit transfers
Hou et al. New Set of Attitude Parameters Dedicated to Relative Attitude Maneuvers
Matsumoto et al. Accurate Real-Time Prediction Guidance Using Numerical Integration for Reentry Spacecraft
CN114383619A (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