КОСМИЧЕСКИЕ ИССЛЕДОВАНИЯ, 2012, том 50, № 3, с. 258-270
УДК 629.78:517.977
МЕТОД ПРОДОЛЖЕНИЯ ДЛЯ ОПТИМИЗАЦИИ МЕЖПЛАНЕТНЫХ
ТРАЕКТОРИЙ С МАЛОЙ ТЯГОЙ © 2012 г. В. Г. Петухов
Научно-исследовательский институт прикладной механики и электродинамики, г. Москва
vgpetukhov@mail.ru Поступила в редакцию 06.07.2011 г.
Рассматривается задача локальной оптимизации межпланетных траекторий с малой тягой с использованием принципа максимума и численных методов продолжения. Анализируются два типа задач: задачи с ограниченной мощностью и задачи с ограниченной тягой. Последняя задача обобщается введением зависимости тяги и удельного импульса от располагаемой электрической мощности. Для сведения задачи оптимального управления к краевой задаче используется принцип максимума Понтрягина, а затем методом продолжения эта краевая задача сводится к задаче Коши. В статье приводятся варианты метода продолжения для оптимизации траекторий с малой тягой, включая новый метод продолжения для задачи с ограниченной тягой, не требующий выбора начального приближения для граничных значений сопряженных переменных.
1. ВВЕДЕНИЕ
Рассматривается оптимизация траектории с малой тягой между двумя планетами в ньютоновском гравитационном поле. Точнее, анализируется задача локальной оптимизации гелиоцентрической траектории между двумя планетами. Эта задача все еще остается достаточно трудной: современные численные методы не гарантируют получения решения типичных задач оптимизации траекторий, особенно задач с переключением тяги.
Существует два основных подхода к оптимизации траекторий с малой тягой: прямой и непрямой. Прямой подход предполагает использование некоторой схемы дискретизации или параметризации управления. В результате, задача оптимального управления (ЗОУ) сводится к задаче нелинейного программирования с большим числом ограничений и неизвестных переменных. Следовательно, задача нелинейного программирования наследует сложность ЗОУ. Основные трудности прямого подхода связаны с выбором схемы дискретизации и эффективностью метода решения задачи нелинейного программирования.
Непрямой подход предполагает использование необходимых или достаточных условий оптимальности. Обычно ЗОУ сводится к краевой задаче (КЗ). Таким образом, в этом случае требуется эффективное средство для решения таких КЗ.
В этой работе рассматривается непрямой подход, основанный на принципе максимума. Соответствующая КЗ решается с использованием ньютоновской гомотопии, реализующей формальную редукцию КЗ к задаче Коши.
Рассматриваются два типа математических моделей малой тяги: задачи с ограниченной мощностью (ОМ) и с ограниченной тягой (ОТ). В рамках допущений ОМ-задачи тяга и удельный импульс ограничены только реактивной мощностью двигателя, тогда как в рамках ОТ-задачи ограничены и тяга, и удельный импульс.
Методы решения ОМ-задачи рассматривались во многих работах, например в [1—5, 9—11]. ОТ-за-дача более реалистична, но более трудна для численного анализа. Целью этой работы является разработка устойчивого и эффективного численного метода для продолжения ОМ-траектории в ОТ-траекторию.
Более того, целью является автоматизация решения ОТ-задачи. Для достижения этой цели необходима либо автоматизация выбора начального приближения для начальных значений сопряженных переменных, либо вычислительная схема решения КЗ без выбора этого начального приближения. Именно последний вариант был реализован в методе продолжения по гравитационному параметру, в котором нулевое начальное приближение для начальных значений сопряженных переменных (соответствующее пассивному движению вдоль начальной орбиты) использовалось для вычисления оптимальной ОМ-траектории с заданной угловой дальностью. Семейство методов продолжения, в частности метод продолжения по гравитационному параметру для решения ОМ-задачи был разработан автором этой статьи в 1994—1995 гг. и с тех пор активно эксплуатировался и неоднократно публиковался (см., например [2—5, 9—11]). Для автоматизации процесса вычисления оптимальной ОТ-траектории требуется эффектив-
ный метод продолжения из ОМ- в ОТ-решение. Эта задача оказалась значительно более сложной: такой метод удалось построить относительно недавно [9].
2. ЗАДАЧА С ОГРАНИЧЕННОЙ МОЩНОСТЬЮ
Рассмотрим ОМ-задачу, т.е. примем, что тяга и удельный импульс ограничены только реактивной мощностью. Реактивная мощность Щ равна половине произведения тяги Р на скорость истечения м>:
N1 = Ръ/ 2,
Р2
2 2 ш а
2Nj 2Nj 0п (х, г)'
(2)
' 1 10
где а — величина реактивного ускорения. Уравнение (2) имеет следующее решение:
ш (г) =
ш.
1 + ш.
N10
I (г)
(3)
где
1 (') = 21:
2* п (х, г)
(4)
— функционал (показатель качества) ОМ-задачи.
Векторы положения х и скорости V КА удовлетворяют дифференциальным уравнениям:
йх = V, ^ = о + а,
. х . (5)
йг йг
где ^ — силовая функция (^ = 1/|х| для центрального ньютоновского поля в безразмерных координатах) и а — реактивное ускорение.
Применяя принцип максимума к ЗОУ (2, 4, 5), имеем следующий гамильтониан:
2 (
н = - а
2П
1 +
2
ш р
N
+ р^О х + Р Уа + Р я
(6)
у0
гдерт, рх, pv — сопряженные переменные к т, х и V соответственно, а = |а|. Следовательно, оптимальное управление имеет вид:
а = ■-—
N,-0 + ш р„
(7)
(1) а уравнения оптимального движения —
С другой стороны, реактивная мощность может быть представлена как произведение к.п.д. электроракетной двигательной установки (ЭРДУ) Цеж на потребляемую ЭРДУ электрическую мощность Ще. И пЕР5-, и Ще могут зависеть от вектора положения КА х и времени I из-за изменения условий освещенности Солнцем, деградации солнечных батарей и двигателей и т.д. Поэтому, в общем случае, Щ также зависит от х и I. Следовательно,
Щ = П(х, 1Що,
где п(х, о = Перу(х, *)Ще(х, 1)/Що, Що = Перу(хо, 10)Ще(х0, 10) — начальная реактивная мощность, х0 — начальный вектор положения КА и 10 — начальное время. Учитывая зависимость Р = шм>, где т — масса КА, дифференциальное уравнение для массы КА примет вид:
л2
а х
йг2
= П х
йш
йг
N-0 + ш р„,
цNjom2p2v 2(( + ш2рш)
Рv,
л2
= П р +
йг
N-0 р.. дп
2 (( 0 + ш2 рш ))
(8)
йрш _ ЦNj0шРшРу
йг
2
((0 + ш2 рш )
Эта система имеет первый интеграл
й ( 2
ш р,
г) = 2шр,
йш , 2 йрш ,--+ ш
йг йг
= 0
ш2 рт
= Сш. (9)
Следовательно, если масса КА свободна на какой-либо границе (например, в конечной точке), то рт = 0 в этой точке вследствие условий трансверсальности и Ст = 0 (и поэтому рт(1) = 0 вследствие (9)). В этом случае дифференциальные уравнения для х и pv не зависят от т ирт, а дифференциальное уравнение для т не зависит от рт. Таким образом, система (8) может быть разделена на независимые динамическую (уравнения для х и р,) и параметрическую (уравнение для массы КА) части если т свободна по крайней мере в одной граничной точке. Это свойство ОМ-задачи было открыто Ирвингом [1] полвека назад (необходимо подчеркнуть, что в экзотическом случае наличия ограничений на массу КА в обеих граничных точках разделение ОМ-задачи на динамическую и параметрическую части невозможно).
Принимая во внимание это разделение, уравнения для динамической части можно записать в виде:
й 2х
ОТ _0 х
йV_0 Р + р^дп
~ "хх!^ "Г"
(10)
ен2
2 дх
2
а
где pv = |pv|. Уравнение для массы имеет вид:
, 2 2 dm _ Пт Pv
dt ~ 2N '
Рассмотрим некоторое начальное приближе-
ние z0:
(11)
jo
(13)
f(z) = 0,
(16)
f(zo) = b Ф 0,
(17)
и масса КА может быть вычислена из (3) с использованием (4, 7) после решения (10). Большинство авторов используют этот подход. Однако, использование гомотопии между ОМ- и ОТ-задачами требует рассмотрения полной системы (8).
Уравнения (8) должны быть дополнены краевыми условиями. В рамках допущений метода гравитационных сфер нулевой протяженности, типичные начальные условия для гелиоцентрического участка межпланетной траектории имеют вид
x(to) = х>(^), = т(^) = т0 (12)
если гиперболический избыток скорости у планеты отправления Ух0 равен нулю или
и погружение задачи (16) в однопараметрическое семейство:
№)] = (1—т)Ь, (18)
где т е [0; 1] — параметр продолжения.
Дифференцирование (18) по т приводит к следующему дифференциальному уравнению:
(19)
x(t0) = Xo(to), v(to) = vo(to) + Vrt Pvo/PvO , m(t 0) = mo
в противном случае. Здесь x0, v0 — векторы положения и скорости планеты отправления соответственно, t0 — заданное время отправления, m0 — начальная масса КА, pv0 = pv(t0), pv0 = |pv0|. Второе уравнение (13) получено из хорошо известного условия трансверсальности.
Типовые конечные условия соответствуют сопровождению или пролету. Для задачи сопровождения эти условия имеют вид:
x(tf) - х/tf) = 0, v(tf) - v/tf) = 0, (14)
а для пролета граничные условия и условия трансверсальности имеют вид:
x(tf) - х/tf) = 0, pv(tf) = 0. (15)
Здесь xf, vf — векторы положения и скорости планеты прибытия соответсвенно, tf — заданное время прибытия.
3. МЕТОД ПРОДОЛЖЕНИЯ ДЛЯ РЕШЕНИЯ ЗАДАЧИ С ОГРАНИЧЕННОЙ МОЩНОСТЬЮ
Основные результаты этого раздела приведены в работе автора [4].
Принцип максимума сводит ОМ-задачу к КЗ (10, 11), (12) или (13) и (14) или (15). Сначала необходимо решить только динамическую часть ОМ-задачи, следовательно уравнения для m и pm могут быть исключены из рассмотрения.
Конечные граничные условия КЗ могут быть интерпретированы как нелинейная система уравнений
где вектор f составлен из компонентов (14) или (15), а z = (pvo, p vo )T.
— = - V1 ь.
йт \дг!
Начальные условия для (19), учитывая (17), имеют вид:
z(0) = 20. (20)
Интегрируя задачу Коши (19), (20) по т от 0 до 1, можно получить решение нелинейной системы (16). Следовательно, нелинейная система (фактически — КЗ) (16) формально редуцирована к задаче Коши (19, 20). Уравнения (19, 20) реализуют ньтоновскую гомотопию между (17) и (16), т.е. продолжение по конечным краевым условиям (ПККУ) ОМ-задачи.
Динамическая часть ОМ-задачи всегда имеет решение, так как реактивное ускорение не ограничено. Более того, существует множество ОМ-решений для любых граничных условий.
Метод гомотопий (19, 20) не работает, если матрица чувствительности становится сингулярной при каком-либо т е [0; 1
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.