基于多等碰撞概率线法的凸多边形航天器安全控制方法

摘要

本发明公开了一种基于多等碰撞概率线法的凸多边形航天器安全控制方法。该方法用于实现目标航天器和跟踪航天器的安全接近控制，包括步骤：建立目标航天器的轨道坐标系；建立跟踪航天器与目标航天器的相对动力学方程；进行相对动力学方程的不确定性传播分析以推导不确定性的协方差矩阵；确定目标航天器中产生作用的作用模块；利用等碰撞概率线法建立多等碰撞概率线法以计算基于作用模块产生的施加于跟踪航天器的避障控制力；确定多等碰撞概率线法的参数的取值范围；确定最终作用于跟踪航天器上的总控制力。本发明的基于多等碰撞概率线法的凸多边形航天器安全控制方法，能够解决目标航天器为凸多边形情况下的跟踪航天器的安全接近控制问题。

权利要求

1.一种基于多等碰撞概率线法的凸多边形航天器安全控制方法，其特征在于，所述方法用于实现目标航天器和跟踪航天器的安全接近控制，所述方法包括如下步骤：

1)建立目标航天器的轨道坐标系；

2)在所述目标航天器的轨道坐标系下，建立跟踪航天器与所述目标航天器的相对动力学方程；

3)进行所述相对动力学方程的不确定性传播分析，采用线性协方差方法推导不确定性的协方差矩阵；

4)将所述目标航天器的几何外形分解为一个最大内包络圆和若干个任意形状部分，分别设定所述最大内包络圆和每个所述任意形状部分的外包络圆为一个作用模块，根据所述跟踪航天器与所述目标航天器的相对位置，确定多个所述作用模块中产生作用的作用模块；

5)利用等碰撞概率线法建立多等碰撞概率线法以计算基于所述作用模块产生的施加于所述跟踪航天器的避障控制力；

6)对所述避障控制力进行分析，确定所述多等碰撞概率线法的参数的取值范围，以使所述跟踪航天器与所述目标航天器不会发生碰撞；

7)确定最终作用于所述跟踪航天器上的总控制力。

2.根据权利要求1所述的基于多等碰撞概率线法的凸多边形航天器安全控制方法，其特征在于，建立目标航天器的轨道坐标系包括：采用o-xyz表示所述目标航天器的轨道坐标系，以所述目标航天器的质心为坐标原点o，x轴由地球地心指向所述目标航天器的质心，z轴指向所述目标航天器的轨道面法线方向，y轴垂直于所述x轴和所述z轴所组成的平面，所述y轴与所述x轴、所述z轴构成右手直角坐标系。

3.根据权利要求2所述的基于多等碰撞概率线法的凸多边形航天器安全控制方法，其特征在于，所述相对动力学方程如式1所示；

其中，r1-t＝[xt,yt]T和v1-t＝[vx,vy]T表示在所述目标航天器的轨道坐标系下所述跟踪航天器在t时刻的相对位置和相对速度， u1-t＝[ux,uy]T表示所述跟踪航天器在t时刻的控制输入，μ表示地球引力常数， a和n分别表示所述目标航天器的轨道半长轴和平均角速度；

假设状态矢量 式1改写为：

其中，A表示状态转移矩阵，B表示控制矩阵，

4.根据权利要求3所述的基于多等碰撞概率线法的凸多边形航天器安全控制方法，其特征在于，采用线性协方差方法推导不确定性的协方差矩阵，包括：利用式5计算所述状态矢量X的解析解；

其中，t0表示初始时间，X0表示与t0对应的状态矢量初始值，Φ(t,t0)表示状态矩阵且满足式6；

Φrr(t,t0)、Φrv(t,t0)、Φvr(t,t0)和Φvv(t,t0)满足式7；

利用式8和式9计算所述状态矢量X的不确定性协方差矩阵；

δX＝X(t)-Ε(X(t))  (8)

其中，N表示所述跟踪航天器所施加的控制脉冲的数量， 和 分别表示在所述目标航天器的轨道坐标系下的初始导航不确定性和控制不确定性的协方差矩阵，CδX表示所述状态矢量的不确定性协方差矩阵，Ε(X(t))表示所述状态矢量的平均值，δX为状态矢量偏差，E(δX)为所述状态矢量偏差的平均值。

5.根据权利要求4所述的基于多等碰撞概率线法的凸多边形航天器安全控制方法，其特征在于，基于所述目标航天器的轨道坐标系，根据所述跟踪航天器所处位置，将所述目标航天器划分为如式10所示的5个作用模块；

其中，[x,y]表示所述跟踪航天器在所述目标航天器的轨道坐标系下的相对位置坐标，Part5表示所述目标航天器中最大内包络圆对应部分，Part1表示位于所述目标航天器的轨道坐标系的xoy第一象限内的目标航天器部分扣除最大内包络圆所属部分后的对应部分，Part2表示位于所述目标航天器的轨道坐标系的xoy第二象限内的目标航天器部分扣除最大内包络圆所属部分后的对应部分，Part3表示位于所述目标航天器的轨道坐标系的xoy第三象限内的目标航天器部分扣除最大内包络圆所属部分后的对应部分，Part4表示位于所述目标航天器的轨道坐标系的xoy第四象限内的目标航天器部分扣除最大内包络圆所属部分后的对应部分。

6.根据权利要求5所述的基于多等碰撞概率线法的凸多边形航天器安全控制方法，其特征在于，根据所述跟踪航天器的最终时刻的相对位置的不同，对施加于所述跟踪航天器的避障控制力进行分类分析。

7.根据权利要求6所述的基于多等碰撞概率线法的凸多边形航天器安全控制方法，其T特征在于，设定所述跟踪航天器的最终时刻的相对位置为rep＝[xep yep]；

