Известия РАН. Механика твердого тела, 2019, № 4, стр. 109-128

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

Задача оптимальной переориентации плоскости орбиты КА решается в общем случае, когда изменения угловых элементов орбиты в процессе управления могут принимать любые конечные значения. Исследуются задачи быстродействия, минимизации импульса тяги реактивного двигателя, характеристической скорости КА, а также задачи, когда функционал, определяющий качество процесса управления, представляет собой линейную свертку с весовыми множителями двух критериев: 1) времени и суммарного импульса величины тяги, затраченных на процесс управления, 2) времени и характеристической скорости КА.

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

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

2. Дифференциальные уравнения ориентации орбиты космического аппарата и орбитальной системы координат и задачи оптимальной переориентации орбиты и плоскости орбиты КА. Будем считать, что вектор a реактивного ускорения центра масс КА от тяги реактивного двигателя КА, а, следовательно, и вектор тяги u* реактивного двигателя во все время управляемого движения КА направлены ортогонально плоскости оскулирующей орбиты КА, т.е. ортогонально радиусу-вектору r и вектору v скорости центра масс КА (коллинеарно вектору с = r × v момента скорости центра масс КА). Тогда дифференциальные уравнения движения центра масс КА в ньютоновском гравитационном поле, описывающие изменение размеров и формы мгновенной орбиты КА, интегрируются, давая уравнение конического сечения. Поэтому управляемое движение центра масс КА в этом случае описывается дифференциальными уравнениями, описывающими изменение мгновенной ориентации орбиты КА или используемой (например, орбитальной) вращающейся системы координат, в которой записываются исходные уравнения движения центра масс КА, и дифференциальным уравнением для истинной аномалии, характеризующей положение центра масс КА на орбите.

Орбита КА в процессе такого управления движением центра масс КА не меняет своей формы и своих размеров, а поворачивается в пространстве под действием управления a или u* как неизменяемая (недеформируемая) фигура.

Движение центра масс КА будем рассматривать в инерциальной системе координат X – геоцентрической экваториальной системе координат OX₁X₂X₃ (X) c началом в центре O притяжения Земли. Ось OX₃ этой системы координат направлена вдоль оси суточного вращения Земли, оси OX₁ и OX₂ лежат в плоскости экватора Земли, ось OX₁ направлена в точку весеннего равноденствия для Земли, ось OX₂ дополняет систему до правой тройки векторов.

Введем также в рассмотрение систему координат ξ, связанную с плоскостью и перицентром орбиты КА. Начало этой системы координат находится в центре O (или в перицентре орбиты), ось ξ₁ направлена вдоль радиуса-вектора перицентра орбиты, ось ξ₃ перпендикулярна плоскости орбиты и имеет направление постоянного по модулю вектора c момента скорости центра масс КА относительно центра O, а ось ξ₂ образует правую тройку с осями ξ₁ и ξ₃. Ориентация системы координат ξ в инерциальной системе координат X характеризует собой ориентацию орбиты КА в инерциальном пространстве и традиционно задается тремя угловыми оскулирующими элементами орбиты [1, 2]: долготой восходящего узла Ω_u, наклоном орбиты I и угловым расстоянием перицентра от узла ω_π.

Дифференциальные уравнения, описывающие мгновенную ориентацию орбиты КА в инерциальной системе координат в угловых элементах орбиты в рассматриваемом случае ортогональности вектора реактивной тяги плоскости оскулирующей орбиты КА, имеют вид [1, 2]

где φ – истинная аномалия (угловая переменная, отсчитываемая в плоскости орбиты от ее перицентра и характеризующая положение центра масс КА на орбите), r = |r| – модуль радиуса-вектора центра масс КА, p и e – параметр и эксцентриситет орбиты, с = |с| = |r × v| – постоянная площадей (модуль вектора момента скорости v центра масс КА), a – проекция вектора ускорения a = u*/m* центра масс КА на направление его вектора момента скорости (алгебраическая величина реактивного ускорения, перпендикулярного плоскости оскулирующей орбиты КА), m* – масса КА.

Задача переориентации орбиты КА в угловых переменных формулируется следующим образом: требуется построить управление a, переводящее орбиту, изменение ориентации которой описывается уравнениями (2.1) из заданного начального положения

в требуемое конечное положение

При этом должен минимизироваться выбранный функционал качества процесса переориентации орбиты КА.

Частный случай этой задачи рассматривался в работах [3–7]. В работе [4] изучался поворот плоскости околоземной круговой орбиты при помощи тяги, нормальной к мгновенной плоскости орбиты, с использованием усредненных уравнений в угловых элементах орбиты. В работе [6] исследовался поворот плоскости круговой орбиты спутника поперечной тягой (тягой, направленной перпендикулярно к мгновенной плоскости орбиты, называемой также в [6] “бинормальной тягой”). Для описания движения использованы уравнения (2.1), записанные в безразмерных переменных. Задача рассматривалась в предположении, что начальная орбита спутника лежит в плоскости экватора и что требуемый наклон орбиты является малым (поэтому в дальнейшем полагалось sin I ≈ I, cos I ≈ 1).

