Децентрализованное управление движением роя малых космических аппаратов для удержания коммуникационной связности

Полный текст

ВВЕДЕНИЕ

В настоящее время бурно развиваются распределённые космические системы, состоящие из большого числа космических аппаратов и решающие единую задачу. Самым распространённым типом таких систем является созвездие, или группировка, когда космические аппараты равномерно распределены на одной или нескольких орбитах. Управление орбитальным движением в таких системах, как правило, осуществляется индивидуально для каждого космического аппарата с помощью команд из центра управления полётом согласно серии манёвров, как, например, в работе [1]. Другой тип распределённых космических систем — это групповой полёт аппаратов на близком расстоянии друг относительно друга. Управление относительным движением такой системы осуществляется с помощью автономных алгоритмов управления, реализованных на бортовых компьютерах космических аппаратов, на основе измерений параметров относительного движения. Относительное движение может оцениваться с помощью специальных систем, например на основе обработки видеоизображения [2–5] или измерений лазерных дальномеров [6], а может и вычисляться на основе переданной по межспутниковому каналу связи информации об орбитальном положении соседних аппаратов, полученной, например, с помощью автономной навигационной системы, как в статье [7]. Однако каждая автономная система определения относительного движения имеет ограниченный радиус работы, равно как и уверенный приём по межспутниковому каналу связи осуществляется лишь в пределах заданного относительного расстояния, определяемого аппаратными особенностями. Таким образом, при построении алгоритмов управления относительным движением прежде всего требуется учитывать коммуникационные ограничения и обеспечивать связность группового полёта космических аппаратов.

Групповой полёт большого числа космических аппаратов называется роем [8]. Применение роя малых аппаратов рассматривается в задачах изучения магнитосферы Земли [9], сборки орбитальных станций [10], построения распределённой системы сенсоров в ионосфере [11], исследования параметров модели гравитационного поля астероидов [12]. Основная особенность управления движением роя космических аппаратов заключается в невозможности учесть в алгоритме управления относительное движение всех аппаратов в группе, в частности из-за коммуникационных ограничений. Поэтому управление строится каждым аппаратом децентрализованно, т.е. независимо от управления других аппаратов, и, как правило, на основе информации об относительных траекториях только ближайших соседей. Это может приводить к случайным хаотическим относительным траекториям, при этом основная задача управления сводится к обеспечению связности роя. В работе [13] был рассмотрен ряд алгоритмов управления с обратной связью с учётом коммуникационных ограничений. Для управления групповым полётом в условиях неопределённостей в статье [14] предлагается использование нейронных адаптивных сетей.

В предыдущей работе авторов [15] был предложен подход к устранению относительного дрейфа траектории с помощью управления на основе расчёта среднего дрейфа аппаратов внутри коммуникационной области. Управляемое относительное движение было исследовано аналитически с использованием теории графов в работе [16]. Однако устранение относительного дрейфа не всегда приводит к обеспечению связности роя, как показано в статье [17], поэтому в настоящей работе рассматривается управление с учётом среднего относительного сдвига траектории, что позволяет задать ограничения на относительное положение аппаратов.

ПОСТАНОВКА ЗАДАЧИ И УРАВНЕНИЯ ДВИЖЕНИЯ

Рассматривается группа малых космических аппаратов (МКА), запущенных на низкую околоземную орбиту. Предполагается, что каждый МКА оснащён либо системой определения относительного движения, либо относительное движение вычисляется на основе измерений автономной навигационной системы, переданных по каналу межспутниковой связи. Эта система позволяет оценивать положение и скорость любого соседнего аппарата в пределах максимальной дальности, которую мы в дальнейшем будем называть зоной видимости или областью коммуникации. Если расстояние до какого-то МКА больше, чем область коммуникации, то информация о его движении отсутствует. В начальный момент времени группа МКА движется в соответствии с заданными условиями отделения, которые зависят от характеристик транспортно-пусковых контейнеров. Такая система запуска имеет определённые погрешности исполнения, которые становятся причиной различных скоростей у МКА после отделения. При отсутствии управления под действием гравитационного поля Земли различия в начальных скоростях приведут к постепенному увеличению относительных расстояний между МКА. Поэтому необходимо разработать алгоритм децентрализованного управления, который будет препятствовать увеличению относительных расстояний. Требуемое управление будет создаваться за счёт разности сил аэродинамического сопротивления, действующих на аппараты. Изменение величины аэродинамического сопротивления, действующего на спутник, достигается корректировкой его ориентации, т.е. изменением площади поперечного сечения относительно набегающего потока. Предполагается, что спутники оснащены трёхосной системой ориентации для управления угловым движением.

Для описания траекторий МКА удобно пользоваться уравнениями движения, записанными в относительных координатах. Общий вид уравнений относительного движения двух произвольных МКА из группы достаточно сложен для аналитического рассмотрения, поэтому в работе используется линеаризованная модель движения, описываемая системой уравнений Хилла – Клохесси — Уилтшира [18, 19]. Модель описывает движение произвольно выбранного МКА из группы относительно орбитальной системы координат (ОрбСК) под действием центрального гравитационного поля Земли. Начало отсчёта этой системы находится в точке O, движущейся по круговой орбите радиусом , с орбитальной угловой скоростью . Ось Oz направлена вдоль прямой, соединяющей центр масс Земли с точкой O, ось Oy направлена по нормали к плоскости орбиты, ось Ox дополняет тройку до правой (рис. 1).

 

