УДК 629.785
ПРИБЛИЖЕННЫЕ МЕТОДЫ РАСЧЕТА ОПТИМАЛЬНЫХ ПЕРЕЛЕТОВ КОСМИЧЕСКИХ АППАРАТОВ С ДВИГАТЕЛЯМИ МАЛОЙ ТЯГИ. ЧАСТЬ I
© 2007 В. В. Салмин, В. В. Васильев, С. А. Ишков, В. А. Романенко, В. О. Соколов, О. Л. Старинова, В. В. Юрин
Самарский государственный аэрокосмический университет
В первой части статьи приводятся математические постановки задач оптимизации перелетов космических аппаратов с двигателями малой тяги и методы их решения. Рассматриваются особенности используемых математических моделей движения для оптимизации управления в рамках различных задач.
Введение
Механика полета с малой тягой (МТ) выделилась в настоящее время в новый раздел механики космического полета, рассматривающий в совокупности проблемы оптимизации траекторий и законов управления движением, а также выбора оптимальных проектных параметров космического аппарата (КА) и его энергодвигательной установки [1, 2].
К настоящему времени в России (а ранее в СССР) и ряде других стран - США, Германии, Франции, Англии, Японии - созданы и испытаны серийные образцы электроракетных двигателей (ЭРД), которые активно используются в космическом пространстве для решения целого ряда практических задач космонавтики. В основном ЭРД использовались в околоземном космосе. Однако уже имеются примеры успешной реализации космических перелетов с малой тягой в дальнем космосе, например, проекты NASA и ESA Deep Space 1 в 1999-2001 гг., SMART 1 в 2003 -2006 гг. (рис. 1.1) и продолжающаяся миссия JAXA Hayabusa 2003 г.
Главным направлением теоретических исследований в течение многих лет в области оптимизации космических перелетов с МТ является развитие аналитических и численных методов решения задач расчета оптимальных траекторий. В последнее время большее значение приобретают вопросы, связанные с учетом дополнительных факторов в математических моделях движения, использованием более полных проектных моделей КА и дополнительных ограничений на воз-
можности управления двигательной установкой (ДУ).
Поэтому вопрос выбора математической модели (или последовательности моделей) для решения вариационных задач приобретает первостепенную важность. Соответственно, проблема оптимизации маневра с МТ не сводится лишь к поиску оптимальных траекторий, а формулируется как проблема совместной оптимизации проектных параметров, траекторий и законов управления движением КА.
В настоящей статье, имеющей две части, представлены в сжатом виде результаты, полученные авторами в ходе многолетних исследований. Более подробное изложение описанных методов и задач содержится в монографиях [16, 25], а также в статьях [6 15, 17 24].
Рис. 1.1. Экспериментальный аппарат NASA Deep Spase 1 с ЭРД
1. Методы и модели оптимизации космических перелетов с малой тягой
1.1. Математические постановки задач оптимизации
Сформулируем проблему совместной оптимизации баллистических параметров и траекторий динамического маневра и проектных параметров КА с двигателем МТ. Под динамическим маневром z из множества маневров Z понимается переход КА из начального состояния х('0) е X0 в конечное x(tK ) е XK . Вектор баллистических параметров маневра Ь включает начальное X0 и конечное XK многообразия в пространстве состояний, внешние условия и ограничения и определяет схему и продолжительность маневра:
г = (......, )Т е Z с Em,
X 0 = X 0( г), ^ = '0 (г), Хк = Хк (г),
'к = 'к( (1.1)
Обозначим символом р вектор проектных параметров р = (р1,р2,...,р1 )Т е Р, соответствующих принятой конструктивно-компоновочной схеме КА. Здесь Р - множество допустимых проектных параметров.
Динамику движения КА будем описывать системой обыкновенных дифференциальных уравнений:
ёх
х = — = /(', х, и, р,г), х(?о) е Хо(г), х(?к) е Хк (г), а
и = и(', х) еи (р), р е Р, £ е Z, г еП(г).
(12)
Здесь х = (х1,х2,...,хп)Т е X - вектор состояния (фазовых координат) системы; и(',х) = (щ ,и2,.....,иг)Т е и(р) - вектор функций управления; и (р) - множество допустимых управлений; иеО(г) - вектор случайных и неопределенных параметров, учитывающий неполноту информации об условиях реализации отдельных маневров; множество О(г) задает априори границы, в которых заключены неопределенности.
Общей задачей совместной оптимизации будем называть задачу отыскания проектных параметров р е Р и совокупности функций (и (', х, г), х (', г)) из множества допустимых Б, обеспечивающих реализацию диапазона динамических маневров Z при минимальном (максимальном) значении заданного критерия эффективности т. Сложность этой задачи состоит в том, что траектории существенно зависят от проектных параметров, и наоборот, параметры КА во многом определяются выбранными траекториями и режимами управления.
Введем интегро-терминальный критерий (функционал) I, зависящий от траектории х('), управления и(', х), баллистических параметров маневра Ь и проектных параметров р, а также неопределенных факторов и:
'к
I[г, р, х('), и (', х), и] = Р [х((0), х('К )] + | /0 (', х, и, и)'
'о
(1.3)
Задачу отыскания экстремума функционала I [г, р, х('), и (', х ),и ] при заданных параметрах Z ир назовем динамической задачей оптимизации. Пусть существует траектория х(') и управление и(',х), доставляющие минимум функционалу
I[г,р,х('),и(', х),и] при фиксированных векторах Z и р для некоторой принятой модели неопределенностей г~ е О:
(х('), и (', х)) = аг§шт I[г, р, х(')и(', х), г~].
х (' )еХ , и (', х)еи (р)
(1.4)
Минимальное значение функционала I, соответствующее этой траектории, будем называть динамической характеристикой
маневра р, г~).
Сформулируем различные постановки задачи оптимизации.
1. Основной задачей оптимизации КА, предназначенного для выполнения единичного маневра с полной информацией
г е Z, назовем задачу отыскания вектора проектных параметров р е Р и вектор функ-
ций (и(t, х, z), x(t, z))e D , доставляющих максимум критерию //:
(х,U,p) = argmax /(z,p,x(t),u(t,х)).(1.5)
(x,u )eD (p), peP
2. Вектор параметров p e P оптимален для диапазона динамических маневров Z, если КА с параметрами p может выполнить любой маневр из заданного диапазона Z и максимальная степень неоптимальности p(z,p) на множестве Z достигает минимального значения при p = p:
p = arg min max p(z, p) .
peP
zeZ
выполняемых независимо друг от друга операций:
1) (x, й) = arg min I[z, x, и], I(z, x, й) = S(z),
(x,u )eD
2) p = argmax/(z,p,S(z)) .
pe P
(1.6) (1.7)
Под степенью неоптимальности р(г, р) понимается мера проигрыша в критерии эффективности /т(г, Р) при отклонении вектора проектных параметров р(г) от оптимального значения р (г) . Вектор р назовем вектором универсальных для множества 2 проектных параметров.
3. Задачу о максимуме критерия оптимальности ¡¡(г, р, х, и) назовем разделяющейся на динамическую и параметрическую, если в критерии 1 удается выделить критерий низшего уровня - функционал I, зависящий только от траекторий и управления и не зависящий от проектных параметров. Минимум I для каждого фиксированного маневра достигается на паре (х (^), и (¿, х)) е Б и обеспечивает локальный максимум критерия ¡ при любом выборе вектора параметров р е Р: 1 = ¡(г, р, I[г, х, и]). Если в соответствии со сказанным выше
min I[z, x(t), u(t, x)] ° S(z) " p e P,
тогда
//(z,p,x,u) ° //(z, p, S(z)).
Очевидно, если < 0, то решение задачи оптимизации реализуется в форме двух
4. Задачу о максимуме критерия m(z, p, x,u) будем называть условно разделяющейся на динамическую и параметрическую части, если минимум функционала I(z, p, x, и) обеспечивает локальный максимум критерия т. Отыскание глобального максимума т реализуется в форме двух последовательных операций:
1) (x, й) = arg min I[z, p, x, и],
( x ,и )eD
I (z, p, x, й) = S (z, p), (1.8)
2) р (г) = а^шах ¡(г, р, £(г, р)) . (1.9)
Решение этой задачи связано с необходимостью иметь зависимость £(г,р), определенную на множестве Р во всем диапазоне маневров Z.
1.2. Методы решения задач оптимального управления
В задачах оптимального управления элементами класса допустимых Б являются управляемые процессы, точнее их математические модели. Рассмотрим систему, которая в каждый момент времени характеризуется
вектором состояния х = (х1,х2,...,хп)Т , являющимся элементом некоторого множества X, называемого пространством состояний. Изменение состояния х во времени называется процессом. Управляемые процессы принято описывать путем указания закономерности перехода от предыдущего состояния к последующему в зависимости от управляющего воздействия, которое характеризуется вектором управлений и = (и1,и2,...,иг)Т , являющимся элементом некоторого множества и (множества управлений).
Математическая модель управляемого процесса представляет собой, как правило, уравнение, связывающее последующее состояние с предыдущим состоянием и управлением.
Следует подчеркнуть, однако, формальный смысл понятий состояния и управления, связанный с принятой формой математической модели, которая для реального объекта может быть не единственной. Например, в задачах оптимизации космических траекторий в качестве состояния принято рассматривать положение и скорость центра масс КА, а в качестве управления - направление вектора тяги. Изменение направления вектора тяги часто достигается путем его поворота с помощью устройств, создающих моменты относительно центра масс (ЦМ), что заставляет учитывать и динамику вращательного движения КА. Это несоответствие отражает тот факт, что математическая модель зависит от тех задач, которые должны решаться с ее помощью.
Пусть на отрезке [t0, tK ] задается множество D как множество кусочно-непрерывных функций n(t) и кусочно-дифференцируемых x(t), удовлетворяющих условиям (1.2) и минимизирующих функционал (1.3). Общий подход к решению задачи оптимизации в постановке (1.2, 1.3) основан на использовании принципа расширения и достаточных условий оптимальности В.Ф. Кротова [3, 4].
Как правило, во многих задачах оптимизации множество допустимых D задается посредством некоторых условий, выделяющих его из более широкого множества E. Принцип расширения множества допустимых состояний и управлений состоит в том, что функционал доопределяется на более широкое множество E', но так, что наименьшее значение он принимает в D. Для того, чтобы функционал (1.3) достигал абсолютного минимума на (х, u) е D, достаточно существования такой непрерывной и дифференцируемой функции j (t, x), чтобы
1. R(x, u, t)= max R(x, u, t),
(x,u )eE (t)
2. G (x (t о), x (tk )) =
min
x (t0 ^x (t0
x (tk FEx (tk
G (x(t о), x(tk)),
(110)
„ n j r dj где R = I f, +,
i = 1
dx.
dt
С(х0,хк) = Р(■) + (р('к,хк) - у('0х). (1.11)
Необходимо задать синтезирующую функцию р(', х), которая доопределяет функционал I на Е^П так, чтобы минимум I принадлежал Б. Один из способов (формализмов) приводит к процедуре принципа максимума Понтрягина, другой - к процедуре динамического программирования, а третий (метод кратных максимумов) пригоден для решения так называемых вырожденных задач.
1.3. Метод усреднения в задачах оптимального управления
При решении некоторых видов динамических задач оптимизации в механике полета с МТ используются методы асимптотического разделения параметров движения на быстрые и медленные компоненты [2]. Это обусловливается, во-первых, наличием в явном виде малого параметра - реактивного ускорения от тяги, которое меньше гравитационного на несколько порядков; во вторых, присутствием циклической переменной - угловой координаты.
Запишем уравнения движения в общем виде, придерживаясь общепринятых в задачах такого рода обозначений:
х = а ■ X (х, р, и (')),
р = т(х,р) + а ■ У (х,р,и (')), (112)
где х - вектор «медленных» переменных размерности п; а - малый параметр; р - быстрая скалярная переменная (фаза); и(')е и -вектор управлений размерности г.
В общем случае в управление могут входить как быстрые, так и медленные составляющие. Поэтому задачу выбора оптимального управления удобно разделить на две: определение управления как функции
быстрой переменной (выбор структуры управления на витке) и определение законов изменения параметров этой программы от витка к витку.
Критерий оптимальности представим в следующем виде:
I = а | (х, и) Ж.
(113)
0 0
а г г I Л Гл л
3 = 2р-» Р° х, и (х,у,р ]рр Ж. (1.17)
Интегралы в правых частях системы (1.16) образуют совокупность так называемых «усредняющих интегралов»:
2р _
ф= |(.)р , ф = (ф . )Т, у = 1,2п +1. (1.18)
Перейдем в системе (1.12) от времени I к быстрой переменной р :
4х = аХ (х ,р, и(( ))• о- (х ,р),
ар
— = а • О— (х ,р)
ар
(1.14)
Здесь т - так называемое «медленное» время, т = а • ^. В соответствии с принципом максимума Понтрягина введем вектор сопряженных переменных ¥ и запишем Гамильтониан системы (1.14):
После усреднения правые части системы (1.16) не содержат циклической переменной р. Поэтому модель «медленной» эволюции вектора состояния и вектора сопряженных переменных может быть представлена в виде системы интегро-дифференциальных уравнений с «медленным» временем в качестве независимой переменной:
ат
ф,
ф
, где д = 1 х,у
Н = уТ (аа 1 (х,р )• Х(х,р )))+ +уГ • аа- (х,р) - аа- (х,р) • Р0 (х,р,и) =
= а¥ (х,р ,у, и ).
(1.15)
Определим локально-оптимальное управление и из условия максимума Гамильтониана на отрезке р е [0,2р]. Проведем затем процедуру усреднения исходной и сопряженной систем по быстро меняющейся переменной р. Усредненная система уравнений будет иметь вид
| а- (х, р) • X | х, р, и [ х,у,р | р
а х _ а ар 2к
Л
р=2Р 0 * [х, *и [ху,р ))ар,
р=2Р Iа -1 [х ,р ]ар,
, = 1,2п +1, ф 2
= а
х,рур. (1.19)
Исходная оптимизационная задача сводится, таким образом, к решению краевой задачи для системы (1.16, 1.17). Однако вычисление усредняющих интегралов (1.18) представляет самостоятельную проблему и требует разработки специальных процедур.
1.4. Приближенный метод решения динамической задачи (методразбиения)
Пусть движение КА описывается системой дифференциальных уравнений (1.2). Откажемся от получения универсального решения для всего пространства переменных и поставим цель определить ряд упрощенных решений для каждой отдельной выделенной области X пространства состояний. Разобьем допустимую область фазового пространства переменных на т подобластей таких, что
X е Х1 и Х2 и Х3 и... и Хт . Заменим исходную двухточечную краевую задачу на многоэтапную последовательность переходов:
а усредненный критерий оптимальности
х0 ® х1 ® х2 ® ••• ® хт-1 ® хК'
0
Л
0
Л
где хк = {х1 ,х2,...,хт-1} - граничные условия промежуточных (нефиксированных) состояний:
х1 е X1 п X2, х2 е X2 п Xз
(120)
хт-1 е Xm-\ П Xm.
Таким образом, сложный динамический маневр представим в виде последовательности маневров с параметрами
= {х0 ,х1} , г2 = {х1 >х2 } '
гт = {хт-1 ,хК } .
(121)
Функция управления для единичного маневра г. определяется из условия минимума функционала
= а^ тт 3.(г. х,и ),
( и ,х )е
(122)
где Б. си 'X- допустимая область управлений и состояний.
В результате получается динамическая характеристика единичного маневра
= тт 3(г^х,и). (и,х) е Д
(123)
Для получения аналитических решений данной задачи по участкам в зависимости от
вида выделенной области X. возможны различные подходы.
1.5. Метод последовательных расширений
Метод последовательных расширений предполагает поэтапную редукцию математической модели задачи оптимизации (временное отбрасывание связей и ограничений; усреднение движений, носящих циклический характер; линеаризацию движения в окрестности опорной орбиты и т. д.); получение последовательности моделей различного уровня сложности, сохраняющих, однако, все важ-
нейшие особенности исходной модели. Затем определяется структура оптимального управления в рамках упрощенной модели, а также приближенных аналитических зависимостей критерия оптимальности (функционала) вариационной задачи от граничных условий динамического маневра; синтезирующая функция и ее частные производные по компонентам вектора состояния в аналитическом виде; проводится построение оценочной функции режимов управления.
В первом приближении выбор проектных и баллистических параметров, а также траекторий и управлений осуществляется с использованием простейших проектных и динамических моделей, например, представляющих КА точкой переменной массы с «бесплатным» управлением. Последующие приближения используют более сложные модели, учитывающие, например, динамику углового движения КА, возмущающие ускорения от гравитационных, аэродинамических и иных сил.
В результате реализуется схема, основанная на использовании последовательности усложняющихся моделей управляемого движения: в качестве первого приближения используются приближенно оптимальные решения для упрощенных моделей, приближенный синтез управления; затем модель динамической задачи последовательно усложняется, а решение, полученное на предыдущей итерации, используется для оценок и сравнения различных режимов управления.
1.6. Итерационная процедура поэтапной оптимизации
Поскольку проектные параметры КА влияют на динамическую характеристику маневра, и наоборот, баллистическая схема и траектории перелета существенным образом влияют на выбор проектных параметров, процесс оптимизации параметров КА и семейства оптимальных траекторий необходимо вести совместно.
Ключевая идея предлагаемого подхода состоит в условном разделении общей проблемы оптимизации на динамическую и параметрическую части с последующим их объединением через динамическую характеристику маневра, являющуюся мерой затрат
на его реализацию, зависящую от граничных условий, проектных параметров и неопределенных факторов, и уточняемую по мере перехода от простых моделей динамической задачи к более сложным.
Зададим последовательность математических моделей {м5}, .=1, 2,... динамической задачи для конкретного маневра 2 из подмножества Z. В рамках каждой из моделей М* определим критерий оптимальности динамической задачи - функционал 'г,р,х,и,у), а также множество допустимых траекторий и управлений П. и получим динамическую характеристику маневра г,р,г>*.). Пусть в результате решения совокупности динамических задач с применением моделей {м5 }
получена последовательность значений глобального критерия оптимальности
{а}={а(*)(2,р,5(*)(2,р,у.*))}, * = 1,2,... и определен вектор оптимальных проектных параметров согласно выражению
= а^шах А(г,р,8(*>(г,р,ь*.))
рер
Процесс оптимального синтеза назовем устойчивым, если сколь угодно малым приращениям вектора проектных параметров соответствуют малые изменения критерия ^. Процесс заканчивается, когда применение модели более высокого уровня не приводит к заметному изменению критерия оптимальности ^ и вектора проектных параметров р .
В качестве первого приближения используются приближенно оптимальные решения для упрощенных моделей, строится приближенный синтез управления. Затем модель динамической задачи последовательно усложняется.
Влияние неконтролируемых факторов приводит к неопределенности динамической характеристики маневра в пределах нижних
и верхних границ ее изменения 5Н , 5В, которые в свою очередь определяются размерами области неопределенности О. Сначала анализируются пределы изменения динамической характеристики
5 (2, р) £ 5(2,ри) £ (2,р)
гг о
и выбираются проектные решения /-го приближения pi -1, p2), соответствующие предельным оценкам SH , SB :
p(,) = arg max ß(z,p,Sjj(z,p))
1 p e P H
p2(,) = arg max m(z,p,SB(z,p)) . (1.24)
2 p e P
Соответственно, будем иметь критерии оптимальности для двух вариантов решений:
m(V = m(z,pг(,)), m2(,) = m(z,p(2,}).
(1.25)
Сравнение компонентов векторов p\' ^,
p2 ^ и критериев m1(,), mi г) позволяет установить влияние неопределенности на облик проектируемого КА и показатель его эффективности.
На следующих итерациях уточнение проектных параметров приводит к необходимости повторного расчета семейства оптимальных траекторий и режимов управления, а также баллистических параметров.
1.7. Проектная модель космического
аппарата с электрореактивной двигательной установкой малой тяги Для выбора оптимальных проектных параметров КА представим его стартовую массу как сумму масс отдельных систем. Анализ работ в области оптимизации КА с ЭРД малой тяги позволяет выделить следующие элементы: 1) энергоустановку, состоящую из источника и преобразователя энергии; 2) ДУ, включающую маршевые и управляющие двигатели вместе с исполнительными органами (кардановым подвесом); 3) рабочее тело, необходимое для осуществления перелета с учетом затрат на управление; 4) систему подачи и хранения рабочего тела (баки, трубопроводы и др.); 5) полезный груз; 6) корпус и прочие элементы конструкции.
Уравнение баланса масс на начальной орбите имеет вид
м 0 = м ПГ + мэ + М Д + М СПХ + МТ + М К,
где М0 - начальная масса, МПГ - масса полезного груза, МЭ - масса источника и преобразователя энергии, Мд - масса ДУ МСПХ -масса системы подачи и хранения рабочего тела, МТ - масса рабочего тела, МК - масса корпуса, прочих элементов и систем.
Определим вектор параметров КА:
Р (Р' РУПР ' N' ГУПР' rM '10 ) ,
(1.27)
где Р - тяга маршевых двигателей, РУПР - тяга управляющих двигателей, N - мощность энергоустановки, гУПР и гМ - векторы, характеризующие расположение относительно центра масс точек крепления управляющих
и маршевых двигателей, ^^ - тензор инерции. Компонентами этого вектора р являются параметры, наиболее полно характеризующие КА, его схему управления, массу, компоновку, энергетическую и двигательную установки и пр.
Массы отдельных компонентов КА зависят от проектных параметров. Обычно применяются следующие зависимости.
Мэ = aэyN, МД = Уду (Р + кРУПР ), (1.28)
МСПХ = УспхМт , МК = /КР + /КN ,
где аэу, Уду , Успх , /К , /К - соответствующие удельные массовые характеристики. Мощность энергоустановки зависит от тяги двигателей и скорости истечения рабочего тела
N =
Ре 1 + ж
2 hrh
(1.29)
T' 1П.Э
где Ж =
Р
УПР
Р
характеризует относительным
расход массы управляющих двигателей, Т]Т -тяговый коэффициент полезного действия, Т]ПЭ - КПД преобразователя энергии.
Введем в рассмотрение вектор управлений: и = ^Т, еУПР, 5УПР
,МУПР} , где е и вупр -направления маршевого и управляющего век-
тора тяги, дУПР и МУПР - функция включения-выключения и величина управляющего момента, соответственно. При жестком креплении двигателей и = {)УПР ,МУПР } е и .
Основная задача оптимизации формулируется следующим образом: определить из допустимого множества Р вектор проектных параметров р и вектор функций управления и (}) е и, доставляющие при заданной массе полезного груза и заданном времени перелета Т минимум начальной массе КА при выполнении граничных условий переходов:
М 0 = ^ М 0 (р,и (),Т,МПГ ) . (1.30)
Для выделения динамической задачи оптимизации введем в рассмотрение «приведенную» характеристическую скорость — меру энергетических затрат на управление движением - динамическую характеристику маневра:
Р
VX = f—
X JM
1 + qynp 15,
Л
q
УПР\
dt,
(1.31)
где д, дУПР - соответственно секундные расходы рабочего тела маршевых и управляющих двигателей.
Проблема оптимизации разделяется на две независимые:
1) динамическую - нахождение оптимальных программ управления:
u
opt1
(t, p) = arg min Vk (u,p,Xo ,xK ),
u
opt 2
(t,p) = arg minVxK 2 (u,p,Xo ,Xk ) (1.32)
и получение динамическом характеристики перелета:
VXK (u opt (t ),p , x0 ,XK) = Vxk (p
, X0 ,XK
); (1.33)
2) параметрическую - нахождение оптимальных проектных параметров КА:
г opt
= argminMo[Vxk (p , X0 , XK ),p,T,M пг J.
peP
ueU
ueU
1.8. Оценки перелетов с малой тягой на расширенном множестве допустимых траекторий и управлений
Рассмотрим задачу перехода КА с двигателем МТ из одной точки пространства
(r0,V0) в другую точку (rK,VK ) за заданное
время T = tK -t0. Пусть требуется обеспечить минимум некоторого функционала I. Полагаем, что уравнения движения заданы в форме векторных дифференциальных уравнений, описывающих движение ЦМ и угловое движение КА, причем граничные условия для ориентации КА и для угловой скорости не заданы.
Нетрудно проследить связь данной постановки задачи с другими, которые уже рассматривались в литературе. Исключив из системы, описывающей движение КА, уравнения, описывающие его угловое движение, получим расширение множества допустимых траекторий и управлений D до некоторого множества E (поскольку исключены некоторые связи). Очевидно, что min I < min I.
ED
Однако ясно, что случай, когда направление тяги является независимым управлением, а не фазовой координатой, соответствует традиционной постановке задач оптимизации траекторий в механике полета, когда не учитывается угловое движение КА, а вектор тяги может произвольно менять свою ориентацию в пространстве.
Таким образом, данное соотношение служит априорной оценкой функционала на расширенном множестве допустимых траекторий и управлений E. Такая оценка может оказаться полезной и содержательной в случае, когда на управляющий момент не накладывается слишком жестких ограничений. В противном случае степень неоптимальности управляемого движения может оказаться завышенной. Однако отыскание mEin I требует
точного решения задачи оптимального управления движением ЦМ, которая сама по себе является достаточно сложной.
При отсутствии ограничений на направление вектора тяги e (t) в связанной системе
координат (СК) ориентация вектора тяги в пространстве не зависит от углового положения КА, и управление траекторией осуществляется независимо от его движения относительно ЦМ.
Если направление тяги фиксировано в связанной СК (двигатель жестко закреплен
относительно корпуса КА), то I — I ° 0 .
dt
СВ
В этом случае изменение направления вектора тяги в пространстве осуществляется только за счет разворота корпуса КА. В этом случае удобно считать, что тяга направлена вдоль одной из связанных осей, например,
ОХ1 (е ° /' ) . Тогда — = о х е .
14 Ж
Если вектор тяги сохраняет неизменное положение в неподвижной СК, получим вы-
ражение
de dt
= -w х e, определяющее ки-
нематику программного разворота тяги относительно корпуса.
1.9. Математическая модель для оптимизации совместного управления траекторным и угловым движением Общей тенденцией развития перспективных КА с ЭРД является увеличение их масс и моментов инерции, что, в свою очередь, создает ряд проблем управления движением. Конструктивные схемы тяжелых КА, как правило, предусматривают жесткое закрепление связки двигателей относительно корпуса. Изменение направления тяги реализуется при этом путем разворота корпуса КА в пространстве с помощью управляющего момента, величина которого ограничена.
Если момент создается самими маршевыми двигателями, то расход рабочего тела на управление отсутствует. Назовем такую схему управления траекторным и угловым движением совместной.
Раздельной схемой управления будем называть такую схему, которая предполагает использование специальных двигателей, создающих только момент относительно центра масс. Это связано с дополнительными затратами рабочего тела.
Для задач совместного управления тра-екторным и угловым движением КА используем модель, учитывающую динамику движения относительно центра масс, ограничения на ориентацию вектора тяги, зависимость тяги двигателя от расстояния КА до Солнца и от ориентации солнечных батарей, влияние несферичности Земли и сопротивления верхних слоев атмосферы, влияние гравитационных полей Солнца и планет и другие факторы:
dr fy dV _ _ 7 Pc__ 7 — = V ; — = a + g + f = — 8 • e + g + f, dt dt m
da dt
Io1 (Mо - a x 10a);
di _ t d/1 _ - dk _ -
— = a xi. • — = a x —1 = a x• (135)
1 1 1
dt
dm dt
dt
(q+qупр)
p
dt
8 + q 8 I
упрупр
где r - вектор положения ЦМ; a , g, f - векторы реактивного, гравитационного и возмущающих ускорений; P = c • q - тяга маршевого двигателя; e - единичный вектор направления тяги; 8 - функция включения-выключения и реверса тяги маршевого двигателя;
qynp и 8упр - секундный расход и функция nY = cos уСБ sinj включения-выключения управляющего двигателя; а - вектор угловой скорости; (ij, j1, k1) - единичные векторы вдоль связанных осей; I0 = I0 [m(t )] = I0 (t) - матрица тензора инерции. Главный момент внешних сил M0 представлен в виде суммы:
1.10. Математическая модель для совместной оптимизации траекторий и ориентации солнечных батарей
Электрическая мощность, вырабатываемая солнечными батареями (СБ) на освещенных участках траектории, зависит от угла / между направлением на Солнце и нормалью к поверхности батарей: N = Мтахео$р. Задачей управления ориентацией СБ является обеспечение максимального значения сояр.
Положение СБ относительно корпуса КА будем характеризовать двумя углами (рис. 1.2): усв - угол крена оси батареи, который составлен осью вращения батареи 01СБ и поперечной осью КА 02р - угол собственного вращения батареи, характеризующий поворот нормали 0УСБ к плоскости батареи вокруг ее собственной оси О1 СБ. Очевидно, с помощью последовательных поворотов на углы усъ и можно добиться постоянного направления нормали УСБ на Солнце. При этом следует учитывать, что КА одновременно осуществляет программу разворота по углу у.
Запишем выражения для проекций единичного вектора ОУнормали к плоскости СБ на оси орбитальной СК 0ХУ1:
nX = cos (СБ cos y - sin уСБ sin (рСБ sin y,
" СБ ■
(1.37)
м о = мУПР + Н,
где МУПР - вектор управляющего момента, Н - вектор момента от гравитационных, аэродинамических и иных внешних сил.
Необходимые условия реализации программных разворотов КА записываются в виде
nZ = cos (СБ sin y + sin уСБ sin (СБ cos y .
Будем считать, что в каждый момент известны компоненты единичного вектора
rS = (rSX, rSY, rSZ) направления на Солнце в орбитальной СК. Максимальная мощность реализуется при cosft = n-r = 1. Для этого случая выражения для программных углов ориентации солнечных батарей имеют вид:
(СБ = arccos(rSX siny + rSZ cosy), arccos r„
уСБ =
' SX /
sin jСБ
(1.38)
- d a _ _
МУПР > I0--H + a x I a.
(1.36)
Однако двухканальное управление ориентацией солнечных батарей в сочетании со
сложной программой изменения положения корпуса КА в пространстве может оказаться трудным для реализации. В этом случае рассматриваются альтернативные варианты од-ноканального управления: 1) батареи вращаются только вокруг оси OZQE, постоянно совпадающей с поперечной осью 0Z1 (gcB = 0); в этом случае угол рСБ обозначим через рр 2) ось вращения батарей совпадает с направлением связанной оси ОУ а значит, и оси
OY; уСБ ° л/2; в этом случае угол собственного вращения обозначим через р Очевидно, для этих двух вариантов cos ¡3j jj < 1.
Получим выражение для cos ¡j и для этого в (1.38) положим ус& = 0. Для обеспечения максимума cos Д найдем оптимальное значение рг
'SY
(139)
jjopt = arctg-
rsx cosy + rsz smy ' При этом
cos ¡¡jmax = V TIy + (rsx cos У + rsZ slf> У )2 .
(1.40)
Для второй схемы управления cos bjj не зависит от углового положения КА. Положим уСБ ° л/2, и тогда
jjj = arctg^-y, cosbn =Vr2Sx + 4 .(1.41) r
'SX
С учетом возможности выключения двигателя при попадании аппарата в тень Земли целесообразно анализировать поведение среднего за виток косинуса угла ¡3. Освещенность КА зависит от взаимного положения Солнца и оскулирующей плоскости орбиты. Поскольку оно меняется в процессе полета, на отдельных этапах перелета может оказаться более выгодной первая схема од-ноканального управления СБ, а на других -вторая. Поскольку целью управления СБ является обеспечение максимальной мощности, можно рассмотреть также комбинированную схему, при которой возможны развороты аппарата по крену на ±90° при смене знака разности (cos ¡j - cos ¡3jj).
1.11. Математическая модель для оптимального управления околоземными орбитами на больших интервалах времени
Характерной особенностью задач уп -равления движением КА на низкой околоземной орбите является наличие возмущающих ускорений, обусловленных нецентральностью гравитационного поля Земли и сопротивлением верхних слоев атмосферы, сравнимых с величиной реактивного ускорения. Модель задачи оптимизации становится при этом достаточно сложной, а эллиптичность орбиты требует аккуратного описания «медленной» эволюции орбиты на больших интервалах времени. В этих задачах на первый план выходит стратегия гарантированного результата как при выборе законов управле-
ния, так и при оптимизации параметров корректирующей ДУ.
Для задачи совместного управления оскулирующими элементами орбиты (А, е, щ О, /), а также относительным угловым положением Аи КА, обеспечивающего минимум характеристической скорости, с использованием метода усреднения была получена структура оптимального управления на отдельном витке орбиты, показанная на рис. 1.3. Здесь ц - эксцентрическая аномалия центра
разгонного участка (ах > 0), X - половина ширины разгонного участка для трансвер-сальной тяги, а - ширина одного пассивного участка трансверсальной тяги; £ - аргумент
широты центра участка с аг > 0, р - половина ширины рабочего участка с аг > 0 , / - ширина одного пассивного участка для бинормальной составляющей тяги.
Рис. 1.3. Полученная структура оптимального управления на витке
Для этой структуры управления получена математическая модель эволюции орбиты в поле земного сфероида с учетом возмущающего влияния атмосферы и малой тяги с оптимальной структурой управления:
dA _ 4a1 dt p \
A3
—(x+)-2oPcpjmA,
— _ — A [- 3e(X + ) + 4 sin(X + f) cos f cosh -dt p\f
- ff sin f(X + f) cos a cos 2h ] - 2e op^J ff,
dw _ a1
— [4 sin(X + -f)cos f sinh -
f
dt pe \
e . „/> a\ „ i e(5cos2 i -1) —sin2(X + f)cosasin2hJ + —-0 5 3 5 ,
2
d&u dt
2f0'5A '5
_4f (A-15 - Ax15)
(1.42)
di a3 A dt p \ f
о ■ 1 b ö b у « f b p ö 2 sin| j + 2 I cos 2 cos z - 3/jl j + 2 - - I-
--3-sin2^ j + -20cosb• (l cos2Z + 1 sin2£)
dW a3
dt p siniy f
_ з 2
O l b ö b ■ r f b p ö
2 sin| j + у I cos у sinz - 3l21 j + 2 - -I-
3sin 2fj + -2j cos b • (l cos2Z + sin 2Z )
ecosi -1
" f0'5A3'5 '
где a1, a3 - составляющие ускорения, направленные по трансверсали и бинормали, соответственно.
1.12. Математическая модель для оптимизации управления относительным движением КА В качестве основной будем рассматривать схему управления относительным движением двух КА. Один из них считается пассивным (КА1), а другой активным (KAII), снабженным ЭРД.
Рассмотрим возмущенное движение КА в цилиндрической системе координат ruz, где r - расстояние от центра Земли до проекции КА на плоскость невозмущенной круговой орбиты, u - угол, отсчитываемый в плоско -сти невозмущенной орбиты от некоторой начальной оси по направлению полета спутника, z - расстояние от плоскости невозмущен-
ной орбиты до КА. Считая, что величина
r
мала, запишем уравнения движения в виде:
z
ёт _ ^ ёи _ Уи ёг _ у
& г' & т ' Ш 2'
У _ у2 _ т+^
& т т2
зу, т
—^ г + Ж Ж Г
УУ
' ж
+ т, (1.43)
Здесь 8, Т, Ж - проекции возмущающих и управляющих ускорений на оси орбитальной
СК, Уг - радиальная скорость, Уи - трансвер-
сальная скорость, У, - нормальная скорость (проекция скорости на перпендикуляр к плоскости невозмущенной орбиты), ¡1 - гравитационный параметр, I - текущее время.
В большинстве практических задач эксцентриситет опорной орбиты невелик, и поэтому уравнения относительного движения записываются в следующем виде:
Ат _ АУг,
АЬ _АУи _ 1Ат,
АУ _ 2ЛАУ, _ 12 Ат + 8,
(144)
АУ, _ _1АУ + Т,
Аг _АУ,
АУ, __Я2Аг + Ж
Здесь 1 _
средняя угловая ско-
рость движения КА1 по опорной орбите;
АЬ _ Аи • т - проекция расстояния между КА на дугу опорной орбиты.
1.13. Математическая модель для оптимизации перелетов между орбитами с большими эксцентриситетами
Для задач оптимизации перелетов между орбитами с большими эксцентриситетами можно рассматривать два варианта ориентации вектора тяги: свободная ориентация и ориентация по трансверсали.
Изменение оскулирующих элементов кеплеровской орбиты описывалось с использованием усредненных уравнений, полученных на основе стандартной процедуры усреднения уравнений в оскулирующих элементах для плоского движения КА:
ёА _ 1
Зу ~ Р
ёе 1
ЗУ ~ 2Р
• ^ +лЯ_е2 • J2),
л/Ат^г
_ е2 х
л/1 _ е2 • J1 + 2 J3 _ е •(J4 _ J2)
ёт 1
ЗУ„ 2рА
_ е2 х
_ J5 + eJ6 + +J7 _
¿.к.
JJ ф} (®,Е)ЗЕ, ]=18.
л/^
(145)
Здесь ¡ - гравитационный параметр Земли,
в - угол, характеризующий ориентацию тяги в плоскости орбиты относительно транс-версали, Е - эксцентрическая аномалия,
J1,...,J8 - усредняющие интегралы - функции
параметров управления.
1.14. Модели для оптимизации межпланетных перелетов с малой тягой
Граничные условия межпланетного перелета определяются его целью и относительными положениями планет старта, финиша и КА. Обычно траектория движения разбивается на участки движения в сферах действия планет и Солнца и оптимальное движение рассчитывается по участкам. На границах участков необходимо осуществлять стыковку траектории по фазовым координатам и массе КА.
Особенностью оптимизации замкнутых межпланетных перелетов (с возвращением КА на планету старта) является дополнительное условие равенства угловых перемещений аппарата и планеты старта в конечный момент времени:
(Т + Т )• тЗ _ (<?2 + )+ Т (тЗ _ тМ ) _ 2Р • П .
(146)
т
х
е
х
2
3
Здесь (р2 и (р4 - угловые дальности прямого и
обратного гелиоцентрических перелетов, тЗ
и тм - средние угловые скорости движения Земли и Марса, п - произвольное целое число. Появляется неоднозначность решения целевой задачи в зависимости от г _ Т2 /Т4 -соотношения длительностей прямого и обратного перелетов и Ц0 - даты старта. Это приводит к необходимости введения и последующей оптимизации дополнительных параметров, описывающих баллистическую
схему перелета Ь _ Ц0, г }Т .
Задача проектно-баллистической оптимизации межпланетного перелета формулируется следующим образом. Требуется определить вектор проектных параметров
р _ {Р0, с} е Р , вектор баллистических параметров Ь _ Ц0, г} е В и вектор функций
управления и ^) _ {е (/), д()}Т е и, доставляющие при заданных массе полезного груза МПГ и длительности перелета Т минимум стартовой массе КА и обеспечивающие выполнение целевой задачи, описываемой множеством допустимых фазовых координат аппарата X:
г; = { c (R )• s dt,
(1.48)
где с (К) - коэффициент, учитывающий падение мощности энергоустановки, а следовательно, тяги двигателей и расхода рабочего тела в зависимости от расстояния от КА до Солнца. Для КА с ядерными источниками
энергии с (К) ° 1, для КА с солнечной энергоустановкой
c (r )= Р® =NfiRLJ
p0
Nn
R
(1.49)
M0 = Min M0(p, b, u(tjMnr = fixe T = fixe>x E
(147)
Задачи оптимизации пилотируемых
экспедиций наиболее сложны, так как множество допустимых фазовых координат кроме граничных условий прямого и обратного перелетов содержат специфические ограничения, связанные с обеспечением безопасности экипажа (ограничения на суммарную длительность экспедиции TS £ ТПРЕд , минимально-допустимое расстояние от КА до Солнца
R| > Rni>m, длительность нахождения в радиационных поясах Земли и др).
Для разделения задачи оптимизации на параметрическую и динамическую части вводится промежуточный критерий оптимизации - приведенное моторное время:
Этот критерий непосредственно определяет суммарные затраты рабочего тела
Р
МТ (р) _ — • Т*(р) на перелет и, следователь-
с
но, является динамической характеристикой маневра.
Баллистическая часть задачи оптимизации состоит в определении вектора функций управления и ()_ {е (), д()}Т и вектора баллистических параметров Ь _ Ц0, г }Т (для замкнутых перелетов), обеспечивающих выполнение целевой задачи с минимальными затратами рабочего тела при фиксированных проектных параметрах КА, и построении зависимости
Т*(р)_ Ь й |р _ Ахе Ахе х е х).
Проектная часть задачи оптимизации состоит в выборе вектора проектных
параметров р _ {Р0, с}Т, обеспечивающих
минимум стартовой массе КА с учетом полученной зависимости.
Баллистическая часть задачи оптимизации решается в соответствии с разработанным подходом, связанным с использованием последовательности усложняющихся моделей.
Модель А описывает движение аппарата в рамках теории сфер действия в центральном поле притяжения Солнца и планет без учета возмущений от других притягива-
0
ющих центров в плоской полярной СК. Орбиты планет считаются круговыми и компланарными. Стыковка плането- и гелиоцентрических участков осуществляется только по массе КА.
В рамках модели Б движение КА разделяется, в соответствии с теорией сфер действия, на гелио- и планетоцентрические участки, и при этом на границах сфер действия проводится точная стыковка траекторий движения по координатам, скорости и массе КА. Орбиты планет считаются эллиптическими и некомпланарными. При расчете движения в сферах действия планет учитываются участки затенения, гравитационные возмущения от других небесных тел и нецентральности гравитационного поля планеты.
Модель В использует уравнения движения в поле притяжения нескольких тел (Солнце и планеты солнечной системы), учитывается эллиптичность и некомпланарность орбит планет, деградация СБ и другие факторы. При решении динамической части задачи учитываются ограничения на проектные и баллистические параметры и вектор управления, траектория рассматривается как непрерывная с оптимальной стыковкой участков.
Список литературы
1. Гродзовский Г. Л., Иванов Ю. Н., То -карев В. В. Механика космического полета (проблемы оптимизации). - М.: Наука, 1975.
2. Лебедев В. Н. Расчет движения космического аппарата с малой тягой.- М.: ВЦ АН СССР, 1968.
3. Гурман В. И. Принцип расширения в задачах управления. - М.: Наука, 1985.
4. Кротов В. Ф., Гурман В. И. Методы и задачи оптимального управления. - М.: Наука, 1973.
5. Гурман В. М., Попов Ю. Б., Салми-н В. В. О возможности реализации траекторий аппаратов с малой тягой с учетом их движения вокруг центра масс //Космические исследования. - 1970. Т.8, № 5. - С. 684-692.
6. Салмин В. В. Оптимизация режимов разгона вращающегося космического аппарата с двигателем малой тяги //Космические исследования. - 1973. Т. 11, № 8. - С. 842 - 853.
7. Брусов В. С., Салмин В. В. Комбинированная двигательная система, универсаль-
ная для диапазона маневров //Космические исследования. - 1974. Т. 12, № 3. - С. 368 - 373.
8. Васильев В. В., Салмин В. В. Оптимальный разгон космического аппарата с электрореактивным двигателем при ограниченной скорости поворота вектора тяги //Космические исследования. - 1976. Т. 14, № 3. -С.336 - 342.
9. Салмин В. В. Многошаговые алгоритмы управления движением космических аппаратов //Космические исследования. - 1979. Т. 17, № 6. - С. 835 - 845.
10. Салмин В. В. Аналитическая оценка оптимальности многошаговых адаптивных алгоритмов управления //Космические исследования. - 1980. Т. 18, № 3. - С. 332 - 342.
11. Васильев В. В. Оптимальное управление эллиптической орбитой спутника Земли с двигателем малой тяги //Космические исследования. - 1980. Т. 18, №5. - С. 707 - 714.
12. Юрин В. В. Оптимальная коррекция параметров орбиты космического аппарата с двигателем малой тяги //Космические исследования. - 1983. Т. 21, №5. - С. 666 - 674.
13. Васильев В. В., Салмин В. В. Многошаговые алгоритмы коррекции орбиты спутника Земли двигателем малой тяги //Космические исследования. - 1984. Т.22, № 4. -С. 507 - 519.
14. Салмин В. В., Ишков С. А. Оптимальные программы управления в задаче межорбитального перелета с непрерывной тягой //Космические исследования. - 1984. Т. 22, № 5. - С. 702 - 711.
15. Васильев В. В., Салмин В. В. Выбор универсальных параметров двигателя малой тяги, предназначенного для поддержания орбиты спутника Земли//Космические исследования. - 1984. Т.22, № 6. - С. 858 - 866.
16. Салмин В. В. Оптимизация космических перелетов с малой тягой: Проблемы совместного управления траекторным и угловым движением.- М.: Машиностроение, 1987.
17. Ишков С. А., Салмин В. В. Оптимизация траекторий и параметров межорбитальных транспортных аппаратов с двигателями малой тяги //Космические исследования. -1989. Т.27, №1. - С. 42-53.
18. Салмин В. В., Соколов В. О. Приближенный расчет маневров формирования
орбиты спутника Земли с двигателем малой тяги //Космические исследования. - 1991. Т. 29, № 6. - С. 872-888.
19. Ишков С. А. Сближение космических аппаратов с малой тягой на околокруговых орбитах //Космические исследования. -1992. Т.30, № 2. - С. 165 - 179.
20. Ишков С. А., Милокумова О. Л., Салмин В. В. Оптимизация замкнутых межпланетных перелетов Земля-Марс-Земля с малой тягой //Космические исследования. -1995. Т.33, №2, - С. 210 - 218.
21. Ишков С. А. Расчет оптимальных межорбитальных перелетов с малой транс-версальной тягой на эллиптическую орбиту //Космические исследования. - 1997. Т.35, № 2. - С. 178 - 188.
22. Ишков С. А., Романенко В. А. Формирование и коррекция высокоэллиптичес-
кой орбиты спутника Земли с двигателем малой тяги //Космические исследования. -1997. Т.35, № 2. - С. 11 - 20.
23. Салмин В. В., Старинова О. Л. Оптимизация межпланетных перелетов КА с двигателями малой тяги с учетом эллиптичности и некомпланарности орбит планет //Космические исследования. - 2001. Т.39, № 1. - С. 51 - 59.
24. Храмов А. А., Ишков С. А. Расчет маневров коррекции слабоэллиптических и круговых орбит с двигателем малой и конечной тяги // Известия Самарского научного центра РАН. - 2002. Т.4, №1. - С. 144-152.
25. Салмин В. В., Ишков С. А., Стари-нова О. Л. Методы решения вариационных задач механики космического полета с малой тягой. - Самара: Издательство Самарского научного центра РАН, 2006.
APPROXIMATE METHODS OF CALCULATING OPTIMAL FLIGHTS OF SPACE VEHICLES WITH LOW-THRUST ENGINES. PART I
© 2007 V. V. Salmin, V. V. Vasilyev, S. A. Ishkov, V. A. Romanenko, V. O. Sokolov, O. L. Starinova, V. V. Yurin
Samara State Aerospace University
The first part of the paper presents mathematical formulations of the tasks of optimizing flights of space vehicles with low-thrust engines and methods of their solution. The peculiarities of mathematical motion models used for control optimization within the frames of different tasks are discussed.