近年来，考虑避免碰撞的编队飞行航天器系统协同控制受到了广泛的研究和关注^[1-7].与鲁棒性较差的单个航天器系统相比较^[8]，多航天器编队飞行系统具有通信成本低、鲁棒性强、效率高等优点^[9]，可以突破单个航天器编队系统的局限，对探测所需要的控制目标有更广阔的前景和更有意义的研究价值.

在航天器编队飞行系统协同控制的实际应用中，要求航天器能获得相互信息以实现对期望状态的跟踪.然而在信息获取的过程中由于距离和传输设备的约束，通常存在通信延迟^[10]、信息无法实时准确获取^[11]、航天器之间存在碰撞等情况.在编队飞行过程中，由于存在燃料损耗，航天器质量等信息无法精确实时获取，且在信息传递时会受到各种外部扰动的影响.根据上述分析的问题，为了提高编队航天器间协同控制的性能和精度，设计相应的自适应协同避碰控制算法具有重要的研究意义.文献[12]中对航天器编队协同控制进行分析，同时考虑避免碰撞提出了3种不同的编队跟踪控制方法.在文献[13]中基于编队卫星姿态动力学模型方程，通过引入速度滤波器设计了能够保证速度误差渐近收敛的自适应协同控制器.文献[14]进一步考虑推力校正和系统增益存在不确定性的情况下，提出了自适应跟踪控制算法.在上述分析的协同控制律中均为全状态反馈，而在实际的航天器编队飞行任务应用中，为节约燃料消耗、防止航天器之间的相对速度测量装置发生故障导致整个系统不稳定，需要保证编队系统在无速度测量情况下的有效协同作用.文献[15]提出了一种仅利用领航-跟踪航天器编队中的相对姿态反馈跟踪相对转动的控制器.文献[16]考虑输入受限的情况下，针对多Euler-Lagrange系统提出了一种无速度测量的分布式有限时间跟踪控制算法.

在编队飞行过程中，航天器之间避免碰撞是实现多航天器编队跟踪控制目标的重要保证.近年来，在多无人机协同控制、机械系统同步控制、多无人车协同控制等方面引起了国内外广泛的研究和关注.碰撞规避的主要研究方法可以分为：博弈论、碰撞概率分析、数学规划、半定规划^[17]和势函数法^[18].而势函数因其简单直观的物理意义成为解决避免碰撞的一种重要方法^[19-20].文献[21]基于势函数的方法对模型不确定的编队飞行卫星系统设计了避免碰撞的自适应协同控制器，同时考虑了参考轨迹作用的问题.文献[22]设计了能保证Lagrange系统协同避碰的控制律，并考虑输入扰动的影响.在此基础上，文献[23]分析了Lagrange系统同时存在外部干扰和参数不确定的情况下的稳定性，并利用势函数的方法设计相应的鲁棒自适应控制算法.

受上述文献的启发，本文针对航天器质量未知和存在通信时延的情况，考虑反馈状态对航天器编队飞行系统的影响，分别对于全状态反馈和无速度测量设计两种自适应协同避碰控制律.首先，本文考虑了存在通信时延的航天器编队飞行系统的自适应协同避碰控制问题，在设计相应控制律时，选取合理满足系统的势函数，实现避免碰撞的控制目标；其次，考虑航天器编队反馈状态的角度，分别提出了全状态反馈和无速度测量的自适应协同避碰控制律，并同时考虑了航天器质量不确定性和通信时延，对于参数不确定有较好的自适应性; 最后通过仿真分析证明了所设计控制律的有效性.

1 基础理论 1.1 相对运动非线性模型

在本文中每个航天器均假设为刚体结构，则n个航天器绕参考航天器编队飞行情况示意图，如图 1所示.其中F_I$ \triangleq ${O_IX_IY_IZ_I}为赤道惯性坐标系，其原点O_I代表地球的中心.F_l$ \triangleq ${o_lx_ly_lz_l}为参考航天器轨道坐标系，用来描述编队航天器间的相互运动，其原点o_l位于参考航天器的质心，x_l轴沿着地球质心指向参考航天器的矢量方向，y_l轴垂直于x_l轴并位于参考航天器轨道平面，z_l轴通过右手定则获得. r_l为参考航天器相对于地球质心的位置矢量^[24].

在不考虑参考航天器所受主动力控制的作用下，假设参考航天器在椭圆轨道上运行.r_i、v_i分别为第i个航天器相对与参考航天器的在编队飞行系统中的位置和速度，则第i个航天器非线性动力学模型的相对运动方程为：

$ {{\mathit{\boldsymbol{\dot r}}}_i} = {\mathit{\boldsymbol{v}}_i}, $

(1)

$ \begin{array}{l} {m_{{f_i}}}{{\mathit{\boldsymbol{\dot v}}}_i} = {\mathit{\boldsymbol{C}}_i}\left( {{{\dot \theta }_{\rm{c}}}} \right){\mathit{\boldsymbol{v}}_i} + {\mathit{\boldsymbol{D}}_i}\left( {{{\ddot \theta }_c},{{\dot \theta }_c},{r_{{f_i}}}} \right){\mathit{\boldsymbol{r}}_i} + \\ {\mathit{\boldsymbol{n}}_i}\left( {{\mathit{\boldsymbol{r}}_{{f_i}}},{\mathit{\boldsymbol{r}}_l}} \right) + {\mathit{\boldsymbol{d}}_i} + {\mathit{\boldsymbol{f}}_i}, \end{array} $

(2)

其中：

$ {\mathit{\boldsymbol{C}}_i}\left( {{{\dot \theta }_{\rm{c}}}} \right) = 2{m_{{f_i}}}{{\dot \theta }_{\rm{c}}}\left[ {\begin{array}{*{20}{c}} 0&1&0\\ { - 1}&0&0\\ 0&0&0 \end{array}} \right], $

(3)