Рис. 1. Орбитальная система координат.

 

Пусть ${\mathbf{r}}_i={\left[x_i,y_i,z_i\right]}^T$ — координаты i-го МКА в ОрбСК, $i=\mathrm{1,...,}N,$ где N — количество МКА в группе. Тогда для координат вектора относительного положения МКА линеаризованная система уравнений движения имеет вид

$\begin{array}{l}{\ddot{x}}_i+2\omega {\dot{z}}_i=\mathrm{0,}\\ {\ddot{y}}_i+{\omega }^2{y}_i=\mathrm{0,}\\ {\ddot{z}}_i-2\omega {\dot{x}}_i-3{\omega }^2{z}_i=0.\end{array}$   (1)

Уравнения справедливы для малых относительных расстояний, поэтому расстояние между i-м МКА и точкой O должно быть на несколько порядков меньше r₀. Система линейных дифференциальных уравнений (1), описывающая движение i-го МКА, имеет следующее решение:

$\begin{array}{l}{x}_i=D_i\left(t\right)+2A_i\cos \left(\omega t+{\psi }_i\right),\\ y_i=B_i\cos \left(\omega t+{\phi }_i\right),\\ z_i=2C_i+A_i\sin \left(\omega t+{\psi }_i\right),\end{array}$   (2)

где $A,B,C,D,\psi ,\phi$ — параметры движения, которые зависят от начальных условий ${\left[x_i\left(0\right),y_i^{}\left(0\right),z_i^{}\left(0\right),{\dot{x}}_i\left(0\right),{\dot{y}}_i^{}\left(0\right),{\dot{z}}_i^{}\left(0\right)\right]}^T$:

$\begin{array}{l}{C}_i=\frac{{\dot{x}}_i\left(0\right)}{\omega }+2z_i\left(0\right),\\ D_i\left(t\right)=-3C_i\omega t+x_i\left(0\right)-\frac{2{\dot{z}}_i\left(0\right)}{\omega },\\ B_i=\sqrt{\frac{{\dot{y}}_i^2\left(0\right)}{{\omega }^2}+y_i^2\left(0\right)},\\ A_i=\frac{1}{\omega }\sqrt{{\dot{z}}_i^2\left(0\right)+4z_i^2\left(0\right){\omega }^2+8z_i\left(0\right){\dot{x}}_i\left(0\right)\omega +4{\dot{x}}_i^2\left(0\right)}.\end{array}$   (3)

Можно заметить, что движение по оси Oy является ограниченным, поэтому в работе будет рассмотрено только движение спутников в плоскости Oxy. По оси Ox относительное расстояние линейно увеличивается со временем пропорционально параметру C_i (рис. 2). Далее этот параметр будем называть параметром дрейфа МКА относительно точки O. Таким образом, чтобы относительные расстояния между МКА в группе не увеличивались со временем, необходимо, чтобы выполнялось условие $C_i-C_j=C_{ij}=0$ для любых двух аппаратов i и j. Приведение к нулю всех относительных параметров дрейфа также приведёт к достижению замкнутых траекторий для всех аппаратов относительно друг друга, поэтому основной задачей управления относительным движением будет устранение относительного параметра дрейфа C_ij между всеми МКА в группе. Однако устранение дрейфа может произойти уже вне области коммуникации, поэтому другим важным параметром движения для поддержания связности группы является относительный сдвиг $D_{ij}=D_i-D_j$, который должен быть меньше радиуса связи.

 

Рис. 2. Пример траектории спутника при ненулевом параметре C_i.

 

Рассмотрим управляемое движение группы МКА. Будем предполагать, что каждый аппарат оснащён системой управления относительным движением. Введём обозначение для вектора ускорения , где u_i — ускорение, реализованное i-м МКА относительно точки O. Тогда параметры движения с учётом управления u_i будут меняться следующим образом [20]:

${\mathbf{u}}_i={\left[u_i^x,u_i^y,u_i^z\right]}^T$, где u_i — ускорение, реализованное i-м МКА относительно точки O. Тогда параметры движения с учётом управления u_i будут меняться следующим образом [20]:

$\begin{array}{l}{\dot{A}}_i=\frac{1}{\omega }\left(u_i^z\cos {\psi }_i-2u_i^x\sin {\psi }_i\right),\\ {\dot{B}}_i=-\frac{1}{\omega }{u}_i^y\sin {\phi }_i,\\ {\dot{C}}_i^{}=\frac{1}{\omega }{u}_i^x,\\ {\dot{D}}_i=-3\omega C_i-\frac{2}{\omega }{u}_i^z,\\ {\dot{\psi }}_i=-\frac{1}{\omega A_i}\left(u_i^z\sin {\psi }_i+2u_i^x\cos {\psi }_i\right),\\ {\dot{\phi }}_i=-\frac{1}{\omega B_i}{u}_i^y\sin {\phi }_i.\end{array}$  (4)