当xep＝0 or yep＝0时，施加于所述跟踪航天器的避障控制力如式48所示；

其中，F50-repel+Fh0-repel,h＝1,2,3,4表示当所述跟踪航天器位于h＝1,2,3或4所对应的区域时，施加于所述跟踪航天器的避障控制力为F50-repel+Fh0-repel，F50-repel,h＝5表示当所述跟踪航天器位于h＝5所对应的区域时，施加于所述跟踪航天器的避障控制力为F50-repel。

8.根据权利要求7所述的基于多等碰撞概率线法的凸多边形航天器安全控制方法，其特征在于，当xep≠0 and yep≠0时，根据相对位置rep＝[xep yep]T建立两个坐标系统O1-ix1-iy1-i(i＝1,2)，并根据相对位置rep＝[xep yep]T的不同分为如式49所示的四个情况进行避障控制力分析；

施加于所述跟踪航天器的避障控制力如式59所示；

其中，F50-repel+Fh1-repel+Fh2-repel,h＝h1表示当所述跟踪航天器位于h＝h1所对应的区域时，施加于所述跟踪航天器的避障控制力为F50-repel+Fh1-repel+Fh2-repel，F50-repel+Fh0-repel,h＝

1,2,3,4&h≠h1表示当所述跟踪航天器位于h＝1,2,3,4&h≠h1所对应的区域时，施加于所述跟踪航天器的避障控制力为F50-repel+Fh0-repel，F50-repel,h＝5表示当所述跟踪航天器位于h＝

5所对应的区域时，施加于所述跟踪航天器的避障控制力为F50-repel。

9.根据权利要求8所述的基于多等碰撞概率线法的凸多边形航天器安全控制方法，其特征在于，根据式70所示条件，确定所述多等碰撞概率线法的参数的取值范围，以使所述跟踪航天器与所述目标航天器不会发生碰撞；

其中，R50、Ri0(i＝1,2,3,4)和Rij(i＝1,2,3,4,j＝0,1,2)表示相对应的作用模块的最小外包络圆的半径，ζ50、ζi0(i＝1,2,3,4)和ζij(i＝1,2,3,4,j＝0,1,2)表示相对应的作用模块的最小外包络圆的圆心到跟踪航天器的最短距离。

10.根据权利要求9所述的基于多等碰撞概率线法的凸多边形航天器安全控制方法，其特征在于，作用于所述跟踪航天器的总控制力利用式81计算；

其中，uTotal表示作用于所述跟踪航天器的总控制力， 表示最优控制，Frepel表示所述避障控制力，m表示所述跟踪航天器的质量。

说明书

基于多等碰撞概率线法的凸多边形航天器安全控制方法

技术领域

[0001] 本发明涉及航天器运动控制技术领域，具体涉及一种基于多等碰撞概率线法的凸多边形航天器安全控制方法。

背景技术

[0002] 近年来，航天器在轨失效事件日渐增多，为了降低在轨失效事件发生的概率，延长航天器工作年限，提高工作性能，越来越多的在轨服务被应用于航天器，而航天器近距离操作作为一项支撑在轨服务的基本技术，航天器近距离操作需满足严格的安全性要求。

[0003] 传统上，采用人工势函数(Artificial Potential Function，APF)法来保证航天器近距离操作的安全性，该方法具有计算简单和理论易于证明等特点；但由于人工势函数未考虑不确定性的影响，因而包括导航和控制等在内的各类不确定性会对航天器近距离操作的安全性产生影响。针对航天器安全接近过程中不确定性的影响，目前采用一种等碰撞概率线法来保证航天器近距离操作的安全性，该等碰撞概率线法由于没有包含超越函数，因此相比于传统的等碰撞概率函数，能在保证有效精度的同时大大降低计算量；且能够通过理论证明等碰撞概率线法的有效性。然而，现有的等碰撞概率线法中，两个航天器的几何外形被简化为球形或椭球形，但是在工程实践中，大多数航天器的几何外形并不是简单的球形或椭球形，在近距离操作中，航天器的不同几何外形对安全控制具有不同的影响，导致现有的等碰撞概率线法在工程实际应用中具有一定的局限性。

发明内容

[0004] 为解决上述现有技术中存在的技术问题，本发明提供一种基于多等碰撞概率线法的凸多边形航天器安全控制方法。

[0005] 为此，本发明公开了一种基于多等碰撞概率线法的凸多边形航天器安全控制方法。所述方法用于实现目标航天器和跟踪航天器的安全接近控制，所述方法包括如下步骤：

[0006] 1)建立目标航天器的轨道坐标系；

[0007] 2)在所述目标航天器的轨道坐标系下，建立跟踪航天器与所述目标航天器的相对动力学方程；

[0008] 3)进行所述相对动力学方程的不确定性传播分析，采用线性协方差方法推导不确定性的协方差矩阵；

[0009] 4)将所述目标航天器的几何外形分解为一个最大内包络圆和若干个任意形状部分，分别设定所述最大内包络圆和每个所述任意形状部分的外包络圆为一个作用模块，根据所述跟踪航天器与所述目标航天器的相对位置，确定多个所述作用模块中产生作用的作用模块；

[0010] 5)利用等碰撞概率线法建立多等碰撞概率线法以计算基于所述作用模块产生的施加于所述跟踪航天器的避障控制力；

[0011] 6)对所述避障控制力进行分析，确定所述多等碰撞概率线法的参数的取值范围，以使所述跟踪航天器与所述目标航天器不会发生碰撞；

[0012] 7)确定最终作用于所述跟踪航天器上的总控制力。

[0013] 进一步地，在所述基于多等碰撞概率线法的凸多边形航天器安全控制方法中，建立目标航天器的轨道坐标系包括：