$ \begin{array}{l} {\mathit{\boldsymbol{D}}_i}\left( {{{\ddot \theta }_c},{{\dot \theta }_c},{r_{{f_i}}}} \right) = - {m_{{f_i}}}\frac{\mu }{{\mathit{\boldsymbol{r}}_{{f_i}}^3}}{\mathit{\boldsymbol{I}}_{3 \times 3}} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{m_{{f_i}}}\left[ {\begin{array}{*{20}{c}} {\dot \theta _c^2}&{{{\ddot \theta }_c}}&0\\ {{{\ddot \theta }_c}}&{\dot \theta _c^2}&0\\ 0&0&0 \end{array}} \right], \end{array} $

(4)

$ {\mathit{\boldsymbol{n}}_i}\left( {{r_{{f_i}}},{r_i}} \right) = \mu {m_{{f_i}}}{\left[ { - \frac{{{r_l}}}{{r_{{f_i}}^3}} + \frac{1}{{r_l^2}}\quad 0\quad 0} \right]^{\rm{T}}}, $

(5)

$ {\mathit{\boldsymbol{d}}_i} = {\mathit{\boldsymbol{d}}_{{f_i}}} - \frac{{{m_{{f_i}}}}}{{{m_l}}}{\mathit{\boldsymbol{d}}_l}, $

(6)

式中：d_l、d_fi分别为参考航天器和第i个航天器所受到的外部干扰总和; θ_c为参考航天器的真近点角; r_l、r_fi分别为参考航天器和第i个航天器到地球中心的距离; I_3×3为适当维数的单位矩阵，μ为地球引力常数.

为能够更好地设计本文所提出的自适应协同避碰控制律，利用Euler-Lagrange方程来描述编队飞行航天器系统的相对运动^[25]：

$ {m_{{f_i}}}{{\mathit{\boldsymbol{\ddot r}}}_i} + \mathit{\boldsymbol{A}}{{\mathit{\boldsymbol{\dot r}}}_i} + \mathit{\boldsymbol{B}}{\mathit{\boldsymbol{r}}_i} + \mathit{\boldsymbol{C}} = {\mathit{\boldsymbol{d}}_i} + {\mathit{\boldsymbol{f}}_i}, $

(7)

式中：

$ \mathit{\boldsymbol{A}} = - {\mathit{\boldsymbol{C}}_i}\left( {{{\dot \theta }_c}} \right),\mathit{\boldsymbol{B}} = - {\mathit{\boldsymbol{D}}_i}\left( {{\theta _c},{\theta _c},{r_{{f_i}}}} \right),\mathit{\boldsymbol{C}} = - {\mathit{\boldsymbol{n}}_i}\left( {{r_{{f_i}}},{r_l}} \right). $

1.2 基础图论

本文中用加权无向图G=(V, E, A)来描述编队航天器之间的信息传递^[26]：V={v₁, …, v_n}为节点集，E⊆V×V为连接节点的边集，A=[a_ij]∈ R^n×n为加权邻接矩阵.若(j, i)∈E，且(i, j)∈E，则称G为无向图.如果第j个节点有指向第i个节点的路径且存在信息传递，则第j个节点就有一条边指向第i个节点，记为(j, i)∈E.在无向图中，如果任意两点之间存在路径连接，则无向图具有连通性^[26].

1.3 势函数设计与控制目标设定 1.3.1 势函数设计

本文利用势函数的方法来完成多航天器编队飞行避免碰撞的控制目标，使系统达到期望的控制状态.势函数的主要优点在于：在期望的状态下具有全局最小值，且具有较高势函数值的区域代表障碍碰撞等限制情况.受文献[27]启发，设计势函数的主要思想如下：首先，定义一个反映系统约束条件的标量势函数.其次，利用势函数的梯度值表示施加在航天器编队上避免碰撞的控制力.最后，设计适当的控制律使系统势函数的值具有减小的趋势，并且应用Lyapunov稳定性理论，保证编队系统的位置和速度均同时收敛于期望状态，满足避免碰撞的控制目标.

R_ij、r_ij分别为编队系统通信和最小安全区域的距离，选取如下的避免碰撞势函数：

$ {\mathit{\boldsymbol{V}}_{ij}}\left( {{r_i},{r_j}} \right) = {\left( {\min \left\{ {0,\frac{{{{\left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\|}^2} - R_{ij}^2}}{{{{\left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\|}^2} - r_{ij}^2}}} \right\}} \right)^2},i \ne j $

式中, V_ij(r_i, r_j)为一个非负连续的可微函数.

V_ij沿向量r_i的梯度表示为

$ \nabla {r_i}{\mathit{\boldsymbol{V}}_{ij}} = \left\{ {\begin{array}{*{20}{l}} {4\frac{{\left( {R_{ij}^2 - r_{ij}^2} \right)\left( {{{\left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\|}^2} - R_{ij}^2} \right)}}{{{{\left( {{{\left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\|}^2} - r_{ij}^2} \right)}^3}}}{{\left( {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right)}^{\rm{T}}},}\\ {{\rm{if}}:{r_{ij}} < \left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\| < {R_{ij}};}\\ 0\\ {{\rm{if}}:\left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\| \ge {R_{ij}},\left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\| \le {r_{ij}}.} \end{array}} \right. $

1.3.2 控制目标设定

在航天器编队飞行控制系统中，假设编队中各航天器均能获取到期望位置r_i^d和期望速度$\mathit{\boldsymbol{\dot r}}_i^d$.定义第i个航天器的位置误差和速度误差如下：

$ {{\mathit{\boldsymbol{\tilde r}}}_i} = {\mathit{\boldsymbol{r}}_i} - \mathit{\boldsymbol{r}}_i^d, $

$ {{\mathit{\boldsymbol{\dot {\tilde r}}}}_i} = {{\mathit{\boldsymbol{\dot r}}}_i} - \mathit{\boldsymbol{\dot r}}_i^d. $

在航天器质量参数不确定性、通信拓扑结构存在时延的情况下，将避免碰撞势函数与自适应协同控制相互结合，实现航天器在编队飞行过程中保持编队稳定不发生碰撞并保证系统的拓扑连通性，最终实现对期望位置和期望速度的跟踪, 即当t→∞时，对于$ \forall i, j = 1, \cdots , n, {\mathit{\boldsymbol{\widetilde r}}_i} \to {\mathit{\boldsymbol{\widetilde r}}_j} \to 0, {\mathit{\boldsymbol{\dot {\tilde r}}}_i} \to {\mathit{\boldsymbol{\dot {\tilde r}}}_j} \to 0 $.

2 控制律设计

在本文中，为实现编队航天器飞行系统的自适应协同避免碰撞控制，分别对于全状态反馈和无速度测量情况，设计相应的控制律.在设计控制律时，为实现编队系统的稳定性应考虑：1)编队航天器对期望位置和期望速度的跟踪；2)编队航天器的协同控制作用使其状态趋于稳定；3)保证编队航天器之间不发生碰撞.