В работе [7] рассматривалась задача поворота плоскости оскулирующей орбиты КА с помощью “бинормальной силы”, создающей “бинормальное ускорение”, с использованием уравнений для угловых оскулирующих элементов (2.1). Рассмотрение ограничивается случаем круговой орбиты, который, по словам авторов этой работы, был исследован в работах [3, 5]. Анализируется оптимальный в смысле минимизации характеристической скорости поворот плоскости круговой орбиты на угол наклонения орбиты ΔI за неограниченное время. Поворот орбиты на угол ΔI рассматривается в работе [7] в приближенной постановке с использованием лишь второго уравнения системы (2.1) для наклонения орбиты.

В статье [8] рассматривается вековое изменение угловых элементов орбиты Ω_u, I, ω_π под действием реактивного ускорения, ортогонального плоскости оскулирующей орбиты КА. Эта задача называется в статье задачей коррекции элементов орбиты Ω_u, I, ω_π “бинормальным реактивным ускорением”. Используемые в статье [8] исходные уравнения движения КА имеют вид уравнений в угловых элементах (2.1). Для решения задачи авторы статьи переходят в этих уравнениях к новой независимой переменной (истинной аномалии φ) по формуле dt = (r²/c)dφ и дополняют их дифференциальным уравнением для характеристической скорости ${{v}_{{ch}}}$. Граничные условия маневра коррекции записываются в виде

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

Решение задачи оптимальной переориентации орбиты КА посредством реактивного ускорения, ортогонального плоскости оскулирующей орбиты, с помощью уравнений (2.1) в угловых элементах орбиты в строгой нелинейной постановке достаточно сложно в силу нелинейности этих уравнений, наличия в них особых точек I = 0, π, а также в силу громоздкости уравнений для сопряженных множителей. Поэтому для решения этой задачи вместо угловых элементов орбиты целесообразно использовать параметры Эйлера (Родрига–Гамильтона).

Дифференциальные уравнения ориентации орбиты КА в параметрах Эйлера имеют вид [9–13]

где Λ_j (j = 0, 1, 2, 3) – параметры Эйлера, характеризующие ориентацию орбиты КА (системы координат ξ) в инерциальной системе координат X; Ω₁, Ω₂, Ω₃ = 0 – проекции вектора Ω мгновенной абсолютной угловой скорости орбиты на связанные с ней координатные оси Oξ_i.

Параметры Эйлера Λ_j связаны с угловыми элементами орбиты соотношениями

Уравнения (2.2) в кватернионной записи принимают вид [9–13]

где Λ = Λ₀ + Λ₁i + Λ₂j + Λ₃k – кватернион ориентации орбиты КА (кватернионный оскулирующий (медленно изменяющийся) элемент орбиты КА); Ω_ξ – отображение вектора Ω на базис ξ (вектор Ω мгновенной абсолютной угловой скорости орбиты направлен вдоль радиуса-вектора r центра масс КА и определяется формулой: Ω = (a/c) r); i, j, k – векторные мнимые единицы Гамильтона, ○ – символ кватернионного умножения.

Уравнения (2.2), (2.3) или (2.5), (2.3) – система пяти нелинейных стационарных дифференциальных уравнений первого порядка относительно параметров Эйлера Λ_j и истинной аномалии φ. Эти уравнения, в отличие от четырех нелинейных дифференциальных уравнений (2.1) ориентации орбиты в угловых элементах орбиты Ω_u, I, ω_π, не имеют особых точек I = 0, π, к тому же при переходе в них от времени t к новой независимой переменной φ в соответствии с дифференциальным соотношением dφ = = (c/r²)dt получается (при a = a(φ)) система четырех линейных нестационарных дифференциальных уравнений относительно параметров Эйлера Λ_j (в то время как дифференциальные уравнения в угловых элементах орбиты остаются существенно нелинейными).

Указанные обстоятельства делают использование уравнений (2.2), (2.3) или (2.5), (2.3) для решения задач переориентации орбиты, плоскости орбиты и коррекции угловых элементов орбиты более удобным и эффективным в сравнении с использованием уравнений в угловых оскулирующих элементах (2.1). Такое решение задачи переориентации орбиты КА в непрерывной постановке (с использованием ограниченной (малой) тяги) рассмотрено в работах [14–17]. В них в качестве функционала качества процесса переориентации орбиты КА рассмотрены комбинированный функционал, равный взвешенной сумме времени переориентации и интеграла от квадрата модуля управления, а также комбинированный функционал, равный взвешенной сумме времени переориентации и импульса управления (характеристической скорости) за время переориентации орбиты КА. В работах [18, 19] с использованием уравнений (2.5), (2.3) построена теория решения задачи оптимальной переориентации орбиты КА в импульсной постановке (с использованием импульсной (большой) реактивной тяги). Приведены алгоритмы решения краевых задач оптимальной двухимпульсной и многоимпульсной переориентации орбиты КА (для нефиксированного числа импульсов реактивной тяги) и примеры численного решения краевых задач оптимальной переориентации орбиты КА с использованием ограниченной (малой) или импульсной (большой) тяги, в которых для описания ориентации орбиты КА используется кватернионный оскулирующий элемент ориентации орбиты.

