CN109062243B - 一种多约束下航天器能量最优的自主碰撞规避方法 - Google Patents

一种多约束下航天器能量最优的自主碰撞规避方法 Download PDF

Info

Publication number
CN109062243B
CN109062243B CN201811285216.1A CN201811285216A CN109062243B CN 109062243 B CN109062243 B CN 109062243B CN 201811285216 A CN201811285216 A CN 201811285216A CN 109062243 B CN109062243 B CN 109062243B
Authority
CN
China
Prior art keywords
satellite
axis
coordinate system
equation
collision avoidance
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.)
Expired - Fee Related
Application number
CN201811285216.1A
Other languages
English (en)
Other versions
CN109062243A (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.)
Changchun University of Technology
Original Assignee
Changchun University of Technology
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 Changchun University of Technology filed Critical Changchun University of Technology
Priority to CN201811285216.1A priority Critical patent/CN109062243B/zh
Publication of CN109062243A publication Critical patent/CN109062243A/zh
Application granted granted Critical
Publication of CN109062243B publication Critical patent/CN109062243B/zh
Expired - Fee Related 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
    • B64G1/244—Spacecraft control systems

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Automation & Control Theory (AREA)
  • Chemical & Material Sciences (AREA)
  • Combustion & Propulsion (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Control Of Position, Course, Altitude, Or Attitude Of Moving Bodies (AREA)

Abstract

一种多约束下航天器能量最优的自主碰撞规避方法，涉及航空航天领域，有效解决了卫星碰撞问题。本发明包括：建立卫星运动学模型；建立坐标变换矩阵，将卫星在质心坐标系下的运动映射到地球惯性坐标系下的运动；制定碰撞规避策略；在约束条件下计算卫星的状态方程；选取航天器能量优化指标和哈密顿函数，计算协调因子；计算航天器能量最优条件下的卫星最优速度增量和方向，通过卫星最优速度增量和方向，规避卫星以最优条件下的运动轨迹进行碰撞规避。本发明控制结果精度高，可以实现两颗或多颗卫星的碰撞规避，具有节省燃料、控制操作简单等优点，普遍适用于微小卫星，同时本发明特别适用于近距离的卫星碰撞规避。

Description

一种多约束下航天器能量最优的自主碰撞规避方法
技术领域
本发明涉及航空航天技术领域，具体涉及一种多约束下航天器能量最优的自主碰撞规避方法。
背景技术
自主碰撞规避是防止低轨小卫星碰撞的有效手段和方法，由于商业小卫星迅速发展，小卫星的数量快速增多，增大了碰撞的风险系数，并且测控站的数量有限，如果几颗卫星同时碰撞，仅通过地面同时规避几乎不可能。
目前，现有的卫星自主碰撞规避方法主要有以下几种方式：
《航天器在轨防碰撞自主规避策略》(国防科技大学学报，姚党鼐，王振国，2012,6(34):100-104)针对航天器与空间目标距离较近时的碰撞问题进行了研究，依靠星上测量设备，设计了近距离自主规避策略，但是未考虑测量设备的测量误差，近距离规避碰撞，仍会存在碰撞的风险。
《编队卫星快速解析碰撞预警方法研究》(第一届中国空天安全会议，王有亮，李明涛，郑建华，2015,83-88)提出一种基于特征点的快速解析碰撞预警方法，降低了对绝对测量信息的依赖，通过割线法和高斯方程求解燃料最省的碰撞规避，但是求的解均是近似解，会存在一定偏差。
《编队卫星碰撞规避方法研究》(上海航天，胡敏，曾国强，2010,3:6-10)将碰撞概率密度函数在危险区域积分，获得编队卫星的碰撞概率，当碰撞概率大于安全阈值时，向卫星施加最小脉冲速度修正量，降低了碰撞概率，但是碰撞规避机动后未保持星座构型。
《基于碰撞概率的交会对接最优碰撞规避机动》(宇航学报，王华，李海阳，唐国金，2008,29(1):220-223)提出了一种求解交会对接过程中最优碰撞规避机动方法，在保证碰撞概率降低到安全值的前提下得到最优的避撞机动冲量，采用机动方向和机动大小分布求解的策略计算冲量，适用于相距几十公里以内的碰撞规避，而不适用于近距离的碰撞规避。
发明内容
为了有效解决卫星碰撞问题，本发明提供了一种多约束下航天器能量最优的自主碰撞规避方法。
本发明为解决技术问题所采用的技术方案如下：
本发明的一种多约束下航天器能量最优的自主碰撞规避方法，包括以下步骤：
步骤一、建立卫星运动学模型；
步骤二、建立坐标变换矩阵，将卫星在质心坐标系下的运动映射到地球惯性坐标系下的运动；
步骤三、制定碰撞规避策略；
步骤四、在约束条件下计算卫星的状态方程；
步骤五、选取航天器能量优化指标和哈密顿函数，计算协调因子；
步骤六、计算航天器能量最优条件下的卫星最优速度增量和方向，通过卫星最优速度增量和方向，规避卫星以最优条件下的运动轨迹进行碰撞规避。
作为优选的实施方式，步骤一的具体实现过程如下：
两颗卫星d和g正常在轨飞行，当t时刻，卫星d与卫星g相交于空间某一点时两颗卫星发生碰撞，卫星d相对于卫星g的位置矢量ρ为：
ρ＝r1-r2 (1)
其中，r1为卫星d的地心矢径，单位为km；r2为卫星g的地心矢径，单位为km；
当航天器的相对位置矢量ρ的模ρ小于目标轨道半径时，满足航天器相对运动方程，表示为：
Figure BDA0001848878430000031
其中，航天器的相对位置矢量ρ＝[x,y,z]，单位为m；x、y、z分别表示ρ在x轴，y轴和z轴的位置；[ux,uy,uz]为速度增量，单位为m/s；n为目标轨道平均角速度，单位为°/s；
Figure BDA0001848878430000032
分别表示ρ在x轴，y轴和z轴的速度；
Figure BDA0001848878430000033
Figure BDA0001848878430000034
分别表示ρ在x轴，y轴和z轴的加速度。
作为优选的实施方式，步骤二的具体实现过程如下：
(1)首先定义卫星质心坐标系、轨道坐标系和地球惯性坐标系3种坐标系；
(2)卫星质心坐标系下给出的相对运动方程，初始状态和目标状态条件均是在卫星质心坐标系下给出的，因此需要将卫星质心坐标系转换到轨道坐标系再转到地球惯性坐标系进行控制；
卫星质心坐标系到轨道坐标系的转换矩阵TD为：
Figure BDA0001848878430000035
其中，
Figure BDA0001848878430000036
表示绕z轴的旋转矩阵，
Figure BDA0001848878430000037
表示绕y轴的旋转矩阵，
Figure BDA0001848878430000038
为绕z轴的旋转角度，
Figure BDA0001848878430000039
为绕y轴的旋转角度；
轨道坐标系转到地球惯性坐标系的转换矩阵TT为：
Figure BDA0001848878430000041
其中，[q1 q2 q3 q4]为轨道坐标系相对地球惯性坐标系的姿态四元数；
通过上述坐标变换矩阵，将卫星在质心坐标系下的运动映射到地球惯性坐标系下的运动。
作为优选的实施方式，步骤三的具体实现过程如下：
对于碰撞规避，通过改变相对半长轴控制量来实现星间的安全飞行。假设两颗卫星，一颗为在轨正常运行的目标卫星，另一颗规避卫星不断靠近目标卫星，采用主动控制策略，假设规避卫星从A点到B点的轨迹路径上的任意一点到目标卫星的距离为r(t)，目标卫星的安全球半径为rs，则定义规避卫星在轨迹路径上与安全球的距离为：
Δr(t)＝r(t)-rs (5)
其中，Δr(t)为规避卫星在轨迹路径上与安全球的距离，单位为m；将安全球半径为rs的最小值Δr(t)min作为判断是否发生碰撞的准则，判断准则如下：
Figure BDA0001848878430000042
其中，rs＝rs1+rΔ，rΔ为预留距离，单位为m；rs1表示规定的安全距离，单位为m；
当Δr(t)min＞0时，规避卫星与目标卫星不会发生碰撞；当Δr(t)min＝0时，规避卫星与目标卫星存在碰撞的风险；当Δr(t)min＜0时，规避卫星与目标卫星会发生碰撞。
作为优选的实施方式，步骤四的具体实现过程如下：
(1)控制约束包括：轨道类型约束、有效载荷约束、测控资源约束、控制策略约束和星座能力约束；
(2)计算运动学的状态方程
状态变量定义为：x＝[x y z vx vy vz]，u＝[ux uy uz]T，x、u分别表示状态方程的状态向量和控制速度增量；
运动方程表示卫星的运动状态，卫星的在轨运动以位置和速度来确定，因此将运动方程变换为状态方程；将C-W方程改写为一阶微分方程组的形式：
Figure BDA0001848878430000051
其中，vx、vy、vz分别表示x轴、y轴、z轴的相对速度；则状态方程
Figure BDA0001848878430000052
为：
Figure BDA0001848878430000053
其中，A为系数矩阵；B为常数矩阵；系数矩阵A和常数矩阵B的表达式为：
Figure BDA0001848878430000054
作为优选的实施方式，步骤五的具体实现过程如下：
控制终端状态为：
X(tf)＝Xtf (9)
其中，X(tf)表示tf时状态，Xtf表示状态值，tf表示时间，给定时间tf，选取能量优化指标J为：
Figure BDA0001848878430000061
根据极小值原理，选取哈密顿函数H为：
Figure BDA0001848878430000062
其中，λ1、λ2、λ3、λ4、λ5、λ6均为协调因子；vx、vy、vz分别表示x轴、y轴、z轴的相对速度；n为目标轨道平均角速度；[ux,uy,uz]为速度增量，单位为m/s；
则协态方程为：
Figure BDA0001848878430000063
上式求解为：
Figure BDA0001848878430000064
其中，μi(i＝1,…6)表示待定常数，t表示时间，nt表示t时间内转动的角度。
作为优选的实施方式，步骤六的具体实现过程如下：
根据极小值原理，最优控制
Figure BDA0001848878430000065
和
Figure BDA0001848878430000066
应使哈密顿函数H达到极小，即
Figure BDA0001848878430000071
将式(13)、式(14)带入式(2)中，可得：
Figure BDA0001848878430000072
规避卫星的初始条件和终端条件如下：
Figure BDA0001848878430000073
其中，X(0)表示0时刻的状态，X(tf)表示tf时状态，x0表示x轴的初始位置，y0表示y轴的初始位置，z0表示z轴的初始位置，vx0表示x轴的初始速度，vy0表示y轴的初始速度，vz0表示z轴的初始速度，Xtf表示x轴的终端位置，Ytf表示y轴的终端位置，Ztf表示z轴的终端位置，vxtf表示x轴的终端速度，vytf表示y轴的终端速度，vztf表示z轴的终端速度；
由于z轴的运动是独立的，因此将式(16)带入式(15)中，得到式(17)：
Figure BDA0001848878430000074
其中，S为sin(nt)缩写，C表示cos(nt)缩写，tf表示终端时间；
由于x和y轴存在耦合，通过解微分方程，利用初始条件和终端条件求得一次代数方程组为：
Figure BDA0001848878430000075
Figure BDA0001848878430000081
Figure BDA0001848878430000082
通过式(18)求得μ1、μ2、μ3和μ4；
则卫星最优速度增量为：
Figure BDA0001848878430000083
通过卫星最优速度增量和方向，规避卫星以最优条件下的运动轨迹进行碰撞规避。
本发明的有益效果是：自主碰撞规避是防止低轨小卫星碰撞的有效手段和方法，由于商业小卫星迅速发展，小卫星的数量快速增多，增大了碰撞的风险系数，并且测控站的数量有限，如果几颗卫星同时碰撞，仅通过地面同时规避几乎不可能。
本发明通过建立卫星运动学模型，建立坐标变换矩阵，将卫星在质心坐标系下的运动映射到地球惯性坐标系下的运动，制定碰撞规避策略，在约束条件下计算卫星的状态方程，选取航天器能量优化指标和哈密顿函数，计算协调因子，计算航天器能量最优条件下的卫星最优速度增量和方向，通过卫星最优速度增量和方向，规避卫星以最优条件下的运动轨迹进行碰撞规避。
本发明采用地面观测预警与自主规避相结合的控制策略，测量设备误差对控制策略没有任何影响，同时采用能量最优方法求解出碰撞规避的解析解，控制结果精度高，不仅可以实现两颗卫星的碰撞规避，而且还可实现多颗卫星的碰撞规避，同时具有节省燃料、控制操作简单等优点，普遍适用于微小卫星，同时本发明特别适用于近距离的卫星碰撞规避。
附图说明
图1为卫星碰撞示意图。
图2为J2000地球惯性坐标系示意图。
图3为本发明的多约束下航天器能量最优的自主碰撞规避方法中燃料最优规避控制策略仿真结果。
图4为本发明的多约束下航天器能量最优的自主碰撞规避方法中的燃料最优规避控制策略仿真结果。
图5为本发明的多约束下航天器能量最优的自主碰撞规避方法中的燃料最优规避控制策略仿真结果。
图6为本发明的多约束下航天器能量最优的自主碰撞规避方法中的燃料最优规避控制策略仿真结果。
图7为本发明的多约束下航天器能量最优的自主碰撞规避方法中的卫星在惯性系下卫星轨道图。
图8为规避卫星控前、控后轨道半长轴图。
图9为碰撞规避图。
图10为本发明的多约束下航天器能量最优的自主碰撞规避方法中的卫星碰撞示意图。
图11为两颗卫星的距离变化曲线。
图12为两颗卫星的距离最近时通过本发明的碰撞规避控制后两颗卫星之间的最小距离。
具体实施方式
下面结合附图和实施例对本发明作进一步详细说明。
本发明的一种多约束下航天器能量最优的自主碰撞规避方法，首先建立运动学模型和坐标变换矩阵，制定碰撞规避策略，其次最优化碰撞规避控制算法的设计。具体包括以下步骤：
步骤一、运动学建模、坐标变换矩阵建立以及碰撞规避策略制定
S11：运动学建模
如图1所示，两颗卫星d和g正常在轨飞行，当t时刻，卫星d与卫星g相交于空间某一点，两颗卫星即发生碰撞。卫星d相对于卫星g的位置矢量ρ为：
ρ＝r1-r2 (1)
其中，r1为卫星d的地心矢径，单位为km；r2为卫星g的地心矢径，单位为km。
当航天器的相对位置矢量ρ的模ρ(ρ＜100km)小于目标轨道半径时，满足航天器相对运动方程(C-W方程)，表示为：
Figure BDA0001848878430000101
其中，航天器的相对位置矢量ρ＝[x,y,z]，单位为m；x、y、z分别表示ρ在x轴，y轴和z轴的位置，[ux,uy,uz]为速度增量，单位为m/s；n为目标轨道平均角速度，单位为°/s；
Figure BDA0001848878430000102
分别表示ρ在x轴，y轴和z轴的速度；
Figure BDA0001848878430000103
Figure BDA0001848878430000104
分别表示ρ在x轴，y轴和z轴的加速度。
S12：坐标变换矩阵建立
首先定义卫星质心坐标系、轨道坐标系和地球惯性坐标系等3种坐标系。
卫星质心坐标系：
卫星质心坐标系是指固连于卫星的Od-xdydzd直角坐标系，坐标原点Od在卫星的质心处，定义xd轴与目标航天器的地心矢径r2重合，yd轴在轨道面内垂直于xd轴，运动方向为正，zd轴由右手法则确定，即轨道动量矩的方向。
轨道坐标系：
轨道坐标系是指Og-xgygzg直角坐标系，其坐标原点为卫星质心Og，+xg轴指向卫星的飞行方向，+zg轴指向地球中心，+yg轴与+xg、+zg轴成右手直角坐标系，此坐标系为对地定向基准。
地球惯性坐标系：
如图2所示，地球惯性坐标系(地球第一赤道惯性坐标系)是指Ce-XeiYeiZei直角坐标系，坐标系原点为地心Ce，Xei轴指向平春分点(2000年1月1日12时)，Zei轴指向平北极(2000年1月1日12时，JD＝2451545.0)，Yei轴和Xei、Zei轴构成右手直角坐标系，图2中示出天球1、黄道平面2、赤道平面3、轨道平面4、轨道倾角5、升交点赤经6、升交点辐角7和黄赤交角8的位置，地球惯性坐标系也为J2000地球惯性坐标系，也为姿态测量基准和轨道测量基准。
卫星质心坐标系下给出的相对运动方程，初始状态和目标状态条件均是在卫星质心坐标系下给出的，因此需要将卫星质心坐标系转换到轨道坐标系再转到地球惯性坐标系进行控制。卫星质心坐标系到轨道坐标系的转换矩阵TD为：
Figure BDA0001848878430000111
其中，
Figure BDA0001848878430000112
表示绕z轴的旋转矩阵，
Figure BDA0001848878430000113
表示绕y轴的旋转矩阵，
Figure BDA0001848878430000114
为绕z轴的旋转角度，
Figure BDA0001848878430000121
为绕y轴的旋转角度；
轨道坐标系转到地球惯性坐标系的转换矩阵TT为：
Figure BDA0001848878430000122
其中，[q1 q2 q3 q4]为轨道坐标系相对地球惯性坐标系的姿态四元数。
通过以上的坐标变换矩阵，可以将卫星质心坐标系下的运动映射到地球惯性坐标系下的运动，实现碰撞规避控制。
S13：碰撞规避策略的制定
对于碰撞规避，通过改变相对半长轴控制量来实现星间的安全飞行。假设两颗卫星，一颗为在轨正常运行的目标卫星，另一颗规避卫星不断靠近目标卫星，采用主动控制策略，假设规避卫星从A点到B点的轨迹路径上的任意一点到目标卫星的距离为r(t)，目标卫星的安全球半径为rs，则定义规避卫星在轨迹路径上与安全球的距离为：
Δr(t)＝r(t)-rs (5)
其中，Δr(t)为规避卫星在路径上与安全球的距离，单位为m。安全球半径为rs的最小值Δr(t)min可以作为用来判断是否发生碰撞的准则，判断准则具体如下：
Figure BDA0001848878430000123
当Δr(t)min＞0时，规避卫星与目标卫星不会发生碰撞，无碰撞风险；当Δr(t)min＝0时，规避卫星与目标卫星存在碰撞的风险；当Δr(t)min＜0时，规避卫星与目标卫星会发生碰撞的风险。其中，rs＝rs1+rΔ，rΔ为预留距离，即准备启控时间，单位为m，rs1表示规定的安全距离，单位为m。
步骤二、最优化碰撞规避控制算法设计
S21：约束条件
控制约束包括：轨道类型约束、有效载荷约束、测控资源约束、控制策略约束和星座能力约束，这些因素直接影响碰撞规避控制策略。
1)轨道类型约束：要求太阳同步轨道，当碰撞规避控制完成后，轨道仍为太阳同步轨道或近太阳同步轨道；
2)有效载荷约束：遥感卫星的载荷与高度密切相关，当轨道高度越高，载荷的分辨率将会下降，不满足用户需求；当轨道高度越低，卫星的轨道衰减越快，使用寿命大大减少；
3)测控资源约束：首先确定可利用的测控站，尽量在测控弧段内进行轨道调整，当碰撞位置和时间确定后，用较近的测控站进行控制指令发送；
4)控制策略约束：遥感小卫星质量轻，携带的燃料少，因此采用最优燃料控制。当规避卫星进入目标卫星的安全范围后，给出规避卫星控制标志，在测控站弧段内制定规避卫星规避策略，设置控制后的目标轨道，控制时刻和控制时长；
5)星座能力约束：对于商业遥感星座，重访和覆盖是客户最为关心的能力，相位分布情况直接影响重访能力，因此，相位分布的保持尤为重要，若同一个轨道面上的卫星高度相差较大，几天时间就会超出允许相位差，因此控后尽量回到原来的高度或百米级的高度差。
S22：计算运动学的状态方程表达
状态变量定义为：x＝[x y z vx vy vz]，u＝[ux uy uz]T，x、u分别表示状态方程的状态向量和控制速度增量。
运动方程表述了卫星的运动状态，卫星的在轨运动主要以位置和速度来确定卫星的运动，因此将运动方程变换为状态方程，方便求解。
将C-W方程改写为一阶微分方程组的形式如下：
Figure BDA0001848878430000141
其中，vx、vy、vz分别表示x轴、y轴、z轴的相对速度；则状态方程
Figure BDA0001848878430000146
可表达为：
Figure BDA0001848878430000142
其中，A为系数矩阵；B为常数矩阵；系数矩阵A和常数矩阵B的表达式为：
Figure BDA0001848878430000143
S23：选取能量优化指标和哈密顿函数，计算协调因子
控制终端状态为：
X(tf)＝Xtf (9)
其中，X(tf)表示tf时状态，Xtf表示状态值，tf表示时间，给定时间tf，则为固定时间和固定端点的最优控制，选取能量优化指标J为：
Figure BDA0001848878430000144
依据极小值原理，选取哈密顿函数H为：
Figure BDA0001848878430000145
其中，λ1、λ2、λ3、λ4、λ5、λ6均为协调因子；vx、vy、vz分别表示x轴、y轴、z轴的相对速度；n为目标轨道平均角速度；[ux,uy,uz]为速度增量，单位为m/s；
则协态方程为：
Figure BDA0001848878430000151
上式求解为：
Figure BDA0001848878430000152
其中，μi(i＝1,…6)表示待定常数，t表示时间，nt表示t时间内转动的角度。
步骤三、计算最优控制率
根据极小值原理，最优控制
Figure BDA0001848878430000153
和
Figure BDA0001848878430000154
应使哈密顿函数H达到极小，即
Figure BDA0001848878430000155
将式(13)、式(14)带入式(2)中，可得：
Figure BDA0001848878430000156
规避卫星的初始条件和终端条件如下：
Figure BDA0001848878430000161
其中，X(0)表示0时刻的状态，X(tf)表示tf时状态，x0表示x轴的初始位置，y0表示y轴的初始位置，z0表示z轴的初始位置，vx0表示x轴的初始速度，vy0表示y轴的初始速度，vz0表示z轴的初始速度，Xtf表示x轴的终端位置，Ytf表示y轴的终端位置，Ztf表示z轴的终端位置，vxtf表示x轴的终端速度，vytf表示y轴的终端速度，vztf表示z轴的终端速度。
由于z轴的运动是独立的，因此将式(16)带入式(15)中，可得：
Figure BDA0001848878430000162
其中，S为sin(nt)缩写，C表示cos(nt)缩写，tf表示终端时间。
由于x和y轴存在耦合，可通过解微分方程，利用初始条件和终端条件求得，一次代数方程组为：
Figure BDA0001848878430000163
Figure BDA0001848878430000164
Figure BDA0001848878430000165
则可求得μ1、μ2、μ3和μ4。
卫星最优速度增量为：
Figure BDA0001848878430000171
步骤四、碰撞规避控制策略
设置好预置时间t、初始状态X(0)和终端状态X(tf)、以及控制时长tf，选择约束条件，当预测到在轨的两颗卫星某一时刻发生碰撞时，规避卫星则到达预置时间时自主开始规避机动，通过最优化碰撞规避控制算法实时给出卫星最优速度增量和方向，规避卫星会以最省能量的运动轨迹进行碰撞规避。
具体实施方式一
采用二体递推模型给出卫星飞行轨迹，最小安全距离设置为300m，小于300m时即发出预警并进行碰撞规避控制，通过两颗卫星轨道碰撞仿真验证所设计的最优控制规避算法的有效性和可行性。
轨道类型约束：要求太阳同步轨道，当规避控制完成后，轨道仍为太阳同步轨道或近太阳同步轨道，若是二体模型轨道，则不做要求；
有效载荷约束：有效载荷分辨率为1m(@500km)，允许分辨率上限为1.002m，则允许调整的轨道高度不高于501km；
控制策略约束：当规避卫星进入目标卫星的安全范围后，以最优最快的路径规避目标卫星的轨道；
星座能力约束：规避卫星轨道有两颗卫星，以180°等相位分布，允许相位偏置7°，当调整50m时，相位差飘出0.059°/天，则118天飘出应用需求允许范围，因此允许调整的轨道高度范围±50m以内；
测控资源约束：可用的测控站位长春站、三亚站和喀什站，启控时刻尽量选择在测控弧段内。
轨道初始值和控制参数如表1和表2所示。
表1卫星初始轨道参数
卫星 目标卫星 规避卫星 共面星
历元时间 2018.7.1022:25:00 2018.7.1022:25:00 2018.7.1022:25:00
X -2084274.526m -2278408.122m 2278408.122m
Y 6465866.521m 6510485.857m -6510485.857m
Z -1719553.117m 467688.828m -467688.828m
V<sub>x</sub> 1151.289m/s 762.625m/s -762.625m/s
V<sub>y</sub> -1973.692m/s 806.49m/s -806.49m/s
V<sub>z</sub> -7197.116m/s -7511.555m/s 7511.555m/s
表2控制策略参数
控制参数 数值
控制时长 300s
启控时刻 2018.7.1022:59:45
初始位置 [0 0 0]<sup>T</sup>m
初始速度 [0 0 0]<sup>T</sup>m/s
目标位置 [300 300 300]<sup>T</sup>m
目标速度 [0 0 0.15]<sup>T</sup>m/s
采用matlab进行数值仿真，仿真时间12小时，设置目标卫星安全距离为300m，两颗卫星2018.7.1103:59:45开始发生碰撞，图3给出了在燃料最优规避控制策略下规避卫星三维坐标系下的运行轨迹，图4给出了X轴速度增量，图5给出了Y轴速度增量，图6给出了Z轴速度增量。
从图3-6可知，规避卫星在质心系下从[0 0 0]T运动到[300 300 300]T，误差在3m以内，将卫星质心系下的位置和速度通过坐标转换到惯性系下的位置和速度。
当两颗卫星距离小于300m时，卫星将要碰撞，提前6小时在测控弧段内开始进行轨道碰撞控制，其中，图7中实线表示规避卫星运行轨道，点画线表示目标卫星运行轨道，虚线表示规避卫星控后运行轨道，图8中实线表示规避卫星控前轨道局部放大曲线，虚线表示规避卫星控后运行轨道局部放大曲线，图9中实线表示规避卫星控后半长轴曲线，虚线表示规避卫星控前半长轴曲线，图10中实线表示半长轴差。从图7-10可知，轨道半长轴控制量最大不超过28.3m，满足星座构形约束和有效载荷约束。
从图11可知，当2018.7.1103:59:45时，两颗卫星的距离变化，从图12可看出两颗卫星的距离最近时小于300m，碰撞规避控制后的两颗卫星最小距离为8000m，以最优燃料控制规避了两星碰撞。
综上，从仿真效果可知，满足各类约束的情况下，采用设计的自主碰撞规避算法很好的实现了卫星的自主碰撞规避，证明了最优碰撞规避算法的有效可行。
以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (1)