[0014] 采用o-xyz表示所述目标航天器的轨道坐标系，以所述目标航天器的质心为坐标原点o，x轴由地球地心指向所述目标航天器的质心，z轴指向所述目标航天器的轨道面法线方向，y轴垂直于所述x轴和所述z轴所组成的平面，所述y轴与所述x轴、所述z轴构成右手直角坐标系。

[0015] 进一步地，在所述基于多等碰撞概率线法的凸多边形航天器安全控制方法中，所述相对动力学方程如式1所示；

[0016]

[0017] 其中，r1-t＝[xt,yt]T和v1-t＝[vx,vy]T表示在所述目标航天器的轨道坐标系下所述跟踪航天器在t时刻的相对位置和相对速度， u1-t＝[ux,uy]T表示所述跟踪航天器在t时刻的控制输入，μ表示地球引力常数， a和n分别表示所述目标航天
器的轨道半长轴和平均角速度；

[0018] 假设状态矢量 式1改写为：

[0019]

[0020] 其中，A表示状态转移矩阵，B表示控制矩阵，

[0021]

[0022] 进一步地，在所述基于多等碰撞概率线法的凸多边形航天器安全控制方法中，采用线性协方差方法推导不确定性的协方差矩阵，包括：

[0023] 利用式5计算所述状态矢量X的解析解；

[0024]

[0025] 其中，t0表示初始时间，X0表示与t0对应的状态矢量初始值，Φ(t,t0)表示状态矩阵且满足式6；

[0026]

[0027] Φrr(t,t0)、Φrv(t,t0)、Φvr(t,t0)和Φvv(t,t0)满足式7；

[0028]

[0029] 利用式8和式9计算所述状态矢量X的不确定性协方差矩阵；

[0030] δX＝X(t)-Ε(X(t))  (8)

[0031]

[0032] 其中，N表示所述跟踪航天器所施加的控制脉冲的数量， 和 分别表示在所述目标航天器的轨道坐标系下的初始导航不确定性和控制不确定性的协方差矩阵，CδX表示所述状态矢量的不确定性协方差矩阵，Ε(X(t))表示所述状态矢量的平均值，δX为状态矢量偏差，E(δX)为所述状态矢量偏差的平均值。

[0033] 进一步地，在所述基于多等碰撞概率线法的凸多边形航天器安全控制方法中，基于所述目标航天器的轨道坐标系，根据所述跟踪航天器所处位置，将所述目标航天器划分为如式10所示的5个作用模块；

[0034]

[0035] 其中，[x,y]表示所述跟踪航天器在所述目标航天器的轨道坐标系下的相对位置坐标。

[0036] 进一步地，在所述基于多等碰撞概率线法的凸多边形航天器安全控制方法中，根据所述跟踪航天器的最终时刻的相对位置的不同，对施加于所述跟踪航天器的避障控制力进行分类分析。

[0037] 进一步地，在所述基于多等碰撞概率线法的凸多边形航天器安全控制方法中，设定所述跟踪航天器的最终时刻的相对位置为rep＝[xep yep]T；

[0038] 当xep＝0oryep＝0时，施加于所述跟踪航天器的避障控制力如式48所示；

[0039]

[0040] 其中，F50-repel+Fh0-repel,h＝1,2,3,4表示当所述跟踪航天器位于h＝1,2,3或4所对应的区域时，施加于所述跟踪航天器的避障控制力为F50-repel+Fh0-repel，F50-repel,h＝5表示当所述跟踪航天器位于h＝5所对应的区域时，施加于所述跟踪航天器的避障控制力为F50-repel。

[0041] 进一步地，在所述基于多等碰撞概率线法的凸多边形航天器安全控制方法中，当xep≠0andyep≠0时，根据相对位置rep＝[xep yep]T建立两个坐标系统O1-ix1-iy1-i(i＝1,2)，并根据相对位置rep＝[xep yep]T的不同分为如式49所示的四个情况进行避障控制力分析；

[0042]

[0043] 施加于所述跟踪航天器的避障控制力如式59所示；

[0044]

[0045] 其中，F50-repel+Fh1-repel+Fh2-repel,h＝h1表示当所述跟踪航天器位于h＝h1所对应的区域时，施加于所述跟踪航天器的避障控制力为F50-repel+Fh1-repel+Fh2-repel，F50-repel+Fh0-repel,h＝1,2,3,4&h≠h1表示当所述跟踪航天器位于h＝1,2,3,4&h≠h1所对应的区域时，施加于所述跟踪航天器的避障控制力为F50-repel+Fh0-repel，F50-repel,h＝5表示当所述跟踪航天器位于h＝5所对应的区域时，施加于所述跟踪航天器的避障控制力为F50-repel。

[0046] 进一步地，在所述基于多等碰撞概率线法的凸多边形航天器安全控制方法中，根据式70所示条件，确定所述多等碰撞概率线法的参数的取值范围，以使所述跟踪航天器与所述目标航天器不会发生碰撞；

[0047]

[0048] 其中，R50、Ri0(i＝1,2,3,4)和Rij(i＝1,2,3,4,j＝0,1,2)表示相对应的作用模块的最小外包络圆的半径，ζ50、ζi0(i＝1,2,3,4)和ζij(i＝1,2,3,4,j＝0,1,2)表示相对应的作用模块的最小外包络圆的圆心到跟踪航天器的最短距离。

[0049] 进一步地，在所述基于多等碰撞概率线法的凸多边形航天器安全控制方法中，作用于所述跟踪航天器的总控制力利用式81计算；

[0050]

[0051] 其中，uTotal表示作用于所述跟踪航天器的总控制力， 表示最优控制，Frepel表示所述避障控制力，m表示所述跟踪航天器的质量。