Отметим, что в работах [14–19] для решения задач переориентации орбиты КА использован принцип максимума, а в качестве управления использовано реактивное ускорение a центра масс КА.

Наряду с дифференциальными уравнениями (2.2), (2.3) или (2.5), (2.3) ориентации орбиты КА в параметрах Эйлера Λ_j для решения задач переориентации орбиты, плоскости орбиты КА и коррекции угловых элементов КА могут быть использованы дифференциальные уравнения ориентации орбитальной системы координат η в параметрах Эйлера λ_j, имеющие вид [11, 12, 20, 21]

где λ_j (j = 0, 1, 2, 3) – параметры Эйлера, характеризующие ориентацию орбитальной системы координат η в инерциальной системе координат X (ось η₁ этой системы координат направлена вдоль радиуса-вектора r центра масс КА, а ось η₃ перпендикулярна плоскости орбиты (параллельна оси ξ₃)); ω₁, ω₂ = 0, ω₃ – проекции вектора ω мгновенной абсолютной угловой скорости орбитальной системы координат на ее же координатные оси.

Отметим, что управление a входит в уравнения (2.6) мультипликативно в одну из двух компонент вектора абсолютной угловой скорости вращения орбитальной системы координат, направленную вдоль радиуса-вектора центра масс КА, а в уравнениях (2.2) это управление входит мультипликативно в две компоненты вектора абсолютной угловой скорости вращения орбиты КА.

Параметры λ_j связаны с угловыми переменными Ω_u, I, ω_π, φ соотношениями, аналогичными соотношениям (2.4):

Уравнения (2.6) в кватернионной записи принимают вид [11, 12, 20, 21]

где λ = λ₀ + λ₁i + λ₂j + λ₃k – кватернион инерциальной ориентации орбитальной системы координат, связанный с кватернионом Λ ориентации орбиты КА соотношением ω_η – отображение вектора ω на базис η.

Отметим также, что кватернионное дифференциальное уравнение ориентации орбитальной системы координат, аналогичное (2.9), также использовалось для описания орбитального движения в работах [22, 23].

Уравнения (2.6), (2.7) или (2.7), (2.9) использовались в работах [12, 24–26] для решения c использованием принципа максимума Понтрягина задачи оптимальной переориентации орбиты КА посредством реактивной тяги, ортогональной плоскости орбиты, в непрерывной постановке (с использованием ограниченной (малой) тяги и с использованием в качестве управления реактивного ускорения a). В качестве критерия оптимальности был использован комбинированный функционал, равный взвешенной сумме времени переориентации и интегрального квадратичного (в отношении управления) функционала качества, или функционал, равный взвешенной сумме времени переориентации и импульса управления (характеристической скорости) за время переориентации орбиты КА. Управление полагалось ограниченным по модулю.

В работах [27, 28] уравнения (2.7), (2.9) использованы для решения задачи оптимальной переориентации орбиты КА посредством импульсной (большой) реактивной тяги, ортогональной плоскости оскулирующей орбиты. В качестве управления использовано ускорение a центра масс КА от тяги реактивного двигателя. Минимизируется комбинированный функционал, равный взвешенной сумме времени переориентации и импульса реактивного ускорения центра масс КА (характеристической скорости) за время переориентации орбиты КА. Решение задачи построено с помощью предельных переходов в уравнениях и соотношениях, полученных в результате решения задачи оптимальной переориентации орбиты КА в непрерывной постановке (с использованием ограниченной (малой) реактивной тяги и принципа максимума Понтрягина). Алгоритмы импульсного решения задачи, построенные в [28], позволяют определять оптимальные моменты включения реактивного двигателя, оптимальные величины импульсов реактивного ускорения КА и их оптимальное число. Приведены примеры численного решения задачи оптимальной импульсной переориентации орбиты КА, демонстрирующие возможности предлагаемого метода. В [28] также показано, что задача оптимальной импульсной переориентации орбиты КА в случае, когда оптимальное управление состоит из двух импульсов реактивного ускорения, сообщаемых КА в начальный и конечный моменты времени движения, решается аналитически.