2.1 全状态反馈自适应协同控制律

将辅助变量定义如下:

$ \mathit{\boldsymbol{\dot z}}_i^r = \mathit{\boldsymbol{\dot r}}_i^d - \eta \left( {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_d}} \right) = \mathit{\boldsymbol{\dot r}}_i^d - \eta {{\mathit{\boldsymbol{\tilde r}}}_i}. $

通过利用滑模变结构理论，第i个航天器可以表示为

$ {\mathit{\boldsymbol{\varpi}} _i} = {{\mathit{\boldsymbol{\dot r}}}_i} - \mathit{\boldsymbol{\dot z}}_i^r = {{\mathit{\boldsymbol{\dot {\tilde r}}}}_i} + \mathit{\boldsymbol{\eta }}{{\mathit{\boldsymbol{\tilde r}}}_i}, $

式中η为正常数.

定理1  考虑航天器编队相对运动模型(1)~(6)，如果通信拓扑结构为无向图，假设第i个航天器所受到的外部扰动d_i=0，在存在航天器质量参数不确定性、通信时延的情况下，同时考虑避免碰撞的控制目标，将航天器编队飞行系统控制律设计如下：

$ \begin{array}{l} {\mathit{\boldsymbol{f}}_i} = {{\hat m}_{{f_i}}}{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i} - {\delta _i}{\mathit{\boldsymbol{\varpi}} _i} - \\ \;\;\;\;\;\;\beta \sum\limits_{j = 1}^n {{a_{ij}}} \left[ {{\mathit{\boldsymbol{\varpi}} _i} - {\mathit{\boldsymbol{\varpi}} _j}\left( {t - {T_{ij}}} \right)} \right] - \sum\limits_{j = 1}^n \nabla {r_i}{\mathit{\boldsymbol{V}}_{ij}}, \end{array} $

(8)

$ {{\dot {\hat m}}_{{f_i}}} = - {\psi _i}\mathit{\boldsymbol{\varpi}} _i^{\rm{T}}{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i}, $

(9)

$ {\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i} = - \frac{1}{{{m_{{f_i}}}}}{\mathit{\boldsymbol{D}}_i}\left( {{{\ddot \theta }_c},{{\dot \theta }_c},{r_{{f_i}}}} \right){\mathit{\boldsymbol{r}}_i} - \frac{1}{{{m_{{f_i}}}}}{\mathit{\boldsymbol{n}}_i}\left( {{r_i},{r_c}} \right). $

(10)

式中：δ_i、ψ_i为正常数；a_ij为加权邻接矩阵元素；T_ij为第i个航天器获取第j个航天器相对信息的定常时延，其中包括位置时延$\tilde{\boldsymbol{r}}_{j}\left(t-T_{i j}\right)$和速度时延$ {\mathit{\boldsymbol{\dot {\tilde r}}}_j}(t - {T_{ij}}) $.

定义第i个航天器的质量估计误差为

$ {{\tilde m}_{{f_i}}} = {m_{{f_i}}} - {{\hat m}_{{f_i}}}, $

式中，${\dot {\tilde m}_{{f_i}}} = - {\dot {\hat m}_{{f_i}}} = {\psi _i}\mathit{\boldsymbol{\varpi}} _i^{\rm{T}}{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i}$.

根据势函数设计与控制目标设定提出的势函数方法和文献[27]所提出的相关内容，本文对势函数作如下假设.

假设1  若拓扑结构为无向图且V_ij(‖r_i-r_j‖)对称，则可以得到▽_riV_ij(‖r_i-r_j‖)=-▽_rjV_ij(‖r_i-r_j‖).

假设2    $\lim \mathop {{\mathit{\boldsymbol{V}}_{ij}}\left( {\left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\|} \right)}\limits_{\left\| {{r_i} - {r_j}} \right\| \to {c_{ij}}} = + \infty $，其中：V_ij为一个连续可微的非负函数，c_ij为第i个航天器与第j个航天器之间避免碰撞的最小距离.

在本文中，假设航天器所受到外部扰动d_i=0，将控制律(8)~(10)代入式(2), 则闭环系统可以写为

$ \begin{array}{l} {m_{{f_i}}}{{\mathit{\boldsymbol{\dot v}}}_i} = {\mathit{\boldsymbol{C}}_i}\left( {{{\dot \theta }_c}} \right){\mathit{\boldsymbol{v}}_i} + \beta \sum\limits_{j = 1}^n {{a_{ij}}} \left[ {{\mathit{\boldsymbol{\varpi}} _i} - {\mathit{\boldsymbol{\varpi}} _j}\left( {t - {T_{ij}}} \right)} \right] - \\ \;\;\;\;\;\;\;\;\;\;\;\sum\limits_{j = 1,j \ne i}^n {{\nabla _{{r_i}}}} {\mathit{\boldsymbol{V}}_{ij}} - {{\tilde m}_{{f_i}}}{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i} - {\delta _i}{\mathit{\boldsymbol{\varpi}} _i}. \end{array} $