[0052] 本发明技术方案的主要优点如下：

[0053] 本发明的基于多等碰撞概率线法的凸多边形航天器安全控制方法，能够解决目标航天器为凸多边形情况下的跟踪航天器的安全接近控制问题，并通过建立一套多等碰撞概率线法的参数选取原则，能便于多等碰撞概率线法的参数的选取；同时，通过改进LQR控制器，能得到更高的控制精度和更高的鲁棒性。

附图说明

[0054] 为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

[0055] 图1为本发明一个实施例的基于多等碰撞概率线法的凸多边形航天器安全控制方法的流程图；

[0056] 图2为本发明一个实施例的坐标系的示意图；

[0057] 图3为本发明一个实施例的目标航天器的作用模块的分布示意图；

[0058] 图4为本发明一个实施例的跟踪航天器在轨道坐标系Ⅰ部分运动时的受力分析示意图，其中xep＝0 or yep＝0；

[0059] 图5为本发明一个实施例的跟踪航天器在轨道坐标系Ⅰ部分运动时的受力分析示意图，其中xep≠0 and yep≠0；

[0060] 图6为本发明一个实施例的避障控制力作用时的跟踪航天器的受力分析示意图，其中xep＝0 or yep＝0；

[0061] 图7为本发明一个实施例的避障控制力作用时的跟踪航天器的受力分析示意图，其中xep≠0 and yep≠0。

具体实施方式

[0062] 为使本发明的目的、技术方案和优点更加清楚，下面将结合本发明具体实施例及相应的附图对本发明技术方案进行清楚、完整地描述。显然，所描述的实施例仅是本发明的一部分实施例，而不是全部的实施例。基于本发明的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

[0063] 以下结合附图，详细说明本发明实施例提供的技术方案。

[0064] 如附图1所示，本发明实施例提供了一种基于多等碰撞概率线法的凸多边形航天器安全控制方法，该方法用于实现目标航天器和跟踪航天器的安全接近控制，该方法包括如下步骤：

[0065] 1)建立目标航天器的轨道坐标系；

[0066] 2)在目标航天器的轨道坐标系下，建立跟踪航天器与目标航天器的相对动力学方程；

[0067] 3)进行相对动力学方程的不确定性传播分析，采用线性协方差方法推导不确定性的协方差矩阵；

[0068] 4)将目标航天器的几何外形分解为一个最大内包络圆和若干个任意形状部分，分别设定最大内包络圆和每个任意形状部分的外包络圆为一个作用模块，根据跟踪航天器与目标航天器的相对位置，确定多个作用模块中产生作用的作用模块；

[0069] 5)利用等碰撞概率线法建立多等碰撞概率线法以计算基于作用模块产生的施加于跟踪航天器的避障控制力；

[0070] 6)对避障控制力进行分析，确定多等碰撞概率线法的参数的取值范围，以使跟踪航天器与目标航天器不会发生碰撞；

[0071] 7)确定最终作用于跟踪航天器上的总控制力。

[0072] 具体地，以下通过具体示例对本发明实施例提供的基于多等碰撞概率线法的凸多边形航天器安全控制方法的各个步骤进行具体说明。

[0073] (1)建立目标航天器的轨道坐标系；

[0074] 具体地，在本发明实施例提供的基于多等碰撞概率线法的凸多边形航天器安全控制方法中，如附图2所示，采用O-XIYIZI表示历元J2000地球惯性坐标系，地球1地心为坐标原点，XI轴指向历元J2000春分点，地球1赤道平面为基本面，ZI轴指向地球1北极，YI轴与XI轴、ZI轴构成右手直角坐标系；

[0075] 相应地，采用o-xyz表示目标航天器的轨道坐标系(Local Vertical Local Horizontal，LVLH)，以目标航天器的质心为坐标原点o，x轴由地球地心指向目标航天器的质心，z轴指向目标航天器的轨道面法线方向，y轴垂直于x轴和z轴所组成的平面，y轴与x轴、z轴构成右手直角坐标系。

[0076] (2)在目标航天器的轨道坐标系下，建立跟踪航天器与目标航天器的相对动力学方程；

[0077] 具体地，设定采用的目标航天器参考轨道为圆轨道，由于轨道面内的运动能与垂直轨道面的运动分离，因此采用对同轨道面内的相对运动进行研究的方式对相对运动模型进行分析计算；由此，在本发明实施例提供的基于多等碰撞概率线法的凸多边形航天器安全控制方法中，跟踪航天器与目标航天器的相对动力学方程如式1所示；

[0078]

[0079] 其中，r1-t＝[xt,yt]T和v1-t＝[vx,vy]T表示在目标航天器的轨道坐标系下跟踪航天器在t时刻的相对位置和相对速度， u1-t＝[ux,uy]T表示跟踪航天器在t时刻的控制输入，μ表示地球引力常数， a和n分别表示目标航天器的轨道半长轴和
平均角速度；

[0080] 进一步地，假设状态矢量 式1可改写为：

[0081]

[0082] 其中，A表示状态转移矩阵，B表示控制矩阵，控制矩阵B表示两个相互独立的控制量。

[0083] (3)进行相对动力学方程的不确定性传播分析，采用线性协方差方法推导不确定性的协方差矩阵；

[0084] 具体地，在上述步骤(2)中，由于设定目标航天器参考轨道为圆轨道，因此，根据上述的状态矢量X的表达式，利用式5可计算状态矢量X的解析解；

[0085]

[0086] 其中，t0表示初始时间，X0表示与t0对应的状态矢量初始值，Φ(t,t0)表示状态矩阵且满足式6；

[0087]