Использование как уравнений (2.2), (2.3) (или (2.5), (2.3)), описывающих собой ориентацию орбиты КА и положение КА на орбите, так и уравнений (2.6), (2.7) (или (2.9), (2.7)), описывающих собой ориентацию орбитальной системы координат и положение КА на орбите, для решения задач оптимальной переориентации орбиты КА имеет свои достоинства. Так, уравнения (2.6), (2.7) являются при r = const (в случае круговой орбиты) и a = const линейными дифференциальными уравнениями с постоянными коэффициентами, в то время как уравнения (2.2), (2.3) в этом случае являются линейными дифференциальными уравнениями с переменными коэффициентами. Поэтому уравнения (2.6), (2.7) (или (2.9), (2.7)) более удобны и эффективны в сравнении с уравнениями (2.2), (2.3) (или (2.5), (2.3)) с аналитической точки зрения. Однако в состав уравнений (2.2), (2.3) входят уравнения (2.2) в скалярных оскулирующих элементах Λ_j, а в состав уравнений (2.5), (2.3) – уравнение (2.5) для кватернионного оскулирующего элемента орбиты Λ. Переменные Λ_j, а также составленный из них кватернион Λ, как уже отмечалось, непосредственно характеризуют собой ориентацию орбиты КА в инерциальном пространстве. При равенстве нулю управления a (реактивного ускорения) эти переменные становятся постоянными величинами, поэтому переменные Λ_j являются скалярными оскулирующими (медленно изменяющимися) элементами орбиты КА, а переменная Λ – кватернионным оскулирующим (медленно изменяющимся) элементом орбиты КА. Такими свойствами скалярные переменные λ_j и эквивалентная им кватернионная переменная λ, фигурирующие в уравнениях (2.6), (2.7) и (2.9), (2.7), не обладают (они являются быстро меняющимися функциями времени). Поэтому изучение задач оптимальной переориентации орбиты, плоскости орбиты и коррекции угловых элементов орбиты КА представляется актуальным как с использованием кватернионного дифференциального уравнения ориентации орбитальной системы координат (2.9), так и с использованием кватернионного дифференциального уравнения ориентации орбиты КА (2.5).

В настоящей статье задача оптимальной переориентации плоскости орбиты КА изучается с использованием уравнений (2.6), (2.7) (или (2.9), (2.7)), более удобных с аналитической точки зрения. При этом в качестве управления принимается не алгебраическая величина a реактивного ускорения центра масс КА, как это делалось во всех цитированных выше работах, а алгебраическая величина u* вектора реактивной тяги u* = m*(t)a двигателя КА, ортогонального плоскости оскулирующей орбиты КА. Масса аппарата m* полагается кусочно дифференцируемой функцией времени t: m* = m*(t). В состав математической модели, описывающей процесс переориентации плоскости орбиты КА, дополнительно включается дифференциальное уравнение, описывающее изменение массы КА в процессе управляемого движения.

3. Постановка задачи оптимального поворота плоскости орбиты КА переменной массы посредством ограниченной (малой) реактивной тяги. Дифференциальные уравнения задачи переориентации орбиты КА переменной массы (и задачи переориентации плоскости орбиты КА) в инерциальной системе координат с помощью реактивной тяги, ортогональной плоскости оскулирующей орбиты КА, в случае использования для решения задачи кватернионного дифференциального уравнения ориентации орбитальной системы координат имеют следующий вид:

где β* – постоянный коэффициент пропорциональности, равный обратной величине скорости истечения рабочего тело из сопла реактивного двигателя [29].

Ориентация плоскости орбиты КА в инерциальной системе координат характеризуется направляющими косинусами n_k единичного вектора n нормали к плоскости орбиты в этой системе координат, которые запишем через параметры Эйлера λ_j и угловые элементы орбиты КА в следующем виде:

Перейдем к безразмерным переменным ρ, τ, m, u, u_m, β по формулам

где $m_{0}^{*}$ – начальная масса аппарата, $u_{m}^{*}$ = $u_{{\max }}^{*}$ – максимальное значение тяги.

Уравнения движения КА (3.1), (3.2) в безразмерных переменных примут вид

Постановка задачи: требуется определить ограниченное по модулю управление u:

ортогональное плоскости оскулирующей орбиты космического аппарата, переводящее КА, движение центра масс которого описывается уравнениями (3.5), (3.6), из заданного начального состояния в конечное состояние, принадлежащее многообразию и минимизирующее комбинированный функционал качества (3.10) или (3.11):

Функционалы (3.10) и (3.11) представляют собой линейные свертки с постоянными весовыми множителями α₁ и α₂ двух критериев: (3.10) – времени и импульса величины тяги реактивного двигателя, затраченных на процесс управления, (3.11) – времени и характеристической скорости КА.

Безразмерное время τ₁ процесса управления заранее не задается. В частном случае, когда ${{\alpha }_{1}} \ne 0,\,\,{{\alpha }_{2}} = 0$, критерием оптимальности является длительность процесса управления и задача оптимального управления становится задачей быстродействия. В других частных случаях, когда ${{\alpha }_{1}} = 0,\,{\kern 1pt} {\kern 1pt} {{\alpha }_{2}} = 1$, критерием оптимальности является импульс величины тяги реактивного двигателя (функционал (3.10)) или характеристическая скорость КА (функционал (3.11)).

Фигурирующие в начальных условиях кватернион Λ(0) и скалярная величина φ₀ характеризуют начальную ориентацию орбиты КА и начальное положение КА на орбите. Величины c, p, e, φ₀, Λ(0) заданы (начальное Λ(0) значение кватерниона Λ ориентации орбиты КА может быть найдено через заданные начальные значения угловых элементов орбиты Ω_u, I, ω_π по формулам (2.4)). Также заданы конечные значения Ω_u(τ₁), I(τ₁) угловых элементов орбиты Ω_u, I, характеризующие требуемую (конечную) ориентацию плоскости орбиты КА и входящие в уравнения (3.9) многообразия, которому принадлежит конечное состояние управляемой системы. Подлежат определению оптимальный закон управления u = u(τ), под действием которого плоскость орбиты КА займет требуемое положение, и величины τ₁, φ₁, характеризующие время процесса управления и конечное положение КА на орбите.