证明1  定义如下Lyapunov函数:

$ \begin{array}{l} V(t) = \frac{1}{2}\sum\limits_{i = 1}^n {\mathit{\boldsymbol{\varpi}} _i^{\rm{T}}} {m_{{f_i}}}{\mathit{\boldsymbol{\varpi}} _i} + \frac{1}{2}\sum\limits_{i = 1}^n {\frac{1}{{{\psi _i}}}\tilde m_{{f_i}}^2} + \\ \;\;\;\;\;\;\;\;\;\frac{\beta }{2}\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{a_{ij}}} } \int_{t - {T_{ij}}}^t {\mathit{\boldsymbol{\varpi}} _j^{\rm{T}}} (\rho ){\mathit{\boldsymbol{\varpi}} _j}(\rho )d\rho + \\ \;\;\;\;\;\;\;\;\;\sum\limits_{i = 1}^n {\sum\limits_{j > i} {{\mathit{\boldsymbol{V}}_{ij}}} } \left( {\left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\|} \right). \end{array} $

(11)

将式(11)对时间求导，可以得到:

$ \begin{array}{l} \dot V\left( t \right) = \sum\limits_{i = 1}^n {\mathit{\boldsymbol{\varpi}} _i^{\rm{T}}} {m_{{f_i}}}{{\dot {\mathit{\boldsymbol{\varpi}}} }_i}\sum\limits_{i = 1}^n {\frac{1}{{{\psi _i}}}} {{\tilde m}_{{f_i}}}\tilde m_{{f_i}}^2 + \\ \;\;\;\;\;\;\;\;\;\;\;\frac{\beta }{2}\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{a_{ij}}} } \left[ {\mathit{\boldsymbol{\varpi}} _j^{\rm{T}}{\mathit{\boldsymbol{\varpi}} _j} - } \right.\\ \;\;\;\;\;\;\;\;\;\;\;\;\left. {\mathit{\boldsymbol{\varpi}} _j^{\rm{T}}\left( {t - {T_{ij}}} \right){\mathit{\boldsymbol{\varpi}} _j}\left( {t - {T_{ij}}} \right)} \right] + \\ \;\;\;\;\;\;\;\;\;\;\;\;\sum\limits_{i = 1}^n {\sum\limits_{j > i} {{{\left( {{\nabla _{{r_i}}}{\mathit{\boldsymbol{V}}_{ij}}} \right)}^{\rm{T}}}} } {{\mathit{\boldsymbol{\dot r}}}_i} + {\left( {{\nabla _{ij}}{\mathit{\boldsymbol{V}}_{ij}}} \right)^{\rm{T}}}{{\mathit{\boldsymbol{\dot r}}}_j}. \end{array} $

(12)

在对式(12)推导过程中，需要分析如下几点：

1) 由于${\mathit{\boldsymbol{C}}_i}\left( {{{\mathop \theta \limits^. }_c}} \right)$是反对称矩阵，则$\sum\limits_{i = 1}^n {\mathit{\boldsymbol{\varpi}} _i^{\rm{T}}} {\mathit{\boldsymbol{C}}_i}\left( {{{\dot \theta }_c}} \right){\mathit{\boldsymbol{\varpi}} _i} = 0$.

2) 若航天器编队系统的通信拓扑结构为无向图，则a_ij=a_ji，进一步可以得到如下等式成立:

$ \sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{a_{ij}}} } \mathit{\boldsymbol{\varpi}} _i^{\rm{T}}{\mathit{\boldsymbol{\varpi}} _i} = \sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{a_{ij}}} } \mathit{\boldsymbol{\varpi}} _j^{\rm{T}}{\mathit{\boldsymbol{\varpi}} _i}. $

3) 若通信拓扑结构为无向图且▽_riV_ij=-▽_rjV_ij，则

$ \begin{array}{*{20}{c}} {\left. {\sum\limits_{i = 1}^n {\left[ {{{\left( {{\nabla _{{r_i}}}{\mathit{\boldsymbol{V}}_{ij}}} \right)}^{\rm{T}}}{{\mathit{\boldsymbol{\dot r}}}_i}} \right]} + {{\left( {{\nabla _{{r_i}}}{\mathit{\boldsymbol{V}}_{ij}}} \right)}^{\rm{T}}}{{\mathit{\boldsymbol{\dot r}}}_j}} \right] - }\\ {\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {\mathit{\boldsymbol{r}}_i^{\rm{T}}} } {\nabla _{{r_i}}}{\mathit{\boldsymbol{V}}_{ij}} = 0.} \end{array} $

4) 时间延迟满足如下不等式:

$ \begin{array}{l} - \frac{\beta }{2}\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{a_{ij}}} } \mathit{\boldsymbol{\varpi}} _i^{\rm{T}}{\mathit{\boldsymbol{\varpi}} _i} + \sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{a_{ij}}} } \mathit{\boldsymbol{\varpi}} _i^{\rm{T}}{\mathit{\boldsymbol{\varpi}} _j}\left( {t - {T_{ij}}} \right) - \\ \;\;\;\frac{\beta }{2}\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{a_{ij}}} } \mathit{\boldsymbol{\varpi}} _j^{\rm{T}}\left( {t - {T_{ij}}} \right){\mathit{\boldsymbol{\varpi}} _j}\left( {t - {T_{ij}}} \right) \le 0. \end{array} $

根据上述分析，式(12)可以写为

$ \dot V(t) \le - \sum\limits_{i = 1}^n {\mathit{\boldsymbol{\varpi}} _i^{\rm{T}}} {\delta _i}{\mathit{\boldsymbol{\varpi}} _i} \le 0. $