[0088] Φrr(t,t0)、Φrv(t,t0)、Φvr(t,t0)和Φvv(t,t0)满足式7；

[0089]

[0090] 进一步地，当考虑导航不确定性和控制输出不确定性时，可以采用线性协方差方法推导不确定性的协方差矩阵，例如可以利用式8和式9计算状态矢量X的不确定性协方差矩阵；

[0091] δX＝X(t)-Ε(X(t))  (8)

[0092]

[0093] 其中，N表示跟踪航天器所施加的控制脉冲的数量， 和 分别表示在目标航天器的轨道坐标系下的初始导航不确定性和控制不确定性的协方差矩阵，CδX表示状态矢量的不确定性协方差矩阵，Ε(X(t))表示状态矢量的平均值，δX为状态矢量偏差，E(δX)为状态矢量偏差的平均值。

[0094] (4)将目标航天器的几何外形分解为一个最大内包络圆和若干个任意形状部分，分别设定最大内包络圆和每个任意形状部分的外包络圆为一个作用模块，根据跟踪航天器与目标航天器的相对位置，确定多个作用模块中产生作用的作用模块；

[0095] 具体地，将目标航天器的几何外形分解为一个最大内包络圆和若干个任意形状的部分，且每个任意形状部分都具有相对应的最小外包络圆，当跟踪航天器围绕目标航天器飞行时，设定最大内包络圆和每个任意形状部分的外包络圆为作用模块，根据跟踪航天器与目标航天器的相对位置，确定多个作用模块中产生作用的作用模块，从而通过作用模块分别采用传统等碰撞概率线法产生避障控制力，并且共同施加在跟踪航天器上。

[0096] 为了便于对本发明实施例提供的基于多碰撞概率线法的凸多边形航天器安全控制方法进行具体说明，以下假设目标航天器的外形为凸多边形且凸多边形为正方形；当然，针对目标航天器的外形为非正方形的其他凸多边形时，本发明实施例提供的基于多碰撞概率线法的凸多边形航天器安全控制方法仍然适用。

[0097] 如附图3所示，将目标航天器的轨道坐标系分为I、Ⅱ、Ⅲ和Ⅳ四部分，将目标航天器分为5个作用模块，当跟踪航天器围绕目标航天器运动时，产生作用的作用模块如式10所示；

[0098]

[0099] 其中，[x,y]表示跟踪航天器在目标航天器的轨道坐标系下的相对位置坐标，Part5表示所述目标航天器中最大内包络圆对应部分，Part1表示在目标航天器的轨道坐标系下目标航天器I部分中除去最大内包络圆部分后的对应部分，Part2表示在目标航天器的轨道坐标系下目标航天器Ⅱ部分中除去最大内包络圆部分后的对应部分，Part3表示在目标航天器的轨道坐标系下目标航天器Ⅲ部分中除去最大内包络圆部分后的对应部分，
Part4表示在目标航天器的轨道坐标系下目标航天器Ⅳ部分中除去最大内包络圆部分后的对应部分。

[0100] (5)利用等碰撞概率线法建立多等碰撞概率线法以计算基于作用模块产生的施加于跟踪航天器的避障控制力；

[0101] 具体地，首先根据跟踪航天器的最终时刻的相对位置的不同，对施加于跟踪航天器的避障控制力进行分类分析。

[0102] 假设跟踪航天器在最终时刻的相对位置位于目标航天器的表面，且设定跟踪航天器的最终时刻的相对位置为rep＝[xep yep]T；为了实现目标航天器与跟踪航天器的安全接近控制，以下根据相对位置的不同，分两类情况对施加于跟踪航天器的避障控制力进行分析。

[0103] 第一类情况：当xep＝0 or yep＝0时

[0104] 如附图4所示，当跟踪航天器在轨道坐标系I部分运动时，仅考虑位于轨道坐标系I部分的目标航天器的作用模块Part1和作用模块Part5与跟踪航天器的影响，此时设定基于作用模块Part5和作用模块Part1产生的作用力分别为F10-repel和F50-repel。同时，建立坐标系O1x1y1，其中，该坐标系的x1轴由LVLH坐标系的原点指向作用模块Part1的最小外包络圆的圆心，y1轴由x1轴绕着LVLH坐标系的oz轴逆时针旋转90度得到；则可利用式11，将跟踪航天器在LVLH坐标系下的相对位置转化为在O1x1y1坐标系下的相对位置。

[0105] r10-t＝W1{r1-t-[x10 y10]T}  (11)

[0106] 其中：

[0107]

[0108]

[0109]

[0110] 其中，r10＝[x10,y10]T表示作用模块Part1的最小外包络圆的圆心在O1x1y1坐标系下的相对位置。

[0111] 假设目标航天器的姿态为固定不动，此时可利用式15，将跟踪航天器在LVLH坐标系中的相对速度转化为在O1x1y1坐标系下的相对速度。

[0112] v10-t＝W1v1-t  (15)

[0113] 同时，由于协方差矩阵 是在LVLH坐标系下定义的，可利用式16将其转化到O1x1y1坐标系。

[0114]

[0115] 而后，可利用对角化处理矩阵G1-tW1-1，通过式17将在O1x1y1坐标系下的相对位置不确定性的协方差矩阵转化为在O1x1dy1d坐标系下的对角化矩阵；

[0116]

[0117] 其中，G1-t表示正交化转移矩阵，O1x1dy1d坐标系由O1x1y1坐标系基于矩阵G1-tW1-1变换得到。

[0118] 相应地，可利用式18和式19得到跟踪航天器在O1x1dy1d坐标系下的相对位置和相对速度；

[0119] r10-td＝G1-t{r1-t-[x10 y10]T}  (18)

[0120] v10-td＝G1-tv1-t  (19)