Поставленная задача – задача оптимального управления с подвижным правым концом траектории. Условия (3.9) на правом конце траектории записаны в соответствии с соотношениями (3.3) для направляющих косинусов n_k единичного вектора n нормали к плоскости орбиты КА. Отметим, что направляющие косинусы n_k связаны условием $n_{1}^{2}$ + $n_{2}^{2}$ + $n_{3}^{2}$ = 1. Поэтому для решения задачи взяты только первые два из трех условий (3.3) на правом конце траектории. Вместо них могут быть взяты первое и третье из условий (3.3) или второе и третье из этих условий.

4. Законы (структура) оптимального управления и условия трансверсальности. Поставленную задачу будем решать с помощью принципа максимума. Для этого введем дополнительные переменные μ, χ и l, сопряженные по отношению к фазовым переменным λ, φ и m. Функция Гамильтона–Понтрягина в случае минимизации функционала (3.10) имеет вид

а в случае минимизации функционала (3.11) – вид

Здесь μ_j – компоненты сопряженной кватернионной переменной μ,

В соответствии с принципом максимума постоянная ψ₀ ≤ 0. В дальнейшем рассматривается случай, когда постоянная ψ₀ < 0. В этом случае в силу однородности функции Гамильтона–Понтрягина относительно сопряженных переменных константу ψ₀ можно выбрать произвольным образом, поэтому в дальнейшем нами полагается ψ₀ = – 1.

Система уравнений для сопряженных переменных в случае минимизации функционала (3.10) имеет вид уравнений (4.4)–(4.6):

В случае минимизации функционала (3.11) система уравнений для сопряженных переменных имеет вид уравнений (4.4), (4.5) и уравнения (4.7):

Введем новую кватернионную переменную

которая является мультипликативной композицией фазовой λ и сопряженной μ кватернионных переменных. В формуле (4.8) и далее верхняя черта – символ кватернионного сопряжения.

Компоненты ${{{\nu }}_{j}}$ кватерниона μ в соответствии с (4.8) определяются формулами

С использованием новых переменных ν_j выражения (4.1), (4.2), (4.3) для функции Гамильтона–Понтрягина и сопряженные уравнения (4.5), (4.6) и (4.7) принимают более простой вид:

Отметим, что кватернионное сопряженное уравнение (4.4) совпадает по своей форме с кватернионным фазовым уравнением (3.5), поскольку кватернионное уравнение (3.5) обладает свойством самосопряженности. Отметим также, что свойство самосопряженности кватернионного кинематического уравнения было впервые установлено В.Н. Бранцем и И.П. Шмыглевским при решении задач оптимальных пространственных разворотов твердого тела [30].

Из (4.5) или (4.12) видно, что в случае круговой орбиты, когда эксцентриситет орбиты e = 0, сопряженная переменная

Закон оптимального управления (т.е. закон управления, удовлетворяющий необходимым условиям оптимальности) находится из условия максимума функции H по переменной u с учетом наложенного ограничения (3.7) и в случае минимизации функционала (3.10) имеет вид (4.16):

а в случае минимизации функционала (3.11) – вид (4.17):

Отметим, что особый режим управления в статье не рассматривается.

Условия трансверсальности, соответствующие многообразию конечного состояния (3.9), после исключения из них двух неопределенных множителей Лагранжа принимают следующий вид:

при τ = τ₁

5. Анализ задачи оптимальной переориентации плоскости орбиты КА. В случае минимизации функционала (3.10) задача сведена к краевой задаче с подвижным правым концом траектории, описываемой системой нелинейных дифференциальных уравнений (3.5), (3.6), (4.4), (4.12), (4.13), (4.16) двенадцатого порядка и восемью краевыми условиями (3.8), (3.9), которые необходимо дополнить четырьмя условиями трансверсальности (4.18), (4.19) и равенством

имеющим место для оптимального управления u и оптимальной траектории на ее правом конце.

В случае минимизации функционала (3.11) задача сведена к краевой задаче с подвижным правым концом траектории, описываемой системой нелинейных дифференциальных уравнений (3.5), (3.6), (4.4), (4.12), (4.14), (4.17) двенадцатого порядка и восемью краевыми условиями (3.8), (3.9), которые необходимо дополнить четырьмя условиями трансверсальности (4.18), (4.19) и равенством

имеющим место для оптимального управления u и оптимальной траектории в конечный момент τ₁ безразмерного времени.

Дифференциальные уравнения задачи имеют первые интегралы

Интегралы (5.3)–(5.6) имеют место для любого управления u. Интегралы (5.3), (5.4) отражают неизменность норм кватернионных переменных λ и μ. Интеграл (5.5) имеет место в силу самосопряженности кватернионного фазового уравнения (3.5).