因此$\dot{V}(t) \leqslant 0$，V(t)有界.由Barbalat引理^[28]分析系统稳定性，根据假设2和Lyapunov函数V(t)的定义，可知：在任意时间∀t≥0时，‖ r_i(t)- r_j(t)‖＞r_ij，则编队航天器不会发生碰撞.进一步考虑到V(t)有界，则结合${\mathit{\boldsymbol{\varpi}} _i} = {\mathit{\boldsymbol{\dot {\tilde r}}}_i} + \eta {\mathit{\boldsymbol{\tilde r}}_i}$，可知${\mathit{\boldsymbol{\tilde r}}_i}$，${\mathit{\boldsymbol{\dot {\tilde r}}}_i}$也是全局一致有界.因此由Barbalat引理可以得到：t→∞时，$ \tilde{\boldsymbol{r}}_{i} \rightarrow \tilde{\boldsymbol{r}}_{j} \rightarrow 0 $，$\dot{\tilde{\boldsymbol{r}}}_{i} \rightarrow \dot{\tilde{\boldsymbol{r}}}_{j} \rightarrow 0$，且$\mathop {\lim }\limits_{t \to \infty } \dot V(t) = 0$.因此航天器编队飞行闭环系统渐近稳定.

2.2 无速度反馈自适应协同控制律

在航天器编队飞行的实际任务中，航天器的相关信息需要先进的通信设备进行实时测量，为降低航天器之间的通信消耗成本，应考虑设计无速度测量的相应控制律.因此在本文中，针对无速度测量中${\mathit{\boldsymbol{\tilde r}}_i}$和$ {\mathit{\boldsymbol{\dot {\tilde r}}}_i} $未知的情况，设计无速度测量的自适应协同避免碰撞控制律，以实现期望的控制目标.

受文献[29]启发，考虑到通信时延的情况，设计一种无源滤波器：

$ {{\mathit{\boldsymbol{\dot {\hat g}}}}_i} = \mathit{\boldsymbol{ \boldsymbol{\varTheta} }}{\mathit{\boldsymbol{{\hat g}}}_i} + \sum\limits_{j = 1}^n {{a_{ij}}} \left[ {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}\left( {\mathit{\boldsymbol{t}} - {T_{ij}}} \right)} \right], $

(13)

$ {\mathit{\boldsymbol{h}}_i} = \mathit{\Lambda }\left\{ {\mathit{\boldsymbol{ \boldsymbol{\varTheta} }}{{\mathit{\boldsymbol{\hat g}}}_i} + \sum\limits_{j = 1}^n {{a_{ij}}} \left[ {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}\left( {t - {T_{ij}}} \right)} \right]} \right\}. $

(14)

式中：i=1, …, n，Θ ∈ R³是Hurwitz矩阵.Λ=Λ^T∈R^3×3为正定矩阵满足如下Lyapunov等式：

$ {\mathit{\boldsymbol{ \boldsymbol{\varTheta} }}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + \mathit{\boldsymbol{ \boldsymbol{\varLambda} \boldsymbol{\varTheta} }} = - \mathit{\boldsymbol{Q}}, $

式中，Λ、Q均为对称正定阵.

定理2  考虑航天器编队相对运动模型(7)，如果通信拓扑结构为无向图，假设第i个航天器所受到的外部扰动d_i=0，在存在航天器质量参数不确定性、通信时延及速度无法测量的情况下，同时考虑避免碰撞的控制目标，将航天器编队飞行系统控制律设计如下：

$ \begin{array}{*{20}{c}} {{\mathit{\boldsymbol{f}}_i} = - \beta \sum\limits_{j = 1}^{n} {{a_{ij}}} \left[ {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}\left( {t - {T_{ij}}} \right)} \right] - }\\ {{\mathit{\boldsymbol{h}}_i} + {{\hat m}_{{f_i}}}{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i} + \sum\limits_{j = 1,j \ne i}^n {{\nabla _{{r_i}}}} {\mathit{\boldsymbol{V}}_{ij}},} \end{array} $

(15)

$ {{\dot {\hat m}}_{{f_i}}} = - {\xi _i}\mathit{\boldsymbol{r}}_i^{\rm{T}}{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i}, $

(16)

式中, ξ_i为自适应更新律的正常数量.则有

$ {{\dot {\hat m}}_{{f_i}}} = {m_{{f_i}}}(0) - {\xi _i}\int_{r_i^{(0)}}^{r_i^{(t)}} {\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i^{\rm{T}}} {\rm{d}}{r_i}. $

引理1  文献[30]假设M∈Q^p×p，N∈Q^q×q，则下式成立：

1)(M⊗N)(A⊗B)= MA⊗NB.

2) 假设M和N为可逆矩阵，则满足(M⊗N)^-1= M^-1⊗N^-1.

3) 如果M和N是对称正定的，则M⊗N同样也为对阵矩阵.

与设计全状态反馈控制律相似，假设航天器编队受到的外部干扰为零，则闭环系统方程可写为

$ \begin{array}{*{20}{c}} {{m_{{f_i}}}{{\mathit{\boldsymbol{\ddot {\tilde r}}}}_i} + A{{\mathit{\boldsymbol{\dot {\tilde r}}}}_i} + B{{\mathit{\boldsymbol{\tilde r}}}_i} + \mathit{\boldsymbol{C}} + \beta \left( {\mathit{\boldsymbol{ \boldsymbol{\varXi} }} \otimes {\mathit{\boldsymbol{I}}_3}} \right){{\mathit{\boldsymbol{\tilde r}}}_i}\left( {t - {T_{ij}}} \right) + }\\ {{s_i} + {{\tilde m}_{{f_i}}}{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i} - \sum\limits_{j = 1,j \ne i}^n {{\nabla _{{r_i}}}} {\mathit{\boldsymbol{V}}_{ij}} = 0.} \end{array} $

证明  定义如下Lyapunov函数:

$ \begin{array}{l} V(t) = \frac{1}{2}\mathit{\boldsymbol{\dot {\hat g}}}_i^{\rm{T}}{\left( {\mathit{\boldsymbol{ \boldsymbol{\varXi} }} \otimes {\mathit{\boldsymbol{I}}_3}} \right)^{ - 1}}\left( {{\mathit{\boldsymbol{I}}_n} \otimes \mathit{\boldsymbol{ \boldsymbol{\varLambda} }}} \right){{\mathit{\boldsymbol{\dot {\hat g}}}}_i} + \frac{1}{2}\sum\limits_{i = 1}^n {\mathit{\boldsymbol{\tilde r}}_i^{\rm{T}}} {m_{{f_i}}}{{\mathit{\boldsymbol{\tilde r}}}_i} + \\ \;\;\;\;\;\;\;\;\;\;\frac{1}{2}\sum\limits_{i = 1}^n {\frac{1}{{{\xi _i}}}} \tilde m_{{f_i}}^2 + \frac{\beta }{2}\int_{t - {T_{ij}}}^t {\mathit{\boldsymbol{\tilde r}}_i^{\rm{T}}} (\rho )(\mathit{\boldsymbol{ \boldsymbol{\varXi} }} \otimes \\ \;\;\;\;\;\;\;\;\;\;\left. {{\mathit{\boldsymbol{I}}_3}} \right){{\mathit{\boldsymbol{\tilde r}}}_i}(\rho ){\rm{d}}\rho + \sum\limits_{i = 1}^n {\sum\limits_{j > i} {{\mathit{\boldsymbol{V}}_{ij}}} } \left( {\left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\|} \right). \end{array} $

(17)

将式(17)对时间求导，可以得到:

$ \begin{array}{l} \dot V(t) = \frac{1}{2}\mathit{\boldsymbol{\dot {\hat g}}}_i^{\rm{T}}{\left( {\mathit{\boldsymbol{ \boldsymbol{\varXi} }} \otimes {\mathit{\boldsymbol{I}}_3}} \right)^{ - 1}}\left( {{\mathit{\boldsymbol{I}}_n} \otimes \mathit{\boldsymbol{ \boldsymbol{\varLambda} }}} \right){{\mathit{\boldsymbol{\ddot {\hat g}}}}_i} + \\ \;\;\;\;\;\;\;\;\;\;\frac{1}{2}\mathit{\boldsymbol{\ddot {\hat g}}}_i^{\rm{T}}{\left( {\mathit{\boldsymbol{ \boldsymbol{\varXi} }} \otimes {\mathit{\boldsymbol{I}}_3}} \right)^{ - 1}}\left( {{\mathit{\boldsymbol{I}}_n} \otimes \mathit{\boldsymbol{ \boldsymbol{\varLambda} }}} \right){{\mathit{\boldsymbol{\dot {\hat g}}}}_i} + \\ \;\;\;\;\;\;\;\;\;\;\sum\limits_{i = 1}^n {\mathit{\boldsymbol{\dot {\tilde r}}}_i^{\rm{T}}} {m_{{f_i}}}{{\mathit{\boldsymbol{\ddot {\tilde r}}}}_i} + \sum\limits_{i = 1}^n {\frac{1}{{{\xi _i}}}} {{\tilde m}_{{f_i}}}{{\dot {\tilde m}}_{{f_i}}} + \\ \;\;\;\;\;\;\;\;\;\;\;\beta \left( {\mathit{\boldsymbol{ \boldsymbol{\varXi} }} \otimes {\mathit{\boldsymbol{I}}_3}} \right)\left[ {\mathit{\boldsymbol{\tilde r}}_i^{\rm{T}}{{\mathit{\boldsymbol{\dot {\tilde r}}}}_i} - \mathit{\boldsymbol{\tilde r}}_i^{\rm{T}}\left( {t - {T_{ij}}} \right){{\mathit{\boldsymbol{\dot {\tilde r}}}}_i}\left( {t - {T_{ij}}} \right)} \right] + \\ \;\;\;\;\;\;\;\;\;\;\;\sum\limits_{i = 1}^n {\sum\limits_{j > i} {{{\left( {{\nabla _{{r_i}}}{\mathit{\boldsymbol{V}}_{ij}}} \right)}^{\rm{T}}}} } {{\mathit{\boldsymbol{\dot r}}}_i} + {\left( {{\nabla _{{r_j}}}{\mathit{\boldsymbol{V}}_{ij}}} \right)^{\rm{T}}}{{\mathit{\boldsymbol{\dot r}}}_j}. \end{array} $

(18)

在对式(18)推导过程中，需要分析如下几点：

1) 根据式(13)可以得到如下:

$ \mathit{\boldsymbol{\ddot {\hat g}}} = \left( {{\mathit{\boldsymbol{I}}_n} \otimes \mathit{\boldsymbol{ \boldsymbol{\varLambda} }}} \right){{\mathit{\boldsymbol{\dot {\hat g}}}}_i} + \left( {\mathit{\boldsymbol{ \boldsymbol{\varXi} }} \otimes {\mathit{\boldsymbol{I}}_3}} \right){{\mathit{\boldsymbol{\dot {\tilde r}}}}_i}. $

(19)

2) 如果编队航天器之间相互作用的通信拓扑结构是无向的，且拓扑图的边指向航天器，则矩阵Ξ =L_A+diag(a₁₀, …, a_n0)是对阵正定的.

3) 根据引理1可知，等式(Ξ⊗I₃)^-1(I_n⊗Λ)= Ξ^-1⊗Λ是对称正定的.

根据上述分析，式(18)可以重新写为