[0121] 其中，r10-td＝[x10-td,y10-td]T和v10-td＝[vx10-td,vy10-td]T分别表示跟踪航天器在O1x1dy1d坐标系下的相对位置和相对速度。

[0122] 在上述设定中，当跟踪航天器在轨道坐标系I部分运动时，基于作用模块Part1在O1x1dy1d坐标系下产生的避障控制力为F10-repel，F10-repel可表示为：

[0123] F10-repel＝F10-paral+F10-perpen  (20)

[0124] 其中，F10-paral和F10-perpen可利用下述公式计算获取；

[0125]

[0126]

[0127]

[0128]

[0129] 同时，根据式11-式14，可得：

[0130]

[0131] 其中：

[0132]

[0133]

[0134] 则利用式25-式27，可得：

[0135]

[0136]

[0137] 式中，r1-td表示跟踪航天器在Oxdyd坐标系下的相对位置， 和 表示与相对位置矢量r1-td方向相反和垂直的单位向量。

[0138] 其中， 矩阵在Oxdyd坐标系下是对角矩阵。

[0139] 基于上述式28和式29，经计算可得：

[0140]

[0141]

[0142] 其中，在上述式20至式31中，F10-repel表示作用模块Part1在O1x1dy1d坐标系下产生的避障控制力，F10-paral和F10-perpen分别表示F10-repel在相对位置r10-td平行反方向上和垂直方向上的分量，D010表示在O1x1dy1d坐标系下的作用模块Part1周围的影响区域的直径；r10-td表示在O1x1dy1d坐标系下的相对位置r10-td的模的大小； 和 表示在O1x1dy1d坐标系下的相对速度在相对位置矢量r10-td上的平行反方向分量与垂直分量； 和 表示与相对位置矢量r10-td方向相反与垂直的单位向量，λ0和d0为确定避障控制力的量级的正常数，amax为跟踪航天器的最大加速度， 和 表示与相对位置矢量r1-t方向相反与垂直的单位
向量。

[0143] 根据上述的跟踪航天器在轨道坐标系I部分运动时的避障控制力的分析，同理可得，当跟踪航天器在轨道坐标系Ⅱ、Ⅲ和Ⅳ部分运动时，基于作用模块Part2、作用模块Part3和作用模块Part4在O1x1dy1d坐标系下产生的避障控制力如式32所示；

[0144] Fi0-repel＝Fi0-paral+Fi0-perpen,(i＝2,3,4)  (32)

[0145] 设定：Fi0-repel(i＝2,3,4)表示作用模块Part2、作用模块Part3和作用模块Part4在O1x1dy1d坐标系下产生的避障控制力，ri0＝[xi0 yi0]T,(i＝2,3,4)为与不同作用模块相对应的最小外包络圆圆心在LVLH坐标系下的相对位置，ri0-td＝[xi0-td,yi0-td]T和vi0-td＝[vxi0-td,vyi0-td]T(i＝2,3,4)分别表示跟踪航天器在O1x1dy1d坐标系下的相对位置和相对速度，D0i0,(i＝2,3,4)表示在O1x1dy1d坐标系下相对应的作用模块周围的影响区域的直径，ri0-td,(i＝2,3,4)为与不同作用模块相对应的最小外包络圆圆心的相对位置ri0-td＝[xi0-td,yi0-td]T的模的大小， 和 表示在O1x1dy1d坐标系下的相对速度在相对位
置矢量ri0-td＝[xi0-td,yi0-td]T上的平行反方向分量与垂直分量，Fi0-paral和Fi0-perpen分别表示Fi0-repel在相对位置ri0-td平行反方向上和垂直方向上的分量。

[0146] ri0-td和vi0-td可通过下述式33和式34计算；

[0147] ri0-td＝G1-t{r1-t-[xi0yi0]T},(i＝2,3,4)  (33)

[0148] vi0-td＝G1-tv1-t,(i＝2,3,4)  (34)

[0149] 同理，参照F10-repel的计算过程，可得：

[0150]

[0151]

[0152]

[0153]

[0154]

[0155]

[0156] 进一步地，根据上述的避障控制力的分析计算，同理可得，基于作用模块Part5在O1x1dy1d坐标系下产生的避障控制力如式41所示；

[0157] F50-repel＝F50-paral+F50-perpen  (41)

[0158] 设定：F50-repel表示作用模块Part5在O1x1dy1d坐标系下产生的避障控制力，r50-td＝[x50-td,y50-td]T和v50-td＝[vx50-td,vy50-td]T分别表示跟踪航天器在O1x1dy1d坐标系下的相对位置和相对速度，D050表示在O1x1dy1d坐标系下作用模块Part5周围的影响区域的直径，r50-td为相对位置r50-td的模的大小， 和 表示在O1x1dy1d坐标系下的相对速度在相对位置矢量r50-td上的平行反方向分量与垂直分量，F50-paral和F50-perpen分别表示F50-repel在相对位置r50-td平行反方向上和垂直方向上的分量。

[0159] r50-td和v50-td可通过下述式42和式43计算；

[0160] r50-td＝G1-tr1-t  (42)

[0161] v50-td＝G1-tv1-t  (43)

[0162] 同理，参照F10-repel的计算过程，可得：

[0163]

[0164]

[0165]

[0166]

[0167] 由此，基于上述的避障控制力的计算分析，当xep＝0 or yep＝0时，施加于跟踪航天器的避障控制力如式48所示；

[0168]

[0169] 其中，F50-repel+Fh0-repel,h＝1,2,3,4表示当跟踪航天器位于h＝1,2,3或4所对应的区域时，施加于跟踪航天器的避障控制力为F50-repel+Fh0-repel，F50-repel,h＝5表示当跟踪航天器位于h＝5所对应的区域时，施加于跟踪航天器的避障控制力为F50-repel。