В силу автономности системы уравнений для фазовых и сопряженных переменных равенства (5.1) и (5.2) будут выполняться в любой момент времени оптимального управляемого движения, т.е. будут являться первыми интегралами этих уравнений для оптимального управления u. В случае круговой орбиты уравнения задачи имеют также первый интеграл (4.15): χ = χ₀ = const.

На активном участке траектории оптимальное управление u является постоянной величиной, поэтому дифференциальное уравнение для массы m (второе уравнение подсистемы (3.6)) интегрируется и дает соотношение

описывающее изменение массы КА на этом участке и показывающее, что на активном участке безразмерная масса m аппарата является линейной функцией безразмерного времени τ (уменьшается линейно с течением времени τ).

Размерная q* и безразмерная q характеристические скорости КА связаны соотношением q* = (c/R)q. Изменение безразмерной характеристической скорости описывается дифференциальным уравнением

Из сопоставления этого уравнения с комбинированным функционалом (3.11) видно, что характеристическая скорость q является частным случаем функционала (3.11), когда коэффициенты этого функционала ${{\alpha }_{1}} = 0,\,{\kern 1pt} {\kern 1pt} {{\alpha }_{2}} = 1$.

Интегрирование уравнения (5.8) с учетом (5.7) на активном участке траектории дает следующий закон изменения характеристической скорости в процессе управляемого движения КА на этом участке траектории:

где m₀ и q₀ – значения m и q в момент времени, соответствующего началу рассматриваемого активного участка траектории.

Оба функционала (3.10) и (3.11) можно выразить через конечное безразмерное время движения τ₁ и конечное значение безразмерной массы аппарата m₁. Имеем:

Таким образом, поставленная выше задача Лагранжа оптимального управления сводится к задаче Майера. Из последних приведенных формул следует, что при α₁ = 0 задачи минимизации импульса тяги и характеристической скорости сводятся к задаче минимизации топлива. Это также видно из ниже приводимых таблиц 4, 5 и 7 численных расчетов.

Из функций Гамильтона–Понтрягина (4.10), (4.11) и законов оптимального управления (4.16), (4.17) видно, что они фактически являются функциями не сопряженных переменных μ_j (j = 0, 1, 2, 3) (компонент кватернионной сопряженной переменной μ), а функциями скалярных переменных ${{{\nu }}_{1}}$ и ${{{\nu }}_{{3~}}}$ (компонент кватерниона ν), которые в свою очередь являются мультипликативными композициями фазовых λ_j и сопряженных μ_j переменных, имеющими вид (4.9). При этом величина ν₀ = const (она равна компоненте $\nu _{0}^{*}$ кватернионной константы ν* управляемого движения, определяемой соотношением (5.5)) и входит она лишь в условие трансверсальности (4.19). Поэтому вместо выше перечисленных двенадцати дифференциальных уравнений краевой задачи оптимального управления относительно фазовых и сопряженных переменных λ_j (j = 0, 1, 2, 3), φ, m, μ_j, χ, l можно использовать для решения задачи оптимальной переориентации плоскости орбиты КА систему одиннадцати дифференциальных уравнений относительно переменных λ_j, φ, m, ν_k (k = 1, 2, 3), χ, l, опуская дифференциальное уравнение для переменной ν₀, переходящее в соотношение ν₀ = const, и условие трансверсальности (4.19), содержащее величину ν₀.

Дифференцируя кватернионное соотношение (4.8) по переменной τ и учитывая уравнения (3.5), (4.4), получаем следующие дифференциальные уравнения для скалярной ${{{\nu }}_{0}}$ и векторной ${{{\mathbf{\nu }}}_{{v}}}$ = ν₁i + ν₂j + ν₃k частей кватернионной переменной ν:

Из уравнения (5.9) имеем: ${{{\nu }}_{{0~}}} = ~{{{\nu }}_{0}}(0) = {\text{const}}.$ Из векторного уравнения (5.10) получаем следующую систему трех скалярных дифференциальных уравнений для компонент ν_k векторной переменной ${{{\mathbf{\nu }}}_{{v}}}$:

Уравнения (5.11) ((5.10)) имеют (для любого управления) первый интеграл

Уравнения (5.11) являются дифференциальными уравнениями функции переключения ${{{\nu }}_{1}}$ = ${{{\nu }}_{1}}({\tau })$, поскольку управление, как это видно из соотношений (4.16) и (4.17), переключается с одной ветви на другую в зависимости от значения координаты ν₁. Переменные ν_k, входящие в уравнения (5.11), связаны с кватернионным первым интегралом (5.5) уравнений краевой задачи преобразованием вращения: ${{{\mathbf{\nu }}}_{{v}}} = ~{\mathbf{\bar {\lambda }}}~ \circ ~{\mathbf{\nu }}_{{v}}^{*}~ \circ ~{\mathbf{\lambda }}.$