$ \begin{array}{l} \dot V(t) = - \sum\limits_{i = 1}^n {\mathit{\boldsymbol{\dot {\tilde r}}}_i^{\rm{T}}} {m_{{f_i}}}{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i} + \mathit{\boldsymbol{\dot r}}_i^{\rm{T}}{\mathit{\boldsymbol{h}}_i} - \sum\limits_{i = 1}^n {\mathit{\boldsymbol{\tilde r}}_i^{\rm{T}}} {m_{{f_i}}}{\mathit{\boldsymbol{ \boldsymbol{\varPhi} }}_i} + \\ \;\;\;\;\;\;\;\;\;\;\frac{1}{2}\mathit{\boldsymbol{\dot {\hat g}}}_i^{\rm{T}}\left( {\mathit{\boldsymbol{ \boldsymbol{\varXi} }} \otimes {\mathit{\boldsymbol{I}}_3}} \right)\left( {{\mathit{\boldsymbol{I}}_n} \otimes \mathit{\boldsymbol{ \boldsymbol{\varLambda} \boldsymbol{\varTheta} }}} \right){{\mathit{\boldsymbol{\dot {\hat g}}}}_i} + \\ \;\;\;\;\;\;\;\;\;\;\frac{1}{2}\mathit{\boldsymbol{\dot {\hat g}}}_i^{\rm{T}}\left( {{\mathit{\boldsymbol{ \boldsymbol{\varXi} }}^{ - 1}} \otimes {\mathit{\boldsymbol{I}}_3}} \right)\left( {{\mathit{\boldsymbol{I}}_n} \otimes {\mathit{\boldsymbol{ \boldsymbol{\varTheta} }}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varLambda} }}} \right){{\mathit{\boldsymbol{\dot {\hat g}}}}_i} + \\ \;\;\;\;\;\;\;\;\;\;{{\mathit{\boldsymbol{\dot {\hat g}}}}_i}\left( {{\mathit{\boldsymbol{ \boldsymbol{\varXi} }}^{ - 1}} \otimes {\mathit{\boldsymbol{I}}_3}} \right)\left( {{\mathit{\boldsymbol{I}}_n} \otimes \mathit{\boldsymbol{ \boldsymbol{\varLambda} }}} \right)\left( {\mathit{\boldsymbol{ \boldsymbol{\varXi} }} \otimes {\mathit{\boldsymbol{I}}_r}} \right){{\mathit{\boldsymbol{\dot {\tilde r}}}}_i} + \\ \;\;\;\;\;\;\;\;\;\;\sum\limits_{i = 1}^n {\sum\limits_{j > i} {{{\left( {{\nabla _{{r_j}}}{\mathit{\boldsymbol{V}}_{ij}}} \right)}^{\rm{T}}}} } {{\mathit{\boldsymbol{\dot r}}}_i} + {\left( {{\nabla _{{r_j}}}{\mathit{\boldsymbol{V}}_{ij}}} \right)^{\rm{T}}}{{\mathit{\boldsymbol{\dot r}}}_j} - \\ \;\;\;\;\;\;\;\;\;\;\sum\limits_{i = 1}^n {\sum\limits_{j = i}^n {\mathit{\boldsymbol{r}}_i^{\rm{T}}} } {\nabla _{{r_i}}}{\mathit{\boldsymbol{V}}_{ij}} = \mathit{\boldsymbol{\dot {\tilde r}}}_i^{\rm{T}}\left( {{\mathit{\boldsymbol{I}}_n} \otimes \mathit{\boldsymbol{ \boldsymbol{\varLambda} }}} \right){{\mathit{\boldsymbol{\dot {\hat g}}}}_i} + \\ \;\;\;\;\;\;\;\;\;\;\mathit{\boldsymbol{\dot {\hat g}}}_i^{\rm{T}}\left( {{\mathit{\boldsymbol{I}}_n} \otimes \mathit{\boldsymbol{ \boldsymbol{\varLambda} }}} \right){{\mathit{\boldsymbol{\dot {\tilde r}}}}_i} + \\ \;\;\;\;\;\;\;\;\;\;\frac{1}{2}\mathit{\boldsymbol{\dot {\hat g}}}_i^{\rm{T}}\left[ {{\mathit{\boldsymbol{ \boldsymbol{\varXi} }}^{ - 1}} \otimes \left( {{\mathit{\boldsymbol{ \boldsymbol{\varTheta} }}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + \mathit{\boldsymbol{ \boldsymbol{\varLambda} \boldsymbol{\varTheta} }}} \right)} \right]{{\mathit{\boldsymbol{\dot {\hat g}}}}_i} = \\ \;\;\;\;\;\;\;\;\;\; - \frac{1}{2}\mathit{\boldsymbol{\dot {\hat g}}}_i^{\rm{T}}\left( {{\mathit{\boldsymbol{ \boldsymbol{\varXi} }}^{ - 1}} \otimes \mathit{\boldsymbol{Q}}} \right){{\mathit{\boldsymbol{\dot {\hat g}}}}_i} \le 0. \end{array} $

由于Ξ^-1∈R^3×3和Q∈R^3×3均为对称正定矩阵，通过引理1可知Ξ^-1⊗Q也为对称正定矩阵.由于V(t)＞0, $ \dot{V}(t) \leqslant 0 $，${\mathit{\boldsymbol{\dot {\hat g}}}_i}, {\mathit{\boldsymbol{\tilde r}}_i}, {\mathit{\boldsymbol{\dot {\tilde r}}}_i}$与${\tilde m_{fi}}$有界，则式(19)中的$\mathit{\boldsymbol{\ddot {\hat g}}}$也为有界.根据Lyapunov函数分析可知，V(t)有界，$\dot{V}(t)$负半定且连续，则t→∞时，$\mathit{\boldsymbol{\dot {\hat g}}} \to 0$.根据上述分析可知，当t→∞时，$\tilde{\boldsymbol{r}}_{i} \rightarrow \tilde{\boldsymbol{r}}_{j}$, $ {\mathit{\boldsymbol{\dot {\tilde r}}}_i} \to {\mathit{\boldsymbol{\dot {\tilde r}}}_j} \to 0 $，编队航天器飞行系统在无速度测量的情况下可以避免发生碰撞.

3 仿真分析

为了验证本文所设计控制律的有效性，以3个航天器编队飞行系统为例，分别对全状态反馈和无速度测量自适应协同控制进行相应的仿真.

编队航天器在仿真初始设定的轨道参数分别如下：半长轴a=125 km，轨道偏心率e=0.07，轨道倾角i=30°，升交点赤经Ω=45°，近地点幅角ω=30°，真近点角为-2.5×10^-4、-5.0×10^-4、2.5×10^-4；仿真初始时刻参数设定如下，参考航天器真近点角0°，航天器质量分别为：110、115、105 kg，航天器质量估计初始值${\hat m_{\rm{i}}}(0) = 100\;{\rm{kg}}$，通信时延T_ij=1.5 s，相对位置分别为：$ \boldsymbol{r}_{1}(0)=\left[\begin{array}{lll} 50 \sqrt{2} & 0 & 50 \end{array}\right]^{\mathrm{T}} $、$\boldsymbol{r}_{2}(0)=\left[\begin{array}{lll} 0 & -100 & 0 \end{array}\right]^{\mathrm{T}}$、$\boldsymbol{r}_{3}(0)=\left[\begin{array}{lll} -50 \sqrt{2} & 0 & -50 \end{array}\right]^{\mathrm{T}}$，m; 相对速度分别为：$ {\mathit{\boldsymbol{\dot r}}_1}(0) = {\left[ {\begin{array}{*{20}{l}} 0&{ - 0.5}&0 \end{array}} \right]^{\rm{T}}} $、${\mathit{\boldsymbol{\dot r}}_2}(0) = {\left[ {\begin{array}{*{20}{l}} { - 0.25}&0&{ - 0.25} \end{array}} \right]^{\rm{T}}}$、$\dot{\boldsymbol{r}}_{3}(0)=\left[\begin{array}{llll} 0 & -0.5 & 0 \end{array}\right]^{\mathrm{T}}$，m/s; 期望位置为$ \boldsymbol{r}_{i}^{d}=\left[\begin{array}{c} 100 \sqrt{2} \cos (n t+\pi / 2 \cdot i) \\ -200 \sin (n t+\pi / 2 \cdot i) \\ 100 \cos (n t+\pi / 2 \cdot i) \end{array}\right] $，m.

将避免碰撞的最小安全距离设置为c_ij=10 m，选取势函数为

$ {\mathit{\boldsymbol{V}}_{ij}}\left( {{\mathit{\boldsymbol{r}}_i},{\mathit{\boldsymbol{r}}_j}} \right) = \sum\limits_{i = 1}^3 {\frac{{211000}}{{{{\left( {{{\left\| {{\mathit{\boldsymbol{r}}_i} - {\mathit{\boldsymbol{r}}_j}} \right\|}^2} - {{20}^2}} \right)}^2}}}} . $

首先，设置全状态反馈自适应协同控制律(14)~(16)中的参数如下：λ_i=50, γ_i=50, α=1，通信拓扑图为无向图，则a_ij=0.5，反之a_ij=0.

在全状态反馈自适应协同避碰控制律(8)~(10)的作用下，图 2、3分别为航天器编队相对位置误差曲线和相对速度误差曲线.由图 2、3可以看出，在跟踪到期望位置的过程中，编队航天器相对于参考航天器速度收敛于零附近，同时编队航天器相对距离在100 s后基本保持不变，在任意时刻相对距离均不小于c=10 m，避免碰撞的发生.

其次，设置无速度测量自适应协同避碰控制律(15)~(16)中的参数如下：β=1, ξ_i=0.1通信拓扑图为无向图，则a_ij=0.3，反之a_ij=0.

在无速度测量自适应协同避碰控制律(15)~(16)的作用下，图 4、5分别为航天器编队相对位置误差曲线和相对速度误差曲线.由图 4、5可以看出，位置跟踪误差曲线和速度跟踪误差曲线均能收敛到零附近，且收敛时间具有同步性并保持编队飞行系统稳定，在图 4中，相对位置跟踪误差以稳态误差约为10^-3 m，而在图 5中相对速度跟踪的稳态误差约为10^-5 m/s.

图 6、图 7分别给出了全状态反馈和无速度测量控制力曲线图.在图 6、图 7中可以看出，在仿真的初始阶段由于存在通信延迟会产生较大的抖振，但随着时间的变化曲线均达到收敛.与全状态反馈自适应协同避碰控制律相比，编队航天器在无速度测量自适应协同避碰控制律作用下，相对振幅较大，是由于航天器质量估计值在自适应更新律中稳态控制不能为零.因此，根据设计的控制律(8)~(10)、(15)~(16)可以看出，当辅助变量收敛到零附近时，控制变量不能达到零.

编队航天器避免碰撞的相对位置变化轨迹如图 8所示，在航天器编队飞行时假设各航天器为刚体结构，但是在实际应用运行过程中会受到姿态展开等因素的影响，这些干扰因素仍然是无法忽略的.由图 8可以看出，各航天器均能到达期望轨迹且不发生碰撞.根据上述分析可知，航天器编队飞行闭环系统稳定性，并完成编队跟踪、避免碰撞、自适应协同控制的控制目标.

4 结论

1) 在航天器编队飞行控制系统存在通信时延、质量参数不确定的情况下，考虑反馈状态信息的不同状态，提出了两种具有较强鲁棒性的自适应协同避碰控制律.

2) 选取合理的势函数包含在自适应协同避碰控制律中，实现编队飞行系统避免碰撞的发生，具有较高的控制精度.

3) 利用Lyapunov稳定性理论对所设计的两种控制律进行理论分析证明，证明航天器编队飞行闭环系统的渐近稳定性.仿真结果进一步验证了所设计控制律的有效性.