1.一种多约束下航天器能量最优的自主碰撞规避方法，其特征在于，包括以下步骤：
步骤一、建立卫星运动学模型；
两颗卫星d和g正常在轨飞行，当t时刻，卫星d与卫星g相交于空间某一点时两颗卫星发生碰撞，卫星d相对于卫星g的位置矢量ρ为：
ρ＝r1-r2 (1)
其中，r1为卫星d的地心矢径，单位为km；r2为卫星g的地心矢径，单位为km；
当航天器的相对位置矢量ρ的模ρ小于目标轨道半径时，满足航天器相对运动方程，表示为：
Figure FDA0002567356380000011
其中，航天器的相对位置矢量ρ＝[x,y,z]，单位为m；x、y、z分别表示ρ在x轴，y轴和z轴的位置；[ux,uy,uz]为速度增量，单位为m/s；n为目标轨道平均角速度，单位为°/s；
Figure FDA0002567356380000012
分别表示ρ在x轴，y轴和z轴的速度；
Figure FDA0002567356380000013
Figure FDA0002567356380000014
分别表示ρ在x轴，y轴和z轴的加速度；
步骤二、建立坐标变换矩阵，将卫星在质心坐标系下的运动映射到地球惯性坐标系下的运动；
(1)首先定义卫星质心坐标系、轨道坐标系和地球惯性坐标系3种坐标系；
(2)卫星质心坐标系下给出的相对运动方程，初始状态和目标状态条件均是在卫星质心坐标系下给出的，因此需要将卫星质心坐标系转换到轨道坐标系再转到地球惯性坐标系进行控制；
卫星质心坐标系到轨道坐标系的转换矩阵TD为：
Figure FDA0002567356380000021
其中，
Figure FDA0002567356380000022
表示绕z轴的旋转矩阵，
Figure FDA0002567356380000023
表示绕y轴的旋转矩阵，
Figure FDA0002567356380000024
为绕z轴的旋转角度，
Figure FDA0002567356380000025
为绕y轴的旋转角度；
轨道坐标系转到地球惯性坐标系的转换矩阵TT为：
Figure FDA0002567356380000026
其中，[q1 q2 q3 q4]为轨道坐标系相对地球惯性坐标系的姿态四元数；
通过上述坐标变换矩阵，将卫星在质心坐标系下的运动映射到地球惯性坐标系下的运动；
步骤三、制定碰撞规避策略；
对于碰撞规避，通过改变相对半长轴控制量来实现星间的安全飞行；假设两颗卫星，一颗为在轨正常运行的目标卫星，另一颗规避卫星不断靠近目标卫星，采用主动控制策略，假设规避卫星从A点到B点的轨迹路径上的任意一点到目标卫星的距离为r(t)，目标卫星的安全球半径为rs，则定义规避卫星在轨迹路径上与安全球的距离为：
Δr(t)＝r(t)-rs (5)
其中，Δr(t)为规避卫星在轨迹路径上与安全球的距离，单位为m；将安全球半径为rs的最小值Δr(t)min作为判断是否发生碰撞的准则，判断准则如下：
Figure FDA0002567356380000031
其中，rs＝rs1+rΔ，rΔ为预留距离，单位为m；rs1表示规定的安全距离，单位为m；
当Δr(t)min＞0时，规避卫星与目标卫星不会发生碰撞；当Δr(t)min＝0时，规避卫星与目标卫星存在碰撞的风险；当Δr(t)min＜0时，规避卫星与目标卫星会发生碰撞；
步骤四、在约束条件下计算卫星的状态方程；
(1)控制约束包括：轨道类型约束、有效载荷约束、测控资源约束、控制策略约束和星座能力约束；
(2)计算运动学的状态方程
状态变量定义为：x＝[x y z vx vy vz]，u＝[ux uy uz]T，x、u分别表示状态方程的状态向量和控制速度增量；
运动方程表示卫星的运动状态，卫星的在轨运动以位置和速度来确定，因此将运动方程变换为状态方程；将C-W方程改写为一阶微分方程组的形式：
Figure FDA0002567356380000032
其中，vx、vy、vz分别表示x轴、y轴、z轴的相对速度；则状态方程
Figure FDA0002567356380000033
为：
Figure FDA0002567356380000034
其中，A为系数矩阵；B为常数矩阵；系数矩阵A和常数矩阵B的表达式为：
Figure FDA0002567356380000041
步骤五、选取航天器能量优化指标和哈密顿函数，计算协调因子；
控制终端状态为：
X(tf)＝Xtf (9)
其中，X(tf)表示tf时状态，Xtf表示状态值，tf表示时间，给定时间tf，选取能量优化指标J为：
Figure FDA0002567356380000042
根据极小值原理，选取哈密顿函数H为：
Figure FDA0002567356380000043
其中，λ1、λ2、λ3、λ4、λ5、λ6均为协调因子；vx、vy、vz分别表示x轴、y轴、z轴的相对速度；n为目标轨道平均角速度；[ux,uy,uz]为速度增量，单位为m/s；
则协态方程为：
Figure FDA0002567356380000044
上式求解为：
Figure FDA0002567356380000051
其中，μi(i＝1,…6)表示待定常数，t表示时间，nt表示t时间内转动的角度；
步骤六、计算航天器能量最优条件下的卫星最优速度增量和方向，通过卫星最优速度增量和方向，规避卫星以最优条件下的运动轨迹进行碰撞规避；
根据极小值原理，最优控制
Figure FDA0002567356380000052
和
Figure FDA0002567356380000053
应使哈密顿函数H达到极小，即
Figure FDA0002567356380000054
将式(13)、式(14)带入式(2)中，可得：
Figure FDA0002567356380000055
规避卫星的初始条件和终端条件如下：
Figure FDA0002567356380000056
其中，X(0)表示0时刻的状态，X(tf)表示tf时状态，x0表示x轴的初始位置，y0表示y轴的初始位置，z0表示z轴的初始位置，vx0表示x轴的初始速度，vy0表示y轴的初始速度，vz0表示z轴的初始速度，Xtf表示x轴的终端位置，Ytf表示y轴的终端位置，Ztf表示z轴的终端位置，vxtf表示x轴的终端速度，vytf表示y轴的终端速度，vztf表示z轴的终端速度；
由于z轴的运动是独立的，因此将式(16)带入式(15)中，得到式(17)：
Figure FDA0002567356380000061
其中，S为sin(nt)缩写，C表示cos(nt)缩写，tf表示终端时间；
由于x和y轴存在耦合，通过解微分方程，利用初始条件和终端条件求得一次代数方程组为：
Figure FDA0002567356380000062
Figure FDA0002567356380000063
Figure FDA0002567356380000064
通过式(18)求得μ1、μ2、μ3和μ4；
则卫星最优速度增量为：
Figure FDA0002567356380000065
CN201811285216.1A 2018-10-31 2018-10-31 一种多约束下航天器能量最优的自主碰撞规避方法 Expired - Fee Related CN109062243B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201811285216.1A CN109062243B (zh) 2018-10-31 2018-10-31 一种多约束下航天器能量最优的自主碰撞规避方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201811285216.1A CN109062243B (zh) 2018-10-31 2018-10-31 一种多约束下航天器能量最优的自主碰撞规避方法

