1.本发明涉及航天器轨道设计与优化技术领域，具体地，涉及一种低轨卫星星座组网电推进变轨策略优化方法和系统。
背景技术：2.近年来，高覆盖率、短重访、甚至是全天时持续覆盖的低轨星座相继被提出，成为当前的研究热点。此类低轨星座的规模可达数十颗到数万颗，占据大量的空间轨位资源，且通常配置电推进系统，采用全电推的方式完成组网机动、构型保持、离轨等轨控任务。低轨卫星星座组网机动是指在星箭分离后，卫星从初始轨道变轨爬升至目标轨道的转移过程。与一般的单星任务相比，星座组网对轨道六根数均提出了设计约束，其转移轨道的设计也存在很大的差异。因此，全电推模式下的组网转移轨道设计与优化成了低轨星座发展的关键技术之一。
3.电推进轨道优化是一个经典的航天动力学设计问题，有很多的学者对其进行了研究，例如小行星探测转移轨道设计、人工太阳同步轨道设计、gto-geo多圈小推力转移、椭圆轨道卫星对地覆盖的电推进控制等，而对于低轨星座组网转移轨道设计的公开文献仍比较少。对于低轨卫星星座，为了降低星座中工作的卫星发生碰撞的概率，卫星的工作轨道常选择小偏心率冻结的近圆轨道。因此，低轨星座组网需要解决多个轨道根数联合控制问题，确保在变轨完成后相位差、半长轴、偏心率及近地点幅角均满足要求，且满足一定的性能指标。针对该问题，本专利考虑电推进推力大小固定，且仅有开和关两种状态的工程约束，构建时间最优的星座组网电推进转移轨道优化问题，提出一种推力参数化方法，将原问题转化为一个控制变量受约束的打靶问题，并用非线性规划算法对其进行求解。该方法可以实现相位差、半长轴、偏心率及近地点幅角四者的联合控制，适用于目标轨道是小偏心率冻结的低轨星座任务。最后，采用地球21阶次引力场模型对该方法的有效性进行了校验。
4.目前，人们对航天器变轨策略进行了一些研究，经检索，与变轨策略设计相关的主要有：
5.专利《一种静止轨道卫星变轨策略优化方法》(公开号：cn102424116a)克服现有技术的不足，提供了一种静止轨道卫星变轨策略优化方法，合理确定变轨策略设计的各个约束条件，以减少变轨策略设计过程中的人工干预以及计算时间和计算量。该专利关注的静止轨道卫星的变轨策略，模型与优化方法与低轨卫星组网存在很大的差异，其相关设计方法无法借用。
6.专利《基于粒子群算法的geo卫星变轨策略计算方法、系统及介质》(公开号：cn108216687a)也是针对静止轨道卫星变轨，将点火时刻和点火方向作为优化变量，将推进剂消耗作为目标函数，通过设定初始粒子种群并按算法进行进化计算，更快速的获得期望解，提高计算效率。该专利同样关注的静止轨道卫星的变轨策略，模型与优化方法与低轨卫星组网存在很大的差异，其相关设计方法无法借用。
7.专利《基于异面变轨策略的卫星星座重构方法、设备及存储设备》(公开号：
cn107885917a)提供了基于异面变轨策略的卫星星座重构方法、设备及存储设备，通过卫星机动变轨的方式，优化卫星的轨道参数，改变卫星星座的空间构型，对卫星星座进行重构，从而满足应急任务对卫星网络组网的性能需求。该专利关注的星座重构方法，而本专利则关注卫星发射入轨初期的组网机动，设计对象不同。
8.专利《一种卫星的自主变轨方法》(公开号：cn101219713)克服现有技术的不足，提出了一种星上自主进行的姿态确定，姿态机动以及自主轨控的自主变轨方法，解决了深空探测卫星的变轨问题，保证变轨准确、可靠、自动地完成。该专利关注深空探测航天器变轨的自主性，与本专利研究对象不同。
9.论文《连续小推力条件下星座轨道机动方法研究》(中国空间科学技术，2020年第01期：1-16)研究了一箭多星星座任务的组网轨道机动问题，采用分时升轨机动实现了相位差和半长轴的联合调整，但并未涉及半长轴、偏心率及近地点幅角三者的联合调整，因而无法直接适用于目标轨道为小偏心率冻结的近圆轨道组网任务。而本专利主要关注的是半长轴、偏心率及近地点幅角三者的联合调整。
10.专利文献cn105799954a(申请号：cn201410845134.3)公开了一种天基分散部署微纳载荷的模块化飞行器及其变轨制导方法，所述飞行器包括：承载与推进一体化单元、测量与控制一体化单元、轻小型级间适配分离装置；可以实现远程自主快速变轨到目标轨道。然而该专利仅对飞行器的变轨进行了研究，与本专利研究对象不同。
技术实现要素：11.针对现有技术中的缺陷，本发明的目的是提供一种低轨卫星星座组网电推进变轨策略优化方法和系统。
12.根据本发明提供的低轨卫星星座组网电推进变轨策略优化方法，包括：
13.步骤1：根据电推进的推力大小，在仅有开和关两种状态的工程约束的基础上，建立时间最优的星座组网电推进转移轨道优化问题；
14.步骤2：对组网过程进行阶段划分，并将推力大小及方向进行分阶段参数化处理；
15.步骤3：确定待求解变量，将星座组网电推进转移轨道优化问题转换为控制变量受约束的打靶问题；
16.步骤4：对打靶问题的初值进行猜测，采用牛顿迭代算法进行求解，最终得到星座组网电推进转移轨道优化问题的最优解。
17.优选的，所述步骤1包括：
18.步骤1.1：在地心惯性参考坐标系下，建立卫星的轨道动力学方程：
[0019][0020][0021][0022]
其中，r和v分表为位置和速度矢量；f
p
为包含j2、大气阻力的摄动加速度；μ为地球引力常数；r为地球半径；ge为地球海平面重力加速度；m为航天器质量；i
sp
为推力器比冲；α为推力方向的单位矢量；t
max
为最大推力；u为实际推力与最大推力的比值；
[0023]
在仅有开和关两种状态的工程约束的基础上，u的取值为0或者1；
[0024]
推力方向的单位矢量在惯性坐标系下进一步表示为：
[0025][0026]
其中，为推力方位角，δ为推力高度角；
[0027]
令卫星的瞬时轨道根数为：
[0028]kosc
＝[a e i ω ω m]
[0029]
其中，a为半长轴，e为偏心率，i为轨道倾角，ω为升交点赤经，ω为近地点幅角，m为平近点角；
[0030]
令表示卫星的平均轨道根数，x＝[r v]为直角坐标参数，则k
osc
、k
ave
及x三者可互相转换；
[0031]
步骤1.2：令初始时刻为t0，末端时刻为tf，则初末时刻卫星面内的组网约束表示为：
[0032]
r(t0)＝r0，v(t0)＝v0，m(t0)＝m0[0033][0034][0035]
式中，及为常值，表示末端时刻期望的半长轴、偏心率及近地点幅角，表示沿迹角，定义为：
[0036][0037]
在末端时刻，期望的根据如下公式确定：
[0038][0039]
式中，表示t0时刻的初始，表示沿迹角变化率；
[0040]
步骤1.3：令控制变量和状态变量：
[0041]
u＝[u α]
t
[0042]
y＝[r v m]
t
[0043]
则轨道动力学方程表示为：
[0044][0045]
相应的边界约束表示为：
[0046]bini
[y(t0)]＝y(t0)-y0＝0
[0047]
[0048]
其中，y0＝[r
0 v
0 m0]
t
[0049]
组网过程电推进转移轨道问题p为寻找最优控制变量：
[0050][0051]
同时满足：
[0052][0053]bini
[y(t0)]＝0，b
ter
[y(tf)]＝0
[0054]
u＝0或者u＝1，α
t
α＝1。
[0055]
优选的，所述步骤2包括：
[0056]
步骤2.1：根据最快入轨的性能指标要求和组网轨道的运动特性，将u和α在时间域上按时刻t1、t2及t3进行划分，在推力为全开的基础上，将转移轨道按时间顺序分为巡航阶段、升轨阶段及偏心率矢量调整阶段；
[0057]
步骤2.2：将各阶段的推力大小和方向设定为：
[0058][0059][0060]
式中，αa和α
eω
分别为升轨阶段和偏心率矢量调整阶段的推力方向单位矢量；
[0061]
步骤2.3：为满足最快入轨的性能指标要求，在升轨阶段和偏心率矢量调整阶段，将推力方向分别设定沿着切向方向和轨道面内垂直拱线方向，表示为：
[0062][0063]
α
eω
＝(cosσ)i
p
+(sinσ)j
p
[0064]
其中，σ表示推力方向与拱线的夹角，且：
[0065][0066][0067]
其中，i
p
和j
p
的计算是基于平均轨道参数进行的，确保推力方向平缓变化。
[0068]
优选的，所述步骤3包括：
[0069]
步骤3.1：对于给定的初始状态t0和y0，定义待求解的变量为：
[0070]
z＝[δt
1 δt
2 δt
3 σ]
t
[0071]
其中，δt1＝t
1-t0，δt2＝t
2-t1，δt3＝t
3-t2；
[0072]
步骤3.2：将问题p转换为如下打靶问题：
[0073][0074]
优选的，所述步骤4包括：
[0075]
步骤4.1：对δt1、δt2、δt3及σ进行猜测，确定迭代变量初值；
[0076]
步骤4.2：采用龙格库塔法计算代价函数f(z)，采用牛顿迭代算法进行迭代计算，直至||f(z)||收敛至小于指定阈值。
[0077]
根据本发明提供的低轨卫星星座组网电推进变轨策略优化系统，包括：
[0078]
模块m1：根据电推进的推力大小，在仅有开和关两种状态的工程约束的基础上，建立时间最优的星座组网电推进转移轨道优化问题；
[0079]
模块m2：对组网过程进行阶段划分，并将推力大小及方向进行分阶段参数化处理；
[0080]
模块m3：确定待求解变量，将星座组网电推进转移轨道优化问题转换为控制变量受约束的打靶问题；
[0081]
模块m4：对打靶问题的初值进行猜测，采用牛顿迭代算法进行求解，最终得到星座组网电推进转移轨道优化问题的最优解。
[0082]
优选的，所述模块m1包括：
[0083]
模块m1.1：在地心惯性参考坐标系下，建立卫星的轨道动力学方程：
[0084][0085][0086][0087]
其中，r和v分表为位置和速度矢量；f
p
为包含j2、大气阻力的摄动加速度；μ为地球引力常数；r为地球半径；ge为地球海平面重力加速度；m为航天器质量；i
sp
为推力器比冲；α为推力方向的单位矢量；t
max
为最大推力；u为实际推力与最大推力的比值；
[0088]
在仅有开和关两种状态的工程约束的基础上，u的取值为0或者1；
[0089]
推力方向的单位矢量在惯性坐标系下进一步表示为：
[0090][0091]
其中，为推力方位角，δ为推力高度角；
[0092]
令卫星的瞬时轨道根数为：
[0093]kosc
＝[a e i ω ω m]
[0094]
其中，a为半长轴，e为偏心率，i为轨道倾角，ω为升交点赤经，ω为近地点幅角，m为平近点角；
[0095]
令表示卫星的平均轨道根数，x＝[r v]为直角坐标参数，则k
osc
、k
ave
及x三者可互相转换；
[0096]
模块m1.2：令初始时刻为t0，末端时刻为tf，则初末时刻卫星面内的组网约束表示为：
[0097]
r(t0)＝r0，v(t0)＝v0，m(t0)＝m0[0098][0099][0100]
式中，及为常值，表示末端时刻期望的半长轴、偏心率及近地点幅角，表示沿迹角，定义为：
[0101][0102]
在末端时刻，期望的根据如下公式确定：
[0103][0104]
式中，表示t0时刻的初始，表示沿迹角变化率；
[0105]
模块m1.3：令控制变量和状态变量：
[0106]
u＝[u α]
t
[0107]
y＝[r v m]
t
[0108]
则轨道动力学方程表示为：
[0109][0110]
相应的边界约束表示为：
[0111]bini
[y(t0)]＝y(t0)-y0＝0
[0112][0113]
其中，y0＝[r
0 v
0 m0]
t
[0114]
组网过程电推进转移轨道问题p为寻找最优控制变量：
[0115][0116]
同时满足：
[0117][0118]bini
[y(t0)]＝0，b
ter
[y(tf)]＝0
[0119]
u＝0或者u＝1，α
t
α＝1。
[0120]
优选的，所述模块m2包括：
[0121]
模块m2.1：根据最快入轨的性能指标要求和组网轨道的运动特性，将u和α在时间域上按时刻t1、t2及t3进行划分，在推力为全开的基础上，将转移轨道按时间顺序分为巡航阶段、升轨阶段及偏心率矢量调整阶段；
[0122]
模块m2.2：将各阶段的推力大小和方向设定为：
[0123][0124][0125]
式中，αa和α
eω
分别为升轨阶段和偏心率矢量调整阶段的推力方向单位矢量；
[0126]
模块m2.3：为满足最快入轨的性能指标要求，在升轨阶段和偏心率矢量调整阶段，将推力方向分别设定沿着切向方向和轨道面内垂直拱线方向，表示为：
[0127][0128]
α
eω
＝(cosσ)i
p
+(sinσ)j
p
[0129]
其中，σ表示推力方向与拱线的夹角，且：
[0130][0131][0132]
其中，i
p
和j
p
的计算是基于平均轨道参数进行的，确保推力方向平缓变化。
[0133]
优选的，所述模块m3包括：
[0134]
模块m3.1：对于给定的初始状态t0和y0，定义待求解的变量为：
[0135]
z＝[δt
1 δt
2 δt
3 σ]
t
[0136]
其中，δt1＝t
1-t0，δt2＝t
2-t1，δt3＝t
3-t2；
[0137]
模块m3.2：将问题p转换为如下打靶问题：
[0138][0139]
优选的，所述模块m4包括：
[0140]
模块m4.1：对δt1、δt2、δt3及σ进行猜测，确定迭代变量初值；
[0141]
模块m4.2：采用龙格库塔法计算代价函数f(z)，采用牛顿迭代算法进行迭代计算，直至||f(z)||收敛至小于指定阈值。
[0142]
与现有技术相比，本发明具有如下的有益效果：
[0143]
本发明中的方法考虑电推进推力大小固定，仅有开和关两种状态，根据最快入轨要求，充分利于组网轨道的运动特性，提出一种推力参数化方法，可以保证转移轨道满足预先设定的组网约束，对姿控系统要求低，且推力参数仅为4个，结构简单，求解效率高，有利
于星上实现，可以有效实现轨道根数的联合控制。
附图说明
[0144]
通过阅读参照以下附图对非限制性实施例所作的详细描述，本发明的其它特征、目的和优点将会变得更明显：
[0145]
图1为本发明方法流程图；
[0146]
图2为具体实施方式中推力剖面；
[0147]
图3为具体实施方式中半长轴及沿迹角偏差；
[0148]
图4为具体实施方式中偏心率矢量及推力大小；
[0149]
图5为具体实施方式中惯性坐标系下推力方向角及燃料消耗。
具体实施方式
[0150]
下面结合具体实施例对本发明进行详细说明。以下实施例将有助于本领域的技术人员进一步理解本发明，但不以任何形式限制本发明。应当指出的是，对本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变化和改进。这些都属于本发明的保护范围。
[0151]
实施例：
[0152]
本发明提供一种低轨卫星星座组网电推进变轨策略优化算法，考虑电推进推力大小固定，仅有开和关两种状态，根据最快入轨要求，充分利于组网轨道的运动特性，提出一种推力参数化方法，可以保证转移轨道满足预先设定的组网约束，对姿控系统要求低，且推力参数仅为4个，结构简单，求解效率高，有利于星上实现，可以有效实现轨道根数的联合控制。
[0153]
如图1，本发明提供了一种低轨卫星星座组网电推进变轨策略优化算法，包括如下步骤：
[0154]
步骤1：考虑电推进推力大小固定，仅有开和关两种状态的工程约束，建立时间最优的星座组网电推进转移轨道优化问题。
[0155]
其中，具体过程包括：
[0156]
步骤101：在地心惯性参考坐标系下，建立卫星的轨道动力学方程为：
[0157][0158][0159]
其中，r和v分表为位置和速度矢量，f
p
为包含j2、大气阻力等摄动加速度，μ为地球引力常数，r为地球半径，ge为地球海平面重力加速度，m为航天器质量，i
sp
为推力器比冲，α为推力方向的单位矢量，t
max
为最大推力，u为实际推力与最大推力的比值。本文考虑发动机推力不可连续调整，仅有开和关两种状态的情况，因此u的取值为0或者1。推力方向的单位矢量在惯性坐标系下可进一步表示为：
[0160][0161]
其中，为推力方位角，δ为推力高度角。
[0162]
令卫星的瞬时轨道根数为：
[0163]kosc
＝[a e i ω ω m]
[0164]
其中，a为半长轴，e为偏心率，i为轨道倾角，ω为升交点赤经，ω为近地点幅角，m为平近点角。令表示卫星的平均轨道根数，x＝[r v]为直角坐标参数，则k
osc
、k
ave
及x三者可互相转换。
[0165]
步骤102：令初始时刻为t0，末端时刻为tf，则初末时刻卫星面内的组网约束可以表示为：
[0166]
r(t0)＝r0，v(t0)＝v0，m(t0)＝m0[0167][0168][0169]
式中，及为常值，表示末端时刻期望的半长轴、偏心率及近地点幅角，表示沿迹角，定义为：
[0170][0171]
在末端时刻，期望的可根据如下式子确定：
[0172][0173]
式中，表示t0时刻的初始，表示沿迹角变化率。
[0174]
步骤103：为便于阐述，令控制变量和状态变量：
[0175]
u＝[u α]
t
[0176]
y＝[r v m]
t
[0177]
则轨道动力学方程可以表示为：
[0178][0179]
相应的边界约束表示为：
[0180]bini
[y(t0)]＝y(t0)-y0＝0
[0181][0182]
其中，y0＝[r
0 v
0 m0]
t
。
[0183]
组网过程电推进转移轨道问题p可以阐述如下：寻找最优控制变量；
[0184][0185]
同时满足：
[0186][0187]bini
[y(t0)]＝0，b
ter
[y(tf)]＝0
[0188]
u＝0或者1，α
t
α＝1
[0189]
步骤2：对组网过程进行阶段划分，并将推力大小及方向分阶段参数化处理。
[0190]
其中，步骤2具体包括：
[0191]
步骤201：根据最快入轨的性能指标要求，充分考虑组网轨道的运动特性，将u和α在时间域上按时刻t1、t2及t3进行划分，在推力基本为全开的基础上，将转移轨道按时间顺序分为巡航阶段、升轨阶段及偏心率矢量调整阶段，如图2。
[0192]
步骤202：将各阶段的推力大小和方向设定为：
[0193][0194][0195]
式中，αa和α
eω
分别为升轨阶段和偏心率矢量调整阶段的推力方向单位矢量。
[0196]
步骤203：为满足最快入轨的性能指标要求，在升轨阶段和偏心率矢量调整阶段，将推力方向分别设定沿着切向方向和轨道面内垂直拱线方向，具体表示为：
[0197][0198]
α
eω
＝(cosσ)i
p
+(sinσ)j
p
[0199]
其中，σ表示推力方向与拱线的夹角，且：
[0200][0201][0202]
需要指出的是，为确保推力方向平缓变化，i
p
和j
p
的计算是基于平均轨道参数进行的。
[0203]
步骤3：确定待求解变量，建立一个控制变量受约束的打靶问题。
[0204]
该步骤的具体过程为：
[0205]
步骤301：经过步骤2的转换，对于给定的初始状态t0和y0，定义待求解的变量为
[0206]
z＝[δt
1 δt
2 δt
3 σ]
t
[0207]
其中，δt1＝t
1-t0，δt2＝t
2-t1，δt3＝t
3-t2。
[0208]
步骤302：将问题p转换为如下打靶问题：
[0209][0210]
可以看出，上述方程组中，未知变量的个数与约束方程数量相等，因此方程的解是完全确定的。
[0211]
步骤4：对打靶问题的初值进行猜测，采用牛顿迭代算法进行求解。
[0212]
该步骤的具体过程为：
[0213]
步骤401：对δt1、δt2、δt3及σ进行猜测，确定迭代变量初值；
[0214]
步骤402：采用龙格-库塔数值计算算法计算代价函数f(z)，采用牛顿迭代算法进行迭代计算，直至||f(z)||收敛至小于指定阈值。
[0215]
以下为一种低轨卫星星座组网电推进变轨策略优化算法的数值仿真验证。
[0216]
假设卫星采用偏低发射入轨，星间分离后卫星的初始轨道高度为1000km，目标轨道高度为1200km。假设卫星的初始质量为170kg，推力大小为20mn，比冲为1600s。在表中，将目标轨道的近地点幅角和偏心率选择为90deg和0.001是为了实现目标轨道的冻结。考虑目标轨道的倾角接近90deg，其轨道冻结受地球高阶引力场影响大，为此在仿真中，将考虑地球21阶次的引力场模型。由于初始轨道偏心率接近0，近地点幅角无法定义，为此在后续分析中采用偏心率矢量的两个分量来描述，即：
[0217]ex
＝ecosω，ey＝esinω
[0218]
在目标轨道，e
x
和ey的目标值为0和0.001。在仿真计算结果和分析中，若无特别说明，则所有的轨道根数均是指平均轨道根数。
[0219]
根据以上计算条件，进行数值仿真计算，结果如图3、图4、图5所示。
[0220]
根据本发明提供的低轨卫星星座组网电推进变轨策略优化系统，包括：模块m1：根据电推进的推力大小，在仅有开和关两种状态的工程约束的基础上，建立时间最优的星座组网电推进转移轨道优化问题；模块m2：对组网过程进行阶段划分，并将推力大小及方向进行分阶段参数化处理；模块m3：确定待求解变量，将星座组网电推进转移轨道优化问题转换为控制变量受约束的打靶问题；模块m4：对打靶问题的初值进行猜测，采用牛顿迭代算法进行求解，最终得到星座组网电推进转移轨道优化问题的最优解。
[0221]
所述模块m1包括：模块m1.1：在地心惯性参考坐标系下，建立卫星的轨道动力学方程：
[0222][0223][0224][0225]
其中，r和v分表为位置和速度矢量；f
p
为包含j2、大气阻力的摄动加速度；μ为地球引力常数；r为地球半径；ge为地球海平面重力加速度；m为航天器质量；i
sp
为推力器比冲；α为推力方向的单位矢量；t
max
为最大推力；u为实际推力与最大推力的比值；在仅有开和关两种状态的工程约束的基础上，u的取值为0或者1；推力方向的单位矢量在惯性坐标系下进一
步表示为：其中，为推力方位角，δ为推力高度角；令卫星的瞬时轨道根数为：k
osc
＝[a e i ω ω m]，其中，a为半长轴，e为偏心率，i为轨道倾角，ω为升交点赤经，ω为近地点幅角，m为平近点角；令表示卫星的平均轨道根数，x＝[r v]为直角坐标参数，则k
osc
、k
ave
及x三者可互相转换；
[0226]
模块m1.2：令初始时刻为t0，末端时刻为tf，则初末时刻卫星面内的组网约束表示为：r(t0)＝r0，v(t0)＝v0，m(t0)＝m0[0227][0228][0229]
式中，及为常值，表示末端时刻期望的半长轴、偏心率及近地点幅角，表示沿迹角，定义为：
[0230]
在末端时刻，期望的根据如下公式确定：
[0231]
式中，表示t0时刻的初始，表示沿迹角变化率；
[0232]
模块m1.3：令控制变量和状态变量：u＝[u α]
t
，y＝[r v m]
t
[0233]
则轨道动力学方程表示为：
[0234]
相应的边界约束表示为：b
ini
[y(t0)]＝y(t0)-y0＝0
[0235][0236]
其中，y0＝[r
0 v
0 m0]
t
[0237]
组网过程电推进转移轨道问题p为寻找最优控制变量：
[0238]
同时满足：b
ini
[y(t0)]＝0，b
ter
[y(tf)]＝0，u＝0或者u＝1，α
t
α＝1。
[0239]
所述模块m2包括：模块m2.1：根据最快入轨的性能指标要求和组网轨道的运动特性，将u和α在时间域上按时刻t1、t2及t3进行划分，在推力为全开的基础上，将转移轨道按时间顺序分为巡航阶段、升轨阶段及偏心率矢量调整阶段；模块m2.2：将各阶段的推力大小和方向设定为：式中，αa和α
eω
分别为升轨阶段和偏心率矢量调整阶段的推力方向单位矢量；模块m2.3：为满足最快入轨的性能指标要求，在升轨阶段和偏心率矢量调整阶段，将推力方向分别设定沿着切向方向和轨道面内垂直拱线方
向，表示为：α
eω
＝(cosσ)i
p
+(sinσ)j
p
，其中，σ表示推力方向与拱线的夹角，且：其中，i
p
和j
p
的计算是基于平均轨道参数进行的，确保推力方向平缓变化。
[0240]
所述模块m3包括：模块m3.1：对于给定的初始状态t0和y0，定义待求解的变量为：z＝[δt
1 δt
2 δt
3 σ]
t
，其中，δt1＝t
1-t0，δt2＝t
2-t1，δt3＝t
3-t2；
[0241]
模块m3.2：将问题p转换为如下打靶问题：
[0242][0243]
所述模块m4包括：模块m4.1：对δt1、δt2、δt3及σ进行猜测，确定迭代变量初值；模块m4.2：采用龙格库塔法计算代价函数f(z)，采用牛顿迭代算法进行迭代计算，直至||f(z)||收敛至小于指定阈值。
[0244]
本领域技术人员知道，除了以纯计算机可读程序代码方式实现本发明提供的系统、装置及其各个模块以外，完全可以通过将方法步骤进行逻辑编程来使得本发明提供的系统、装置及其各个模块以逻辑门、开关、专用集成电路、可编程逻辑控制器以及嵌入式微控制器等的形式来实现相同程序。所以，本发明提供的系统、装置及其各个模块可以被认为是一种硬件部件，而对其内包括的用于实现各种程序的模块也可以视为硬件部件内的结构；也可以将用于实现各种功能的模块视为既可以是实现方法的软件程序又可以是硬件部件内的结构。
[0245]
以上对本发明的具体实施例进行了描述。需要理解的是，本发明并不局限于上述特定实施方式，本领域技术人员可以在权利要求的范围内做出各种变化或修改，这并不影响本发明的实质内容。在不冲突的情况下，本技术的实施例和实施例中的特征可以任意相互组合。
技术特征：1.一种低轨卫星星座组网电推进变轨策略优化方法，其特征在于，包括：步骤1：根据电推进的推力大小，在仅有开和关两种状态的工程约束的基础上，建立时间最优的星座组网电推进转移轨道优化问题；步骤2：对组网过程进行阶段划分，并将推力大小及方向进行分阶段参数化处理；步骤3：确定待求解变量，将星座组网电推进转移轨道优化问题转换为控制变量受约束的打靶问题；步骤4：对打靶问题的初值进行猜测，采用牛顿迭代算法进行求解，最终得到星座组网电推进转移轨道优化问题的最优解。2.根据权利要求1所述的低轨卫星星座组网电推进变轨策略优化方法，其特征在于，所述步骤1包括：步骤1.1：在地心惯性参考坐标系下，建立卫星的轨道动力学方程：步骤1.1：在地心惯性参考坐标系下，建立卫星的轨道动力学方程：步骤1.1：在地心惯性参考坐标系下，建立卫星的轨道动力学方程：其中，r和v分表为位置和速度矢量；f
p
为包含j2、大气阻力的摄动加速度；μ为地球引力常数；r为地球半径；g
e
为地球海平面重力加速度；m为航天器质量；i
sp
为推力器比冲；α为推力方向的单位矢量；t
max
为最大推力；u为实际推力与最大推力的比值；在仅有开和关两种状态的工程约束的基础上，u的取值为0或者1；推力方向的单位矢量在惯性坐标系下进一步表示为：其中，为推力方位角，δ为推力高度角；令卫星的瞬时轨道根数为：k
osc
＝[a e i ω ω m]其中，a为半长轴，e为偏心率，i为轨道倾角，ω为升交点赤经，ω为近地点幅角，m为平近点角；令表示卫星的平均轨道根数，x＝[r v]为直角坐标参数，则k
osc
、k
ave
及x三者可互相转换；步骤1.2：令初始时刻为t0，末端时刻为t
f
，则初末时刻卫星面内的组网约束表示为：r(t0)＝r0，v(t0)＝v0，m(t0)＝m
00
式中，及为常值，表示末端时刻期望的半长轴、偏心率及近地点幅角，表示沿迹角，定义为：在末端时刻，期望的根据如下公式确定：
式中，表示t0时刻的初始，表示沿迹角变化率；步骤1.3：令控制变量和状态变量：u＝[u α]
t
y＝[r v m]
t
则轨道动力学方程表示为：相应的边界约束表示为：b
ini
[y(t0)]＝y(t0)-y0＝0其中，y0＝[r
0 v
0 m0]
t
组网过程电推进转移轨道问题p为寻找最优控制变量：同时满足：b
ini
[y(t0)]＝0，b
ter
[y(t
f
)]＝0u＝0或者u＝1，α
t
α＝1。3.根据权利要求2所述的低轨卫星星座组网电推进变轨策略优化方法，其特征在于，所述步骤2包括：步骤2.1：根据最快入轨的性能指标要求和组网轨道的运动特性，将u和α在时间域上按时刻t1、t2及t3进行划分，在推力为全开的基础上，将转移轨道按时间顺序分为巡航阶段、升轨阶段及偏心率矢量调整阶段；步骤2.2：将各阶段的推力大小和方向设定为：步骤2.2：将各阶段的推力大小和方向设定为：式中，α
a
和α
eω
分别为升轨阶段和偏心率矢量调整阶段的推力方向单位矢量；
步骤2.3：为满足最快入轨的性能指标要求，在升轨阶段和偏心率矢量调整阶段，将推力方向分别设定沿着切向方向和轨道面内垂直拱线方向，表示为：α
eω
＝(cosσ)i
p
+(sinσ)j
p
其中，σ表示推力方向与拱线的夹角，且：其中，σ表示推力方向与拱线的夹角，且：其中，i
p
和j
p
的计算是基于平均轨道参数进行的，确保推力方向平缓变化。4.根据权利要求3所述的低轨卫星星座组网电推进变轨策略优化方法，其特征在于，所述步骤3包括：步骤3.1：对于给定的初始状态t0和y0，定义待求解的变量为：z＝[δt
1 δt
2 δt
3 σ]
t
其中，δt1＝t
1-t0，δt2＝t
2-t1，δt3＝t
3-t2；步骤3.2：将问题p转换为如下打靶问题：5.根据权利要求4所述的低轨卫星星座组网电推进变轨策略优化方法，其特征在于，所述步骤4包括：步骤4.1：对δt1、δt2、δt3及σ进行猜测，确定迭代变量初值；步骤4.2：采用龙格库塔法计算代价函数f(z)，采用牛顿迭代算法进行迭代计算，直至||f(z)||收敛至小于指定阈值。6.一种低轨卫星星座组网电推进变轨策略优化系统，其特征在于，包括：模块m1：根据电推进的推力大小，在仅有开和关两种状态的工程约束的基础上，建立时间最优的星座组网电推进转移轨道优化问题；模块m2：对组网过程进行阶段划分，并将推力大小及方向进行分阶段参数化处理；模块m3：确定待求解变量，将星座组网电推进转移轨道优化问题转换为控制变量受约束的打靶问题；模块m4：对打靶问题的初值进行猜测，采用牛顿迭代算法进行求解，最终得到星座组网电推进转移轨道优化问题的最优解。7.根据权利要求6所述的低轨卫星星座组网电推进变轨策略优化系统，其特征在于，所述模块m1包括：
模块m1.1：在地心惯性参考坐标系下，建立卫星的轨道动力学方程：模块m1.1：在地心惯性参考坐标系下，建立卫星的轨道动力学方程：模块m1.1：在地心惯性参考坐标系下，建立卫星的轨道动力学方程：其中，r和v分表为位置和速度矢量；f
p
为包含j2、大气阻力的摄动加速度；μ为地球引力常数；r为地球半径；g
e
为地球海平面重力加速度；m为航天器质量；i
sp
为推力器比冲；α为推力方向的单位矢量；t
max
为最大推力；u为实际推力与最大推力的比值；在仅有开和关两种状态的工程约束的基础上，u的取值为0或者1；推力方向的单位矢量在惯性坐标系下进一步表示为：其中，为推力方位角，δ为推力高度角；令卫星的瞬时轨道根数为：k
osc
＝[a e i ω ω m]其中，a为半长轴，e为偏心率，i为轨道倾角，ω为升交点赤经，ω为近地点幅角，m为平近点角；令表示卫星的平均轨道根数，x＝[r v]为直角坐标参数，则k
osc
、k
ave
及x三者可互相转换；模块m1.2：令初始时刻为t0，末端时刻为t
f
，则初末时刻卫星面内的组网约束表示为：r(t0)＝r0，v(t0)＝v0，m(t0)＝m
00
式中，及为常值，表示末端时刻期望的半长轴、偏心率及近地点幅角，表示沿迹角，定义为：在末端时刻，期望的根据如下公式确定：式中，表示t0时刻的初始，表示沿迹角变化率；模块m1.3：令控制变量和状态变量：u＝[u α]
t
y＝[r v m]
t
则轨道动力学方程表示为：
相应的边界约束表示为：b
ini
[y(t0)]＝y(t0)-y0＝0其中，y0＝[r
0 v
0 m0]
t
组网过程电推进转移轨道问题p为寻找最优控制变量：同时满足：b
ini
[y(t0)]＝0，b
ter
[y(t
f
)]＝0u＝0或者u＝1，α
t
α＝1。8.根据权利要求7所述的低轨卫星星座组网电推进变轨策略优化系统，其特征在于，所述模块m2包括：模块m2.1：根据最快入轨的性能指标要求和组网轨道的运动特性，将u和α在时间域上按时刻t1、t2及t3进行划分，在推力为全开的基础上，将转移轨道按时间顺序分为巡航阶段、升轨阶段及偏心率矢量调整阶段；模块m2.2：将各阶段的推力大小和方向设定为：模块m2.2：将各阶段的推力大小和方向设定为：式中，α
a
和α
eω
分别为升轨阶段和偏心率矢量调整阶段的推力方向单位矢量；模块m2.3：为满足最快入轨的性能指标要求，在升轨阶段和偏心率矢量调整阶段，将推力方向分别设定沿着切向方向和轨道面内垂直拱线方向，表示为：α
eω
＝(cosσ)i
p
+(sinσ)j
p
其中，σ表示推力方向与拱线的夹角，且：
其中，i
p
和j
p
的计算是基于平均轨道参数进行的，确保推力方向平缓变化。9.根据权利要求8所述的低轨卫星星座组网电推进变轨策略优化系统，其特征在于，所述模块m3包括：模块m3.1：对于给定的初始状态t0和y0，定义待求解的变量为：z＝[δt
1 δt
2 δt
3 σ]
t
其中，δt1＝t
1-t0，δt2＝t
2-t1，δt3＝t
3-t2；模块m3.2：将问题p转换为如下打靶问题：10.根据权利要求9所述的低轨卫星星座组网电推进变轨策略优化系统，其特征在于，所述模块m4包括：模块m4.1：对δt1、δt2、δt3及σ进行猜测，确定迭代变量初值；模块m4.2：采用龙格库塔法计算代价函数f(z)，采用牛顿迭代算法进行迭代计算，直至||f(z)||收敛至小于指定阈值。
技术总结本发明提供了一种低轨卫星星座组网电推进变轨策略优化方法和系统，包括：步骤1：根据电推进的推力大小，在仅有开和关两种状态的工程约束的基础上，建立时间最优的星座组网电推进转移轨道优化问题；步骤2：对组网过程进行阶段划分，并将推力大小及方向进行分阶段参数化处理；步骤3：确定待求解变量，将星座组网电推进转移轨道优化问题转换为控制变量受约束的打靶问题；步骤4：对打靶问题的初值进行猜测，采用牛顿迭代算法进行求解，最终得到星座组网电推进转移轨道优化问题的最优解。本发明可以保证转移轨道满足预先设定的组网约束，对姿控系统要求低，结构简单，求解效率高，有利于星上实现，可以有效实现轨道根数的联合控制。可以有效实现轨道根数的联合控制。可以有效实现轨道根数的联合控制。
技术研发人员：郑艺裕 邹兴 袁双 杨先睿 刘培 周必磊
受保护的技术使用者：上海卫星工程研究所
技术研发日：2022.07.20
技术公布日：2022/11/1