Так как на величину C_i влияет только проекция ускорения на ось Ox, то далее будет рассматриваться управление вида ${\mathbf{u}}_i={\left[u_i^x\mathrm{,0,0}\right]}^T$. Далее такой подход позволит использовать силу аэродинамического сопротивления в качестве управляющего воздействия.

В настоящей работе предполагается, что начальные условия отделения аппаратов таковы, что параметры относительных траекторий A_ij и B_ij малы по сравнению с размерами зоны видимости. Это соответствует малой разнице в наклонениях орбит аппаратов и малым значениям эксцентриситетов околокруговых траекторий. В случае больших значений амплитуд A_ij и B_ij возможен вариант, когда даже при малых значениях параметров C_ij и D_ij аппараты будут вне зоны видимости. Управление движением для обеспечения коммуникационной связности в этом случае не может быть реализовано только с помощью разницы сил аэродинамического сопротивления и будет рассмотрено авторами в последующих работах.

АЛГОРИТМ УПРАВЛЕНИЯ

Так как основной целью управления является устранение относительного параметра дрейфа C_ij, рассмотрим управление по оси Ox, построенное на основе средних значений параметров дрейфов аппаратов внутри группы. Для i-го МКА известно относительное движение не всех аппаратов в группе, а только тех, которые попадают в его зону видимости. Зона видимости представляет собой сферу радиуса R_comm с центром в центре масс i-го МКА. Количество аппаратов, попавших в такую зону, обозначается $N_{comm}^i$ (рис. 3).

 

Рис. 3. Зона видимости i-го МКА.

 

Построим алгоритм управления на основе прямого метода Ляпунова сначала для двух аппаратов. Рассмотрим кандидат-функцию Ляпунова в следующем виде:

$V=\frac{1}{2}\left(C_{ij}^2+D_{ij}^2\right).$

Полная производная этой функции по времени в силу уравнений (4) и управления только вдоль оси Ox имеет следующий вид:

$\dot{V}=C_{ij}{\dot{C}}_{ij}+D_{ij}{\dot{D}}_{ij}=C_{ij}\left(\frac{1}{\omega }{u}_{ij}^x\right)+D_{ij}\left(-3\omega C_{ij}\right)$

где $u_{ij}^x=u_i^x-u_j^x$ — разница управляющих ускорений, действующих на спутники. Для обеспечения асимптотической устойчивости, согласно теореме Барбашина —Красовского [21], производная функции Ляпунова должна быть неотрицательна, а область, где она равна нулю, не должна содержать целых траекторий. Это достигается при управлении следующего вида:

$u_{ij}^x=-k{C}_{ij}+3{\omega }^2{D}_{ij},$   (5)

где $k>0$ — параметр управления. Управление (5) стремится к устранению относительного параметра дрейфа C_ij, а также к достижению нулевого относительного сдвига D_ij мгновенного центра относительно эллиптической траектории. При реализации децентрализованного управления $u_{ij}^x$ предполагается, что каждый аппарат действует независимо от других аппаратов и прикладывает управление $u_i^x$, считая $u_j^x=0$.

Теперь рассмотрим случай, когда внутри области видимости i-го спутника находится $N_{comm}^i$ аппаратов. Для решения задачи устранения всех относительных дрейфов и удержания их в области коммуникации предлагается следующий децентрализованный закон управления, схожий с законом (5):

$u_i^x=-k{\overline{C}}_i+3{\omega }^2{\overline{D}}_i,$   (6)

где ${\overline{C}}_i$ — это средний параметр дрейфа, а ${\overline{D}}_i$ — средний сдвиг траекторий, которые вычисляются следующим образом:

${\overline{C}}_i=\sum _{j=1}^{N_{comm}^i}{C}_{ij}/N_{comm}^i$, ${\overline{D}}_i=\sum _{j=1}^{N_{comm}^i}{D}_{ij}/N_{comm}^i.$   (7)

АНАЛИТИЧЕСКОЕ ИССЛЕДОВАНИЕ УПРАВЛЯЕМОГО ДВИЖЕНИЯ

В разделе проведено аналитическое исследование закона управления относительным движением.

Исследуем характеристики движения группы спутников под действием предложенного управления (6). Подставляя его в уравнения (4), описывающие изменения параметров D_i и C_i, можно получить следующие линейные дифференциальные уравнения, записанные в векторном виде:

$\left(\begin{array}{c}\dot{C}\\ \dot{D}\end{array}\right)=3\omega \left(\begin{array}{cc}aS& S\\ -E& 0\end{array}\right)\cdot \left(\begin{array}{c}C\\ D\end{array}\right),$   (8)

где S — матрица размера $N\times N$ с элементами

$s_{ij}=\left\{\begin{array}{l}\mathrm{1,}\text{\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}если }i=j,\\ \mathrm{0,}\text{\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}если }i\ne j\text{ и }i\text{-й спутник не видит }j\text{-й,}\\ -1/N_{comm}^i,\text{\hspace{0.17em}\hspace{0.17em}если }i\ne j\text{ и }i\text{-й спутник видит }j\text{-й,}\end{array}\right.$   (9)