Таким образом, задача оптимальной переориентации плоскости орбиты КА в непрерывной постановке (с использованием ограниченной (малой) реактивной тяги) сводится к решению краевой задачи с подвижным правым концом траектории, описываемой в случае минимизации функционала (3.10) системой одиннадцати нелинейных стационарных дифференциальных уравнений (3.5), (3.6), (5.11), (4.12), (4.13), (4.16), а в случае минимизации функционала (3.11) – системой одиннадцати нелинейных стационарных дифференциальных уравнений (3.5), (3.6), (5.11), (4.12), (4.14), (4.17) относительно переменных λ_j, φ, m, ν_k, χ, l и восемью краевыми условиями (3.8), (3.9), которые необходимо дополнить тремя условиями трансверсальности (4.18) и равенством (5.1) в первом случае или равенством (5.2) во втором случае.

Оптимальное управление, как это видно из (4.16) и (4.17), состоит из двух возможных режимов: активного участка движения с постоянной максимально возможной по модулю реактивной тягой u = ±u_m, направление (знак) которой определяется знаком переменной ν₁, и пассивного участка движения с нулевой тягой (напомним, что особый режим управления в статье не рассматривается).

Как видно из условий трансверсальности (4.18), переменные ν₃ и χ на правом конце траектории равны нулю, поэтому условие (5.1) (в случае минимизации функционала (3.10)) и условие (5.2) (в случае минимизации функционала (3.11)), имеющие место для оптимального управления u и оптимальной траектории на этом конце траектории принимают соответственно вид (5.12) и (5.13):

Сопряженная переменная χ и второе условие трансверсальности (4.18) исключаются из рассмотрения, а краевая задача будет описываться системой нелинейных стационарных дифференциальных уравнений (3.5), (3.6), (5.11), (4.13), (4.16) или (3.5), (3.6), (5.11), (4.14), (4.17) десятого порядка относительно переменных λ_j, φ, m, ν_k, l и восемью краевыми условиями (3.8), (3.9), которые необходимо дополнить двумя условиями трансверсальности ν₃ = 0, l = 0 и равенством (5.12) или (5.13) на правом конце траектории.

Отметим, что возможно понижение размерности исходной краевой задачи, описываемой дифференциальными уравнениями для переменных λ_j, φ, m, ν_k, χ, l одиннадцатого порядка, на две единицы за счет перехода к новой независимой переменной – истинной аномалии φ (при таком переходе исключаются из рассмотрения дифференциальные уравнения для истинной аномалии φ и сопряженной к ней переменной χ, а дифференциальные уравнения краевой задачи становятся нестационарными).

Отметим также, что дифференциальные уравнения (5.11), описывающие функцию ν₁ = ν₁(φ) переключения управления, и второе дифференциальное уравнение подсистемы (3.6) для массы m при переходе в них к новой независимой переменной φ (истинной аномалии) выделяются для оптимального управления u = ±u_m в замкнутую подсистему четырех нелинейных нестационарных дифференциальных уравнений первого порядка относительно переменных ν_k, m, имеющую вид

6. Результаты численного решения задачи оптимальной переориентации плоскости орбиты КА и их анализ. Приводимые ниже решения задачи для одних и тех же исходных данных, но для разных функционалов качества процесса управления, представлены в безразмерных переменных, угловые величины выражены в градусах. Эксцентриситет орбиты е = 0.1, коэффициент, характеризующий расход топлива на создание единицы тяги β = 4.0, начальная масса аппарата m₀ = 1.0. Начальная ориентация орбиты КА определяется угловыми элементами орбиты I_n = 65.8°, Ω_un = 217.0°, ${{\omega }_{{\pi n}}} = 0.0^\circ ;$ положение КА на орбите в начальный момент τ = 0 определяется истинной аномалией ${{\varphi }_{n}}$ = = 30.0° . Ориентация плоскости, на которую необходимо перевести орбиту КА, определяется угловыми элементами I_k = 64.8°, Ω_uk = 215.25°.

6.1.1. Решение задачи быстродействия для максимальной тяги ${{u}_{m}}$ = 0.004. Оптимальное управление состоит из 4-х этапов с постоянными управлениями ${{u}_{m}}$ или $ - {{u}_{m}}$. За время процесса управления КА совершает один полный виток. В табл. 1 описаны этапы управления: приводятся моменты ${{\tau }_{i}}$ окончания этапов, истинная аномалия φ_i, угловые элементы орбиты, масса m_i КА в конце i-го этапа, значение u_i управления на этом этапе. Из табл. 1 видно, что длительность оптимального процесса составляет 10.998431 единиц безразмерного времени и за это время израсходовано 0.175975 топлива (от начальной массы аппарата).

6.1.2. Решение задачи быстродействия для максимальной тяги u_m = 0.002. Оптимальное управление состоит из 8 этапов с постоянными управлениями u_m или –u_m. За время процесса управления КА совершает 3 полных витка. В табл. 2 описаны этапы управления. Видно, что длительность оптимального процесса составляет 22.948044 единиц безразмерного времени и за это время израсходовано 0.183584 топлива (от начальной массы аппарата). На каждом внутреннем этапе (т.е. за исключением первого и последнего этапов) с управлением u = 0.002 истинная аномалия изменяется на 180.0004°, а на этапе с управлением u = –0.002 на 180.0002°. Длительности внутренних этапов с управлением u = 0.002 постепенно уменьшались на 0.000004, а длительности этапов управлением u = – 0.002 увеличивалась на 0.000002 единиц безразмерного времени.