[0170] 第二类情况：当xep≠0 and yep≠0时

[0171] 如附图5所示，当跟踪航天器在轨道坐标系I部分运动时，仅考虑位于轨道坐标系I部分的目标航天器的作用模块Part1和作用模块Part5与跟踪航天器的影响，此时设定基于作用模块Part5产生的作用力分别为F50-repel；同时根据相对位置rep＝[xep yep]T，建立两个T坐标系统O1-ix1-iy1-i(i＝1,2)，并根据相对位置rep＝[xep yep]的不同分为如式49所示的四个情况进行避障控制力分析；

[0172]

[0173] 具体地，如附图5所示，例如当满足条件h1＝1时，此时跟踪航天器在LVLH坐标系I部分运动时，此时根据相对位置rep＝[xep yep]T将作用模块Part1划分为作用模块1-1和作用模块1-2，并分别建立O1-1x1-1y1-1坐标系和O1-2x1-2y1-2坐标系。

[0174] 同理，当满足h1其他条件，可参照上述的方式进行作用模块Part2、Part3和Part4的模块划分。

[0175] 参照当xep＝0 or yep＝0时，跟踪航天器在轨道坐标系I部分运动时的避障控制力的分析，同理可得，当跟踪航天器在轨道坐标系I、Ⅱ、Ⅲ和Ⅳ部分运动时，基于相对应的作用模块产生的避障控制力如式50所示；

[0176] Fij-repel＝Fij-paral+Fij-perpen(i＝1,2,3,4,j＝1,2)  (50)

[0177] 设定：Fij-repel(i＝1,2,3,4,j＝1,2)表示相对应的作用模块在O1-ix1-iy1-i(i＝1,2)坐标系下产生的避障控制力，rij＝[xij yij]T(i＝1，2,3,4,j＝1,2)为与不同作用模块相对应的最小外包络圆圆心在LVLH坐标系下的相对位置，rij-td＝[xij-td,yij-td]T和vij-td＝T[vxij-td,vyij-td] (i＝1，2,3,4,j＝1,2)分别表示跟踪航天器在O1x1dy1d坐标系下的相对位置和相对速度，D0ij(i＝1,2,3,4,j＝1,2)表示在O1x1dy1d坐标系下相对应的作用模块周围的影响区域的直径，rij-td(i＝1,2,3,4,j＝1,2)为与不同作用模块相对应的最小外包络圆圆心的相对位置rij-td的模的大小， 和 表示在
O1x1dy1d坐标系下的相对速度在相对位置矢量rij-td上的平行反方向分量与垂直分量，Fij-paral和Fij-perpen分别表示Fij-repel在相对位置rij-td平行反方向上和垂直方向上的分量。

[0178] rij-td和vij-td可通过下述式51和式52计算；

[0179] rij-td＝G1-t{r1-t-[xijyij]T},(i＝1,2,3,4,j＝1,2)  (51)

[0180] vij-td＝G1-tv1-t,(i＝1,2,3,4,j＝1,2)  (52)

[0181] 同理，参照F10-repel的计算过程，可得：

[0182]

[0183]

[0184]

[0185]

[0186]

[0187]

[0188] 由此，基于上述的避障控制力的计算分析，当xep≠0 and yep≠0时，施加于跟踪航天器的避障控制力如式59所示；

[0189]

[0190] 其中，F50-repel+Fh1-repel+Fh2-repel,h＝h1表示当跟踪航天器位于h＝h1所对应的区域时，施加于跟踪航天器的避障控制力为F50-repel+Fh1-repel+Fh2-repel，F50-repel+Fh0-repel,h＝1,2,3,4&h≠h1表示当跟踪航天器位于h＝1,2,3,4&h≠h1所对应的区域时，施加于跟踪航天器的避障控制力为F50-repel+Fh0-repel，F50-repel,h＝5表示当跟踪航天器位于h＝5所对应的区域时，施加于跟踪航天器的避障控制力为F50-repel。

[0191] (6)对避障控制力进行分析，确定多等碰撞概率线法的参数的取值范围，以使跟踪航天器与目标航天器不会发生碰撞；

[0192] 基于上述的避障控制力的计算分析，在t时刻，当和 满足时，需要在跟踪
航天器上施加避障控制力。

[0193] 基于上述的式11至式59进行计算分析，可以获取下述条件：

[0194] M1ijd＞M2ijd＞0，M150d＞M250d＞0，λ0＞0，D0ij＞R0＞0，D050＞R0＞0，d0＞0，amax＞0，rij-td＞0，r50-td＞0， 其中，(i＝1，2,3,4,j＝0,1,2)。

[0195] 由于避障控制力在相对位置矢量的平行方向上的分量远远大于其在相对位置矢量的垂直方向上的分量，即：

[0196]

[0197]

[0198] 进一步地，根据式48和式59，总避障控制力可改写为如式62所示形式；

[0199] FTotal-repel＝FTotal-paral+FTotal-perpen  (62)

[0200] 其中：

[0201]

[0202]

[0203] 在上述公式中，避障控制力FTotal-repel在相对位置矢量的平行方向上的分量FTotal-paral与单位矢量 的方向相反，该分量的作用是阻止跟踪航天器指向目标航天器的运动；避障控制力FTotal-repel在相对位置矢量的垂直方向上的分量FTotal-perpen与单位矢量的方向相反，该分量的作用是控制跟踪航天器绕开障碍物。

[0204] 进一步地，由于避障控制力在相对位置矢量的平行方向上的分量远远大于其在相对位置矢量的垂直方向上的分量，在本发明实施例提供的基于多等碰撞概率线法的凸多边形航天器安全控制方法中，仅对相对位置矢量的平行方向上的分量进行受力分析。