$C={\left[C_1,C_2,\mathrm{...},C_N\right]}^T,D={\left[D_1,D_2,\mathrm{...},D_N\right]}^T$ — векторы параметров дрейфов и сдвигов МКА относительно точки O соответственно; N — количество аппаратов в группе; E — единичная матрица размером $N\times N$; параметр $a=-\frac{k}{3{\omega }^2}$.

Рассмотрим следующую матрицу:

$T=\left(\begin{array}{cc}aS& S\\ -{\rm E}& 0\end{array}\right).$   (10)

Собственные и присоединённые векторы и собственные значения матрицы полностью определяют решение уравнений (8). Поскольку T — блочная матрица, будем искать её собственные векторы в блочном виде. Пусть ${\left(\begin{array}{cc}v& w\end{array}\right)}^T$ является собственным вектором матрицы T (здесь v и w — N-мерные векторы); пусть  является соответствующим собственным значением, тогда

$T\cdot \left(\begin{array}{c}v\\ w\end{array}\right)=\mu \cdot \left(\begin{array}{c}v\\ w\end{array}\right).$

Блочная форма матрицы T (10) позволяет привести это к системе

$\left\{\begin{array}{l}aSv+Sw=\mu v,\\ -v=\mu w.\end{array}\right.$   (11)

Второе уравнение также означает, что каждый из векторов v, w ненулевой. Подставив второе уравнение в первое, получим квадратное уравнение на $\mu$:

$\left(1-\mu a\right)Sw=-{\mu }^{2}w.$   (12)

Здесь $\mu =1/a$ не даёт ненулевых решений, а если $(1-\mu a)\ne 0$, то w является собственным вектором матрицы S, а значит, для соответствующего ему собственного значения $\lambda$ выполняется

$Sw=\lambda w.$   (13)

Полученное соотношение подставим в (12), тогда

$\left(1-\mu a\right)\lambda w=-{\mu }^{2}w.$

Отсюда можно получить связь между собственными значениями матрицы T и матрицы S:

${\mu }_{\pm }=\frac{a\lambda \pm \sqrt{a^2{\lambda }^2-4\lambda }}{2}.$   (14)

Если $\lambda \ne 0\text{\hspace{0.17em} или }\frac{4}{a^2},$ то каждому собственному значению S соответствует два различных собственных значения T.

Введём некоторые понятия из теории графов для дальнейшего исследования закона управления. Матрица смежности A для графа с числом вершин N — это квадратная матрица $N\times N$ такая, что элемент этой матрицы a_ij равен числу рёбер, исходящих из i-й вершины в j-ю. Матрица степеней графа D_g — матрица, на главной диагонали которой степени вершин графа, а остальные элементы нули. Также будем полагать, что граф является связным, т.е. существует путь из любой вершины в любую, и матрица D_g является невырожденной. Матрица Кирхгофа $\mathcal{L}$ и нормированная матрица Кирхгофа  получаются следующим образом:

$L=D_g-A,\mathcal{L}=D_g^{-1/2}L{D}_g^{-1/2}.$

Можно заметить, что матрица S связана с нормированной матрицей Кирхгофа $\mathcal{L}$ следующим образом:

$S=E-D_g^{-1}A=D_g^{-1/2}\mathcal{L}{D}_g^{1/2}.$   (15)

Это означает, что собственные числа матрицы S совпадают с собственными числами матрицы $\mathcal{L}$. Для рассматриваемого алгоритма управления следующие свойства спектра $\mathcal{L}$ представляют особый интерес. $\mathcal{L}$ — симметричная матрица, так что её собственные числа вещественны. Нуль всегда является собственным числом. Кратность нулевого собственного числа равна количеству компонент связности графа. Компонента связности графа — некоторое множество вершин графа такое, что для любых двух вершин из этого множества существует путь из одной в другую и не существует пути из вершины этого множества в вершину не из этого множества. В нашей постановке задачи компонента связности будет равна количеству подгрупп, на которые разделилась группа МКА. Второе важное свойство нормированной матрицы Кирхгофа — все ненулевые собственные числа строго положительные. Единственное нулевое собственное число в таком случае будет соответствовать величине параметров дрейфа всей группы относительно точки O.

Так как для неориентированных графов матрица $\mathcal{L}$ является симметричной, то собственные векторы $\mathcal{L}$ образуют полный набор векторов в пространстве. Из сказанного выше следует, что собственные векторы матрицы S также образуют полный набор. Рассмотрим все возможные варианты значений собственных чисел.

Случай $\lambda \ne 0\text{\hspace{0.17em} и }\lambda \ne \frac{4}{a^2}$

Пусть $Sw=\lambda w,$ причём $\lambda \ne 0\text{\hspace{0.17em} и }\lambda \ne \frac{4}{a^2},$ в этом случае

$T\left(\begin{array}{c}-{\mu }_{\pm }w\\ w\end{array}\right)={\mu }_{\pm }\left(\begin{array}{c}-{\mu }_{\pm }w\\ w\end{array}\right).$

Следовательно, каждый такой собственный вектор порождает пару собственных векторов матрицы T с собственными значениями ${\mu }_{\pm }$.

Случай $\lambda =0$

Собственное значение $\lambda =0$ матрицы S в этом случае порождает собственное значение  кратности два у матрицы T. Пусть $w_0$ — это собственный вектор матрицы S, соответствующий $\lambda =0$. Из уравнений (13) и (11) при подстановке $\mu =0$ и $\lambda =0$ получим, что ${\left(\begin{array}{cc}0& w_0\end{array}\right)}^T$ является собственным вектором T. Других линейно независимых собственных векторов нет, так что будем искать присоединённый вектор ${\left({\begin{array}{cc}{p}_0& q\end{array}}_0\right)}^T$ в виде

$\left(\begin{array}{cc}aS& S\\ -{\rm E}& 0\end{array}\right)\cdot \left(\begin{array}{c}{p}_0\\ q_0\end{array}\right)=\left(\begin{array}{c}0\\ w_0\end{array}\right).$   (16)

Получим систему уравнений

$\begin{array}{l}aS{p}_0+S{q}_0=\mathrm{0,}\\ -p_0=w_0.\end{array}$

Подставим второе уравнение в первое и, учитывая, что w₀ является собственным вектором матрицы S для $\lambda =0$, получим уравнение для q₀:

$-aS{w}_0+S{q}_0=0\Rightarrow S{q}_0=0.$

Замечая, что q₀ является собственным вектором S для $\lambda =0$, находим присоединённый вектор ${\left(\begin{array}{cc}-w_0& \alpha w_0\end{array}\right)}^T$ при произвольном $\alpha$. Таким образом, для связного графа корневое подпространство матрицы T, отвечающее собственному значению 0, двумерно.

Случай $\lambda =\frac{4}{a^2}$

Рассмотрим теперь случай, когда $\lambda =\frac{4}{a^2},$ что означает, что $S{w}_1=\frac{4}{a^2}{w}_1.$ Из (11) при подстановке согласно $\mu =\frac{2}{a}$ получаем, что ${\left(\begin{array}{cc}-\frac{2}{a}{w}_1& w_1\end{array}\right)}^T$ является собственным вектором для матрицы T. Для присоединённого вектора ${\left({\begin{array}{cc}{p}_1& q\end{array}}_1\right)}^T$ должно выполняться

$\left[\left(\begin{array}{cc}aS& S\\ -{\rm E}& 0\end{array}\right)-\frac{2}{a}\left(\begin{array}{cc}E& 0\\ 0& E\end{array}\right)\right]\cdot \left(\begin{array}{c}{p}_1\\ q_1\end{array}\right)=\left(\begin{array}{c}-\frac{2}{a}{w}_1\\ w_1\end{array}\right).$

Преобразуя, получим систему уравнений

$\begin{array}{l}aS{p}_1+S{q}_1-\frac{2}{a}{p}_1=-\frac{2}{a}{w}_1,\\ -p_1-\frac{2}{a}{q}_1=w_1.\end{array}$   (17)

Выражая q₁ из второго уравнения и подставляя в первое, получаем следующее:

$S{p}_1=\frac{4}{a^2}{p}_1.$

Заметим, что p₁ является собственным вектором, соответствующим $\lambda =\frac{4}{a^2}$. Подставив $p_1=w_1$ в (17) и выразив q₁, получим присоединённый вектор ${\left(\begin{array}{cc}{w}_1& -a{w}_1\end{array}\right)}^T$. Таким образом, каждый линейно независимый собственный вектор матрицы S с собственным значением $\lambda =\frac{4}{a^2}$ порождает двумерное корневое подпространство у матрицы T.

Суммируя все рассмотренные случаи, можно сделать вывод, что N − 1 положительным значениям $\lambda$, будут соответствовать 2N − 2 значения  (с учётом кратности) с отрицательной действительной частью. Значение $\mu =0$ будет иметь кратность равную двум и будет порождать двумерное корневое подпространство у матрицы T. Зная собственные значения матрицы T и соответствующие им собственные и присоединённые векторы, можно выписать общее решение для (8):

$\left(\begin{array}{c}C\\ D\end{array}\right)=b_0\left(\begin{array}{c}0\\ w_0\end{array}\right)+b_1\left(3\omega t\left(\begin{array}{c}0\\ w_0\end{array}\right)+\left(\begin{array}{c}-w_0\\ 0\end{array}\right)\right)+\sum _{i=1}^{2N-1}{p}_i\left(t\right)e^{3\omega {\mu }_{i}t}\left(\begin{array}{c}{v}_i\\ w_i\end{array}\right).$

Здесь p_i — это многочлен, кратность которого зависит от кратности соответствующего собственного числа и количества собственных векторов у него; b₀ и b₁ — константы, зависящие от начальных условий. Поскольку все ненулевые $\mu$ имеют отрицательную вещественную часть, все слагаемые в сумме экспоненциально затухают. Собственный вектор w₀, который удовлетворяет $S{w}_0=0$, для связного графа пропорционален ${\left(\begin{array}{cccc}1& 1& \dots & 1\end{array}\right)}^T$. Следовательно, на больших промежутках времени движение происходит по закону

$\begin{array}{l}{C}_i=-b_1,\\ D_i=b_0+3b_1\omega t.\end{array}$

Это означает, что на больших промежутках времени построенное управление (6) обеспечивает движение группы как целого с общим для каждого из спутников дрейфом, так что относительные сдвиги стремятся к нулю. Таким образом, предложенный децентрализованный закон управления обеспечивает сходимость к нулю всех относительных параметров дрейфов и сдвигов в случае, когда группа образует связный граф. В случае нарушения связности графа система (8) имеет больше одного нулевого собственного числа, что соответствует разделению роя на независимые подгруппы.

ЧИСЛЕННОЕ ИССЛЕДОВАНИЕ УПРАВЛЯЕМОГО ДВИЖЕНИЯ

Проведём математическое моделирование управляемого движения группы спутников под действием предложенного управления. При моделировании движения спутников учитывалось влияние центрального поля Земли. Уравнения движения в инерциальной системе координат (ИСК) имеют следующий вид:

${\ddot{R}}_i=-\mu \frac{R_i}{{R_i}^3}+u_i^{},$

где R_i — радиус-вектор i-го МКА в инерциальной системе координат, начало отсчета которой совпадает с центром Земли, $R_i=\left|R_i\right|$; $\mu$ — гравитационный параметр Земли; u_i — требуемое управление, рассчитанное в ИСК. Так как основной целью численного моделирования является верификация результатов аналитического исследования, в работе предполагается, что требуемое управление реализуется точно, ошибки знания параметров орбиты и ошибки реализации управления не рассматриваются. Так как требуемое управление направлено вдоль трансверсали, оно может быть реализовано с помощью разницы действующих на аппараты аэродинамических сил, но особенности этой реализации (неточность знания плотности атмосферы, ошибки ориентации аппаратов для обеспечения требуемой площади поперечного сечения относительно набегающего потока) в этой работе не рассматриваются и становятся предметом дальнейших исследований.

Для достижения связности группы КА ограничим с помощью управления средний сдвиг $\left|{\overline{D}}_i\right|$ каждого спутника некоторой величиной меньше радиуса видимости $R_{comm}$. Для этого закон управления приложенного вдоль оси Ox ОрбСК предлагается сделать зависимым от среднего сдвига траектории i-го МКА в его зоне видимости:

$u_i^x=\left\{\begin{array}{l}-k{\overline{C}}_i,\left|{\overline{D}}_i\right|<0.85R_{comm},\\ -k{\overline{C}}_i+3{\omega }^2{\overline{D}}_i,0.85R_{comm}\le \left|{\overline{D}}_i\right|\le R_{comm}.\end{array}\right.$   (18)

В этом случае, как только $\left|{\overline{D}}_i\right|$ подходит близко к границе зоны видимости, происходит переключение управления на то, которое уменьшает относительные сдвиги траекторий аппаратов $\left|D_{ij}\right|$, что в свою очередь будет уменьшать и относительные расстояния между аппаратами. Таким образом будет сохраняться связность всей группы МКА. В случае же, когда $\left|{\overline{D}}_i\right|<0.85R_{comm}$, управление стремится только к устранению среднего относительного параметра дрейфа. Закон управления (18) позволяет задать некоторые границы для относительного движения, тогда как закон вида (6) устремляет все относительные сдвиги D_ij к нулю, что будет соответствовать набору относительных замкнутых эллиптических траекторий с единым центром, что в свою очередь повышает вероятность столкновения МКА.

Оценка величины необходимого R_comm для сохранения связности группы при управлении вида $u_i^x=-k{\overline{C}}_i$ представлена в работе авторов [16]. Требуемый размер зоны видимости в случае кластерного запуска N аппаратов вдоль вектора трансверсали с нормально распределёнными случайными скоростями отделения (математическое ожидание ${\mu }_V$, дисперсия скорости отделения ${\sigma }_V$) определяется следующим образом:

$R_{comm}={\mu }_D+3{\sigma }_D=3\Delta t{\mu }_V+3{\sigma }_V\sqrt{9\Delta t^2(2N^2-2N+1)+\frac{8}{{\omega }^2}+\frac{18{\omega }^2}{k^2{\lambda }_{\min }^2}}$,   (19)

где ${\lambda }_{\min }^{}$ — наименьшее положительное собственное число матрицы S, определяемой в (9); $\Delta t$ — время между запусками; ${\sigma }_D$ — среднеквадратический сдвиг относительных траекторий при свободном относительном движении после кластерного запуска с заданными случайными ошибками отделения. Все параметры, используемые для моделирования управляемого движения группы МКА, представлены в таблице. Минимальное собственное значение ${\lambda }_{\min }$ матрицы S было посчитано для полного графа с 20 вершинами. Интервал между пересчётом управления означает, что каждые 300 с данные о положениях МКА обновляются. По этим данным рассчитываются средние относительные параметры дрейфа и сдвиги, и получаются новые величины требуемого управления для каждого МКА. Так как движение вдоль оси Oy ОрбСК является ограниченным согласно (2), траектории движения МКА представлены в проекции на плоскость Oxz.

 

Параметры моделирования

Количество спутников, N	20
Масса спутника, m	3 кг
Математическое ожидание скорости отделения, ${\mu }_V$	0.05 м/с
Дисперсия скорости отделения, ${\sigma }_V$	0.01 м/с
Минимальное собственное значение, ${\lambda }_{\min }$	−1.053
Интервал времени между запусками, $\Delta t$	3 с
Высота орбиты, h	350 км
Наклонение орбиты, i	${51.7}^{\text{o}}$
Интервал времени между пересчётом управления	300 с
Коэффициент управления, k	$1.85\cdot {10}^{-7}$

 

Пример успешной работы алгоритма управления

Для начальных условий, представленных в табл. 1, необходимый R_comm, согласно (19), составляет 802 м. На рис. 4 представлены траектории движения всех спутников относительно первого запущенного аппарата, которые, как видно, постепенно становятся замкнутыми, и таким образом аппараты образуют связную группу. Изменение величин относительных параметров дрейфа $C_{1j},j=\mathrm{2,...,}N,$ вычисленных также относительно первого запущенного аппарата, представлено на рис. 5. Для рассматриваемого примера относительные параметры дрейфов спутников сходятся к нулю, что ещё раз подтверждает, что относительные траектории спутников в группе становятся ограниченными. Также можно заметить, что устранение относительного параметра дрейфа происходит за первые 15 ч после запуска. После этого величины параметров дрейфов по модулю не выходят за границы 0.01 м.

 

Рис. 4. Относительные траектории движения МКА при реализации управления.

 

Рис. 5. Величины относительных параметров дрейфов.

 

На рис. 6 представлен график приложенного ко всем спутникам в рое управления в зависимости от времени. Предполагается, что на все аппараты действуют силы аэродинамического сопротивления, поэтому величина управления вдоль оси Ox ОрбСК отрицательна. Коэффициент управления k подобран таким образом, чтобы величина управления не превышала ускорение, которое возможно получить с помощью силы аэродинамического сопротивления путём поворота 3U кубсата (параллелепипед с размерами 10 × 10 × 34 см) относительно набегающего потока. При этом учитывались характерная плотность воздуха на данной высоте орбиты и возможные средние значения параметров траекторий C и D при заданных условиях отделения. Если для некоторых аппаратов в начальный момент времени требуемое управляющее ускорение превышало максимально допустимое значение, в численном моделировании реализовывалось это максимальное значение управляющего ускорения. Отметим также, что по мере уменьшения относительных параметров дрейфа разница приложенного управления уменьшается. Резкие изменения величины приложенного управления соответствуют моментам времени, когда количество аппаратов в зоне видимости изменяется, или моментам переключения между законами управления (18). Изменение $D_{1j},j=\mathrm{2,...,}N,$ во время моделирования представлено на рис. 7. 

Значения параметров сходятся к константным значениям, так как относительные параметров дрейфа сходятся к нулю, что соответствует уравнениям (4).

 

Рис. 6. Реализованная величина управления.

 

Рис. 7. Величины относительных сдвигов.

 

На рис. 8 представлено изменение степеней вершин графа, который образуют МКА в течение моделирования. В течение первых 2.25 ч степень всех вершин в графе равна 19, а значит, МКА образуют полный граф. Далее степень вершин варьируется от 8 до 19, что означает, что некоторые аппараты всё ещё могут отслеживать движение всех спутников в группе, а некоторые только какой-то подгруппы спутников. Несмотря на это, МКА образуют связный граф на протяжении всего моделирования. Так как вид графа меняется во время моделирования, то меняется и скорость сходимости. На рис. 9 показана величина $-\frac{k}{\omega }{\lambda }_{\min }$, определяющая скорость сходимости, что является наибольшим собственным значением, за исключением нуля для матрицы системы (8) в случае устранения средних параметров дрейфа.

 

Рис. 8. Степени вершин в графе.

 

Рис. 9. Величина $\mathbf{-}\left(k/\omega \right){\mathit{\lambda }}_{\mathbf{m}\mathbf{i}\mathbf{n}}$.

 

Пример разделения группы спутников 

Рассмотренный выше пример является демонстрацией успешного применения предложенных правил для управления в задаче устранения относительных дрейфов МКА и удержания относительных средних сдвигов в заданных интервалах. Однако успешность работы алгоритма управления зависит от ряда параметров. Рассмотрим пример, когда применение управления приводит к разделению группы — образованию нескольких независимых подгрупп спутников. Пусть все параметры моделирования останутся такими же, как и для предыдущего примера, кроме R_comm. Возьмём меньший радиус связи, а именно $R_{comm}={\mu }_D+0.75{\sigma }_D$, что соответствует 200 м для выбранных параметров моделирования. На рис. 10 представлены относительные траектории аппаратов, из которого видно, что спутники удаляются друг от друга. На рис. 11 изображены параметры дрейфов МКА относительно первого запущенного аппарата. В данном примере не все значения относительных параметров дрейфа сходятся к нулевому значению. Это происходит вследствие того, что R_comm был выбран меньше необходимого значения и спутники рассчитывают управление на основе движения лишь некоторых соседних аппаратов. В итоге это привело к тому, что относительные параметры дрейфов сошлись к четырём различным значениям, что соответствует количеству подгрупп, на которые разделился рой.

 

Рис. 10. Относительные траектории движения.

 

Рис. 11. Величины относительных параметров дрейфов.

 

Количество подгрупп, на которые разделяются аппараты, равно количеству компонент связности графа, который образуют спутники с их коммуникационными связями. Изменение числа компонент связности за время моделирования можно увидеть на рис. 12. В течение первых 2 ч 10 мин управления компонента связности одна, и все аппараты образуют связный граф. После этого граф становится несвязным, а к концу моделирования число компонент связности становится равным четырём, что подтверждает разделение МКА на четыре независимые подгруппы. Так как относительные параметры дрейфов сходятся к ненулевым значениям, то из уравнений на параметр D_ij (4) следует, что относительные параметры D_ij для всех спутников линейно растут от времени, что можно увидеть на рис. 13.

 

Рис. 12. Число компонент связности.

 

Рис. 13. Величины относительных сдвигов.

 

Связность группы при различных значениях радиуса связи

Большой интерес представляет исследование эффекта разделения спутников на подгруппы в зависимости от величины R_comm. Так как ошибка в скорости отделения спутников является случайной, то и результат работы алгоритмов по формированию связной группы тоже будет случайной величиной. Исследуем работу алгоритма управления с помощью метода Монте-Карло, варьируя число , которое учитывается при подсчёте R_comm согласно (19). При фиксированном радиусе связи проведём ряд одинаковых численных экспериментов с различающимися скоростями отделения. После каждого эксперимента выполняется проверка, сошлись ли все относительные параметры дрейфов к нулю. В случае разделения группы вычисляется количество МКА в каждой подгруппе с одинаковыми относительными параметрами дрейфов. Обозначим число спутников в самой многочисленной группе как N_group. Будем измерять эффект разделения на подгруппы с использованием отношения числа самой многочисленной группы к общему числу спутников $N_{group}/N$. Если $N_{group}/N=1$, то это означает, что ни один МКА из группы не отделился. Если же отношение $N_{group}/N$ близко к единице, то это соответствует случаю, когда небольшая часть аппаратов покинула группу, но подавляющее большинство осталось в ней.

Для каждого значения параметра $\alpha$, где $R_{comm}={\mu }_D+\alpha \cdot {\sigma }_D$, было проведено сто численных экспериментов c различными скоростями отделения. Значение $N_{group}/N$ и количество подгрупп вычислялись по результатам каждого численного моделирования. Результаты исследования представлены на рис. 14 и 15 в виде диаграммы размаха, где 50% результатов находятся внутри прямоугольника, горизонтальная черта в прямоугольнике соответствует медиане, на отрезки сверху и снизу прямоугольника приходится по 25% результатов, крестики соответствуют выбросам. Из рис. 14 и 15 видно, что начиная с величины $R_{comm}={\mu }_D+2\cdot {\sigma }_D=536\text{\hspace{0.17em}м}$ эффект разделения роя практически не наблюдается, что означает, что предложенный алгоритм успешно справляется с задачей удержания коммуникационной связности роя МКА.

 

Рис. 14. Относительное количество аппаратов в наибольшей подгруппе.

 

Рис. 15. Количество подгрупп.

 

Таким образом, параметрическое численное исследование управляемого движения выявило размеры коммуникационной области, которые зависят от ошибок в скорости отделения космических аппаратов, при которых задача удержания связности роя успешно решается с помощью предложенного алгоритма управления.

ЗАКЛЮЧЕНИЕ

Аналитическое и численное исследование управляемого относительного движения показало, что предложенный децентрализованный закон управления успешно справляется с задачей удержания связности роя малых космических аппаратов с учётом коммуникационных ограничений. Расчётное управление может быть реализовано как с помощью бортовых двигателей малой тяги, так и с использованием разности сил аэродинамического сопротивления на низких околоземных орбитах. Эффект разделения роя на независимые подгруппы проявляется при величинах радиуса коммуникационной сферы меньших, чем удвоенная величина среднеквадратического сдвига относительных траекторий при свободном относительном движении после кластерного запуска с заданными случайными ошибками отделения.

ФИНАНСИРОВАНИЕ РАБОТЫ

Исследование поддержано грантом РНФ № 22-21-00845.

Об авторах

У. В. Монахова

Институт прикладной математики им. М. В. Келдыша РАН

Автор, ответственный за переписку.
Email: danilivanovs@gmail.com
Россия, Москва

С. А. Шестаков

Институт прикладной математики им. М. В. Келдыша РАН

Email: danilivanovs@gmail.com
Россия, Москва

Я. В. Маштаков

Институт прикладной математики им. М. В. Келдыша РАН

Email: danilivanovs@gmail.com
Россия, Москва

Д. С. Иванов

Институт прикладной математики им. М. В. Келдыша РАН

Email: danilivanovs@gmail.com
Россия, Москва

Список литературы

Дополнительные файлы