Publications (2)

Publication Number Publication Date
CN109062243A CN109062243A (zh) 2018-12-21
CN109062243B true CN109062243B (zh) 2021-04-20

Family

ID=64789048

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201811285216.1A Expired - Fee Related CN109062243B (zh) 2018-10-31 2018-10-31 一种多约束下航天器能量最优的自主碰撞规避方法

Country Status (1)

Country Link
CN (1) CN109062243B (zh)

Families Citing this family (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109765919B (zh) * 2019-02-27 2020-04-14 中国人民解放军军事科学院国防科技创新研究院 基于等碰撞概率线法的航天器近距离安全操作控制方法
CN110109476A (zh) * 2019-05-06 2019-08-09 中国人民解放军军事科学院国防科技创新研究院 基于碰撞概率的异构无人飞行器自主碰撞预警方法
CN110377046B (zh) * 2019-09-04 2022-03-11 哈尔滨工业大学 无人机在舰船上着陆的控制方法
CN110595486B (zh) * 2019-09-05 2021-04-23 上海航天控制技术研究所 基于双星在轨遥测数据的高精度半长轴偏差计算方法
CN111221254B (zh) * 2020-03-26 2021-08-10 大连理工大学 一种基于割线法确定终端时间的航天器姿态调整能量-时间混合最优控制方法
CN114326774B (zh) * 2022-03-14 2022-05-24 北京航天驭星科技有限公司 航天器碰撞规避策略生成的方法及系统
CN115416877B (zh) * 2022-09-06 2023-05-02 宁波天巡科技有限公司 一种基于小推力的卫星碰撞高度规避方法及智能化系统
CN115373425B (zh) * 2022-10-26 2023-01-20 中国西安卫星测控中心 一种基于编队绕飞模型的轨控碰撞安全性分析方法
CN116331518B (zh) * 2023-01-13 2023-09-26 天津大学 一种基于安全自适应动态规划的星群智能编队避碰控制方法
CN116166049B (zh) * 2023-04-25 2023-07-07 中国西安卫星测控中心 不稳定多星串联编队系统星间距离保持控制方法
CN117332684B (zh) * 2023-09-25 2024-04-26 同济大学 一种基于强化学习的多航天器追逃博弈下的最优捕获方法

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105242680B (zh) * 2015-10-28 2018-07-06 哈尔滨工业大学 一种相对非合作目标的航天器相对轨道有限时间抗饱和控制方法
CN105549606B (zh) * 2015-12-21 2018-04-13 北京理工大学 针对失效卫星的超近距离最优防撞接近方法
CN106707751B (zh) * 2016-12-20 2019-05-17 哈尔滨工业大学 航天器终端接近的有限时间饱和避碰控制方法

Also Published As

Publication number Publication date
CN109062243A (zh) 2018-12-21

Similar Documents

Publication Publication Date Title
CN109062243B (zh) 一种多约束下航天器能量最优的自主碰撞规避方法
CN107380485B (zh) 一种微小卫星大面阵广域多模凝视成像控制方法
Vavrina et al. Safe rendezvous trajectory design for the restore-l mission
CN110733667A (zh) 地月平动点轨道间转移设计方法
CN109774974B (zh) 一种用于空间碎片接近的轨道设计方法
CN110450991B (zh) 微纳卫星集群捕获空间非合作目标的方法
CN109190158B (zh) 一种考虑非合作目标禁飞区约束的最优轨道设计方法
CN108663052B (zh) 一种星上自主空间非合作目标相对导航相机指向控制方法
Wolf et al. Toward improved landing precision on Mars
CN112498746B (zh) 一种自主规划卫星沿经线推扫时间及姿态的方法
Somov et al. Nonlinear digital and pulse-width control at approaching a space robot with a geostationary satellite.
Martin et al. Astrone–GNC for Enhanced Surface Mobility on Small Solar System Bodies
Kos et al. Altair descent and ascent reference trajectory design and initial dispersion analyses
Shen et al. Point return orbit design and characteristics analysis for manned lunar mission
Xie et al. Guidance, navigation, and control for spacecraft rendezvous and docking: theory and methods
Zhang et al. Hybrid ground-space target visit problem with a coplanar impulse
Wood The Evolution of Deep Space Navigation: 2004–2006
CN112329202B (zh) 一种火星车对环绕器天线指向算法的优化实现方法
O’Shaughnessy et al. Fire Sail: MESSENGER’s use of solar radiation pressure for accurate Mercury flybys
Mimasu et al. Hayabusa2 operation for MASCOT delivery to Ryugu surface
Vaughan et al. Navigation flight operations for Mars Pathfinder
Treder Space station GN&C overview for payloads
Chen et al. Near-Earth Orbit Satellite Collision Probability Estimation and Collision Avoidance
Xie et al. Autonomous guidance, navigation, and control of spacecraft
Somov et al. Space Robot Motion Control during Rendezvous and Visual Inspecting a Geostationary Satellite State

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
CB03 Change of inventor or designer information

Inventor after: Wang Guogang

Inventor after: Liu Keping

Inventor after: Dai Lu

Inventor after: Fan Lindong

Inventor after: Li Feng

Inventor after: Hu Jianlong

Inventor after: Meng Deli

Inventor before: Wang Guogang

Inventor before: Dai Lu

Inventor before: Xu Kai

Inventor before: Fan Lindong

Inventor before: Li Feng

Inventor before: Hu Jianlong

Inventor before: Meng Deli

CB03 Change of inventor or designer information
TA01 Transfer of patent application right

Effective date of registration: 20210317

Address after: 130000 No. 17, Yanan Road, Changchun, Jilin

Applicant after: CHANGCHUN University OF TECHNOLOGY

Address before: No.1759 Mingxi Road, Gaoxin North District, Changchun City, Jilin Province

Applicant before: CHANG GUANG SATELLITE TECHNOLOGY Co.,Ltd.

TA01 Transfer of patent application right
GR01 Patent grant
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20210420

Termination date: 20211031

CF01 Termination of patent right due to non-payment of annual fee