[0205] 具体地，由于任意的总避障控制力FTotal-repel都包含了避障力分量F50-paral，因此，以下通过分析F50-paral为例，来说明本发明实施例提供的多等碰撞概率线法的有效性。

[0206] 首先，作用模块Part5在相对位置矢量的平行方向上的分量F50-paral对相对位置r50-td求偏导，可得：

[0207]

[0208] 而后，定义以相对位置r50-td为自变量的辅助函数h(r50-td)，辅助函数h(r50-td)为：

[0209]

[0210] 由式66可知，辅助函数是变量为r50-td的二阶函数，且其二阶参数为负；由此，可以得到辅助函数h(r50-td)的极大值为：

[0211]

[0212] 同时，由式66可知，该二次方程对称轴为0，且该二次方程含有两个零点；由于相对位置r50-td＞0，因此，假设x3-2＞0为该二次方程的零点，则相应可得：

[0213]

[0214] 基于上述假设，可得：当r50-td＞x3-2时，h(r50-td)＜0；当0＜r50-td＜x3-2时，h(r50-td)＞0。

[0215] 进一步地，假设x3-3为相对位置r50-td的一个极大值，即：

[0216]

[0217] 则可得：当x3-2＜r50-td＜x3-3时， Fparal是个递减函数；当0＜r50-td＜x3-2时， Fparal是个递增函数。

[0218] 如图6和图7所示，图6和图7分别给出了两种不同的最终时刻的相对位置时避障控制力作用下的跟踪航天器的受力分析。当跟踪航天器运动到目标航天器的影响区域的边界时，设定跟踪航天器的相对平行速度为 当跟踪航天器在影响区域内运动且时，避障控制力FTotal-repel的分量FTotal-paral施加于跟踪航天器，以阻止跟踪航天器飞向目标航天器。

[0219] 定义临界值ζ50为跟踪航天器与目标航天器的最短相对距离，即作用模块Part5的圆心到跟踪航天器的最短距离，如此，当跟踪航天器到达临界点ζ50时，跟踪航天器的平行相对速度 需降低至0，跟踪航天器才不会与目标航天器发生碰撞。

[0220] 因此，基于上述的分析，通过选取合适的多等碰撞概率线法的参数值，使式70所示条件满足时，跟踪航天器将不会与目标航天器发生碰撞；

[0221]

[0222] 其中，R50、Ri0(i＝1,2,3,4)和Rij(i＝1,2,3,4,j＝0,1,2)表示相对应的作用模块的最小外包络圆的半径，ζ50、ζi0(i＝1,2,3,4)和ζij(i＝1,2,3,4,j＝0,1,2)表示相对应的作用模块的最小外包络圆的圆心到跟踪航天器的最短距离。

[0223] (7)确定最终作用于跟踪航天器上的总控制力；

[0224] 本发明实施例中，基于上述的多等碰撞概率线法，通过改进LQR控制器对跟踪航天器进行跟踪预定轨迹，具体地，改进LQR控制器的过程如下所示：

[0225] 根据式30-31、式39-40、式46-47、式57-58以及式62-64，定义F3-repel为：

[0226]

[0227] 其中，b3-1和b3-2表示辅助参数。

[0228] 进一步地，将式71转变为：

[0229] F3-repel＝K3-1X  (72)

[0230] 其中， X表示上述的状态矢量。

[0231] 根据式73推导，可得：

[0232]

[0233] 其中， I2×2表示2×2阶单位矩阵。

[0234] 进一步地，定义矩阵K3-3为矩阵 的上界值，且K3-3为：

[0235]

[0236] 其中，

[0237] 进一步地，本发明实施例中，基于上述计算分析，针对本发明实施例提供的如式2所示的相对动力学方程，提出了一种与多等碰撞概率线法相对应的增强LQR控制器，该增强LQR控制器的能量函数J3-1如式78所示；

[0238]

[0239] 其中，X表示上述状态矢量，m表示跟踪航天器的质量，u表示跟踪航天器的控制输入，Q和R分别表示LQR控制器权矩阵。

[0240] 基于式78，根据最小值原则，可得到相应的最优控制

[0241]

[0242] 其中，K3-4＝-R-1BTS3-2表示反馈矩阵，B表示上述控制矩阵，S3-2表示矩阵，且S3-2满足黎卡提方程，可得：

[0243]

[0244] 其中，A表示上述状态转移矩阵。

[0245] 显而易见，K3-4/m2+Q>Q，因此，增强LQR控制器能得到更高的控制精度和更高的鲁棒性。

[0246] 进一步地，根据式2和式78，通过解黎卡提方程，可以得到最优控制

[0247] 由此，基于得到的最优控制 可以得到作用于跟踪航天器的总控制力，总控制力可以利用式81计算；

[0248]

[0249] 其中，uTotal表示作用于跟踪航天器的总控制力， 表示最优控制，Frepel表示避障控制力，m表示跟踪航天器的质量。

[0250] 可见，本发明实施例提供的基于多等碰撞概率线法的凸多边形航天器安全控制方法，能够解决目标航天器为凸多边形情况下的跟踪航天器的安全接近控制问题，并通过建立一套多等碰撞概率线法的参数选取原则，能便于多等碰撞概率线法的参数的选取；同时，通过改进LQR控制器，能得到更高的控制精度和更高的鲁棒性。

[0251] 需要说明的是，在本文中，诸如“第一”和“第二”等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。此外，本文中“前”、“后”、“左”、“右”、“上”、“下”均以附图中表示的放置状态为参照。

[0252] 最后应说明的是：以上实施例仅用于说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。