В табл. 3 приводятся результаты решения задачи быстродействия для различных значений максимальной тяги. Нижний индекс k у величин τ, φ, ${{{\omega }}_{{\pi }}}$ означает, что приводятся их значения для конечного момента времени τ_k. Видно, что с уменьшением максимальной величины управления u_m увеличивается время τ_k процесса управления, число этапов управления и полных витков траектории КА. Так, для ${{u}_{m}} = ~$ 0.0005 число этапов управления и полных витков траектории КА становится равным 30 и 14 соответственно.

6.2. Решение задачи в случае минимизации импульса величины тяги реактивного двигателя КА. В этом случае коэффициенты функционалов оптимизации ${{\alpha }_{1}} = 0,\,{\kern 1pt} {\kern 1pt} {{\alpha }_{2}} = 1$. Установлено, что краевая задача оптимизации допускает в этом случае при одних и тех же граничных условиях бесчисленное множество решений. В полученных решениях для максимальной тяги u_m = 0.001 с увеличением длительности процесса управления уменьшается значение импульса величины тяги и по мере увеличения длительности оно стремится к некоторому пределу. В табл. 4 приведены решения задачи оптимального управления, в которых в одном полном витке траектории укладываются четыре этапа с управлениями u = 0.001, u = 0.0, u = –0.001, u = 0.0.

С увеличением длительности процесса управления изменяется структура оптимального управления. В один полный виток укладываются только два этапа с управлениями u = 0.0, u = 0.001. При длительности τ_k = 214.875157 значение функционала J = = 0.029935, а при τ_k = 297.674094,  J = 0.029557. Результаты для этих длительностей процесса управления представлены в табл. 5.

6.3. Решение задачи в случае минимизации комбинированного функционала, включающего в себя время и импульс величины тяги реактивного двигателя КА. Установлено, что краевая задача оптимизации допускает в этом случае при одних и тех же исходных данных задачи несколько решений, из которых необходимо выбрать оптимальное (с меньшим значением функционала). В табл. 6 (в ней $\Delta {{m}_{k}}$ – безразмерный расход топлива, ${{\varphi }_{k}}$ – конечная истинная аномалия) приводятся два решения такой задачи для максимальной тяги u_m = 0.001 в случае минимизации функционала

Из таблицы видно, что оптимальным является решение, представленное в первой строке. Решения для оптимальных управлений на каждом полном витке траектории содержат четыре этапа: u = 0.001, u = 0.0, u = –0.001, u = 0.0.

6.4. Решение задачи в случае минимизации характеристической скорости КА. Установлено, что краевая задача оптимизации допускает в этом случае много решений при одних и тех же исходных данных задачи. В табл. 7 приводится несколько таких решений для случая ${{u}_{m}} = 0.001$. Решения для оптимальных управлений, представленные в двух верхних строках, в каждом полном витке содержат четыре этапа: u = 0.001, u = 0.0, u = –0.001, u = 0.0. В решениях, представленных в третьей, четвертой и пятой строках, при указанном в таблице времени τ* происходит переход от режима с четырьмя этапами в одном полном витке траектории к режиму с двумя этапами в одном полном витке: u = 0.0, u = –0.001.

6.5. Решение задачи в случае минимизации комбинированного функционала

Решение для оптимального управления, приведенное в табл. 8, в каждом полном витке содержит четыре этапа: u = 0.001, u = 0.0, u = –0.001, u = 0.0.

7. Заключение. Решена с использованием кватернионного дифференциального уравнения ориентации орбитальной системы координат и принципа максимума Понтрягина задача оптимального поворота плоскости орбиты КА переменной массы в инерциальной системе координат. Управление поворотом плоскости орбиты КА на любые по величине углы производится с помощью ограниченной по модулю реактивной тяги, ортогональной плоскости оскулирующей орбиты КА. Под действием такой тяги форма и размеры орбиты КА остаются неизменными во все время процесса управления, а сама орбита поворачивается в инерциальном пространстве как неизменяемая (недеформируемая) фигура (окружность, эллипс). Учитывается изменение массы аппарата за счет расхода рабочего тела на процесс управления. Частным случаем изучаемой задачи является задача оптимальной коррекции угловых элементов орбиты КА, имеющая важное значение в механике космического полета.

Рассмотрены задачи быстродействия, минимизации импульса тяги, характеристической скорости КА, а также задачи, когда функционал, определяющий качество процесса управления, представляет собой линейную свертку с весовыми множителями двух критериев: 1) времени и суммарного импульса величины тяги, затраченных на процесс управления, 2) времени и характеристической скорости КА. Структура оптимального управления имеет в зависимости от выбранного функционала качества процесса управления различный характер и включает в себя пассивные участки траектории, на которых тяга отсутствует и плоскость орбиты остается неизменной, и активные участки траектории, на которых реактивная тяга принимает максимальное значение и происходит поворот плоскости орбиты КА.

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