УДК 629.785
ВЫБОР ЗАКОНОВ УПРАВЛЕНИЯ ТРАЕКТОРНЫМ И УГЛОВЫМ ДВИЖЕНИЕМ КОСМИЧЕСКОГО АППАРАТА С ЯДЕРНОЙ ЭЛЕКТРОРЕАКТИВНОЙ ДВИГАТЕЛЬНОЙ УСТАНОВКОЙ ПРИ НЕКОМПЛАНАРНЫХ МЕЖОРБИТАЛЬНЫХ ПЕРЕЛЁТАХ
© 2013 В.В. Салмин, А.С. Четвериков
Самарский государственный аэрокосмический университет имени академика С.П. Королёва (национальный исследовательский университет)
Поступила в редакцию 25.11.2013
Исследуется задача выбора законов управления траекторным и угловым движением космического аппарата с электрореактивной двигательной установкой и ядерным источником энергии при некомпланарных перелётах. Предлагается методика решения динамической задачи (выбор траекторий и законов управления движением), основанная на итерационной схеме с последовательным усложнением модели движения. На начальной итерации космический аппарат представляется точкой переменной массы. На последующих итерациях используются более сложные модели, учитывающие угловое движение космического аппарата, ограничения на траектории и режимы управления, возмущающие силы и моменты. На заключительной итерации формируется закон терминального управления периодом обращения и долготой стояния на геостационарной орбите. Учитываются основные возмущения, ошибки реализации тяги и ошибки исполнения программы управления. Предлагается алгоритм решения задачи управления, которая формулируется как минимаксная задача с интегро-терминальным функционалом. Приведены результаты математического моделирования процесса перевода космического аппарата в заданную точку геостационарной орбиты с учётом возмущающих факторов. Ключевые слова: межорбитальный транспортный аппарат, электрореактивная двигательная установка, ядерная энергоустановка, геостационарная орбита, угловое движение, закон управления
ФОРМУЛИРОВКА ПРОБЛЕМЫ
Проблема повышения эффективности транспортных операций в космосе в настоящее время приобретает особую актуальность. Одним из возможных путей решения этой проблемы является использование двигательных систем, основанных на использовании электрореактивных двигателей малой тяги. Высокая скорость истечения рабочего тела (15 - 40 км/с) обеспечивает значительно меньший расход рабочего тела по сравнению с двигателем на химическом топливе. Однако перелёты с малой тягой (ускорение от тяги составляет 0,1 - 1,0 мм/с2) в "сильных" гравитационных полях достаточно продолжительны и требуют от нескольких недель до десятков месяцев. Особый интерес для космонавтики представляет геостационарная орбита.
В условиях действия различных возмущающих факторов для геостационарной орбиты необходимо решать задачу "фазирования" (попадания в заданную точку на орбите).
При оптимизации баллистических схем подобных перелётов необходимо искать компромисс между массой полезной нагрузки и продол-
Салмин Вадим Викторович, доктор технических наук, профессор, заведующий кафедрой летательных аппаратов. E-mail: [email protected]
Четвериков Алексей Сергеевич, инженер кафедры летательных аппаратов. E-mail: [email protected]
жительностью перелёта - основными критериями эффективности.
При решении прикладных задач необходимо рассматривать в совокупности проблемы оптимизации траекторий и законов управления движением, а также выбора оптимальных соотношений масс основных компонентов КА с электрореактивным двигателем (ЭРД) и солнечной или ядерной энергетической установкой (ЯЭУ) [1].
На сегодняшний день существуют различные проекты межорбитальных транспортных аппаратов (МТА), причём предполагается многоразовое их использование, т.е. перелёт включает в себя выведение полезного груза на целевую орбиту и возвращение на исходную орбиту. Недостатком МТА с солнечной энергоустановкой является большая потребная площадь солнечных батарей (800 - 1200 м2). Кроме того, КА с солнечной электрореактивной двигательной установкой (ЭРДУ) при движении в околоземном пространстве периодически попадает в тень Земли, где двигательная установка выключается, что приводит ещё большему увеличению продолжительности перелёта.
Альтернативой является создание МТА с ядерной энергоустановкой. Подобные проектные проработки проводятся в РКК "Энергия" (рис. 1), исследовательском центре имени М.В. Келдыша (рис. 2) [2, 3, 4].
Аппараты такого типа имеют значительные массово-инерционные характеристики, что суще-
Рис. 1. Многоразовый межорбитальный буксир: 1 - термоэмиссионная ЯЭУ; 2 - ферменная конструкция; 3 - приборный отсек; 4 - антенна системы связи; 5 - контейнер с ПН; 6 - поворотные блоки ЭРДУ
ственно затрудняет процесс управления движением. Возникает проблема совместной оптимизации траекторного и углового движений. В процессе полёта на аппарат действуют различные возмущающие силы и моменты, которые будут приводить к значительным отклонениям фактической траектории полёта от номинальной.
Таким образом, задача выбора законов управления движением межорбитального транспортного космического аппарата с двигательной установкой малой тяги требует дополнительных исследований.
1. ОБЩАЯ ЗАДАЧА ОПТИМИЗАЦИИ.
ВЫДЕЛЕНИЕ ДИНАМИЧЕСКОЙ ЗАДАЧИ
1.1. Баллистическая схема замкнутого перелёта
Рассматривается баллистическая схема перелёта, которая предполагает выведение аппарата с исходной круговой орбиты на рабочую орбиту (геостационарная орбита) с приведением в заданную точку на орбите и возвращение аппарата обратно на исходную орбиту (рис. 3).
Рис. 2. Многоразовый межорбитальный буксир
с капельным холодильником-излучателем: 1 - ЯЭУ; 2 - модуль ЭРД; 3 - полезная нагрузка;
4 - капельный холодильник излучатель
Задача оптимизации в общей постановке формулируется как задача совместного выбора проектных параметров р , баллистических параметров Ь и совокупности функций ы(Ь, х), х(Ь) из множества допустимых Д обеспечивающих перевод КА из начального состояния х(Хд ) = Хд в конечное многообразие х(^к ) € Xк при максимальном значении критерия оптимальности ¡Л .
В качестве основного критерия оптимальности /и принимается относительная масса полезной нагрузки (отношение массы полезной нагрузки МПН к стартовой массе аппарата М0):
МПН *
¡Л = —; другим критерием выступает общая
М о
продолжительность перелёта Т, складывающаяся из моторного времени ТМ и времени навигационных измерений ТН.
Граничные условия прямого и обратного перелётов будут выглядеть следующим образом:
при прямом перелёте
г(0) = го, V (0) = Уд, М (0) = Мо,
г(Тх) = гк,У(ТХ) = Ук,М(ТХ) = Мд -Мть (1)
при обратном перелёте
. Исходная орбита
Рис. 3. Баллистическая схема межорбитального перелёта с возвращением
г (0) = гк ,У (0) = Ук, М (0) =
= М0 - Мт 1 - Мпн , г (Т2) = г0,У (Т2) = У0, М (Т2) =
= М 0 - Мт 1 - Мпн - МТ 2-
(2)
N = ■
2Т '
где т^т - тяговый КПД.
Масса рабочего тела на прямой и обратный перелёты определяется как
Мт 1 + Мт2 = М0 (+ ^2 - ^2 ) - Мпн г2 ,(6)
где
21 = 1 - ехр| -
У
ХК1
2 2 = 1 - ехр I -
У
ХК 2 С
VХК1, VХК2 - затраты характеристической скорости на прямой и обратный перелёты.
1.3. Выделение динамической задачи
Критерий ^ можно представить в виде [5] / = тах / (р, Ь, и (г), Т, М 0 ) =
и (г )еи ре Р
Здесь Т1 и Т2 - продолжительность прямого и обратного перелётов, г0, гк - радиусы начальной и конечной орбит, У0, Ук - скорости МТА на начальной и конечной орбитах, М - текущая масса аппарата, М0 - стартовая масса МТА, МТ1, МТ2 -соответственно масса рабочего тела, расходуемого на прямой и обратный перелёты, МПН - масса полезной нагрузки.
1.2. Проектная модель массы КА
Стартовую массу МТА с ядерной электрореактивной энергодвигательной установкой представим в виде суммы масс отдельных его частей, которые определяются по функциональным признакам [1]:
М0 = Мэу + МДУ + МТ1 + МТ 2 + МСПХ+МПН + Мк ,(3) где МЭУ - масса ядерной энергоустановки; МДУ -масса двигательной установки; МСПХ - масса системы подачи и хранения рабочего тела; МК - масса конструкции МТА.
Здесь массы отдельных компонентов определяются следующими зависимостями:
МЭУ = (*ЭУN, МДУ = УД •Р, МТ1 = (Р/с)Т1,
МТ2 = (Р/с)Т2 , МСПХ = ГСПХ (МТ1 + МТ2) ,
Мк = / КМ 0, (4)
где (%эу - удельная массовая характеристика энергоустановки; уд - удельная масса двигательной установки; /ик - относительная масса конструкции; Успх - отношение массы системы подачи и хранения рабочего тела к массе рабочего тела; Р - тяга маршевых и управляющих двигателей; с - скорость истечения рабочего тела.
Потребная мощность энергоустановки определяется выражением:
Р • С
, (5)
= тах
и (г )еи ре Р
Г 1 - А •( + 2 2 - 2, • 2 2 )-/к Л 1 - А • 2 2
(7)
(
где а = 1 + Успх +
Т
аэу
Л
2 •ЛТ
уд
Задача оптимизации в такой постановке разделяется на две задачи: динамическую (выбор оптимальных траекторий и законов управления) и параметрическую (выбор оптимальных проектных параметров).
В рамках данной работы проектные параметры определяются только на начальной итерации. Оптимизация баллистических схем, траекторий и режимов управления в качестве конечной цели сводится к получению зависимостей УХК1(и X Р, х0, ХК ) и УХК2(и(гX P, x0, ХК ).
1.4. Модель управляемого движения КА
Введём следующие системы координат: неподвижную экваториальную геоцентрическую систему координат ОХИУИ2И, подвижную орбитальную систему координат ОТ8Ж, оси которой направлены по трансверсали ОТ, по радиус-вектору ОБ и по нормали к плоскости орбиты OW и связанную с корпусом КА ОХУ2 с осями, совпадающими с главными центральными осями инерции КА (рис. 4).
В модель движения МТА кроме уравнений движения центра масс аппарата необходимо включить уравнения движения относительно центра масс, а также ограничения на управление вектором тяги [6]:
Рис. 4. Положение орбиты и схема управления вектором тяги
С
с
+
С
dr - dV _ _ - P _ _ _ ? — = V ; — = a + g + f = -8-e + g + f ■ dt dt m
dm t \ P c c
— = ~\q+qупр )=~ c8+дупр8
f P c
упр
dm dt
dk dt
= I o1 t )(
= шхk.
+ M.
вн
mx
I o ()m);
(8)
Тогда условие реализации движения с заданной программой управления в скалярном виде записывается следующим образом
М расп М X
> |М
X
М расп
> м
Y
М |асп
М,
> Мzl, (9)
располагаемый управляющий
Здесь г - вектор положения центра масс; а, Я, / - векторы реактивного, гравитационного и возмущающих ускорений; Р = с ■ q - тяга маршевого двигателя; q - секундный расход маршевого двигателя, с - скорость истечения рабочего тела; е - единичный вектор направления тяги; к - единичный вектор, направленный вдоль продольной оси КА; 8 - функция включения-выключения и реверса тяги маршевого двигателя; qупр и 8 упр - секундный расход и функция включения-выключения управляющих двигателей; Ш- вектор угловой скорости КА; 10 = 10[т(Щ = 10(Ь)- матрица тензора инерции; М - вектор момента, необходимого для поддержания заданной ориентации КА в орбитальной системе координат; Мт - главный возмущающий момент внешних сил (основной вклад даёт гравитационный момент и момент от тяги ЭРДУ, возникающий в результате смещения центра масс).
1.5. Условие реализации заданной программы управления
Угловая скорость Ш в выражении (6) складывается из угловой скорости Ш>1, которая определяет вращательное движение осей орбитальной системы координат вокруг центра притяжения вследствие движения КА по оскулирующей орбите, и угловой скорости <Х>2 , обусловленной программными разворотами связанной системы координат относительно осей орбитальной системы координат.
Примем, что оси системы координат ОХУ2 и постоянно совпадают, а управление вектором тяги осуществляется путём поворота блока ЭРД при постоянной ориентации корпуса КА (<2 = 0 ). Тогда вращательное движение КА определяется только угловой скоростью <<. Момент, необходимый для поддержания постоянной ориентации корпуса КА, получим из уравнений системы (8) для углового движения
М = 10 - Мвн + Ш1 х 10 (< .
м
где +расп момент.
При этих условиях управление вектором тяги можно рассматривать независимо от углового движения аппарата.
Выражение (9) позволяет оценить максимальный уровень "стабилизирующих" моментов по связанным осям КА и проверить возможность их создания путём доворотов отдельных блоков ЭРД и форсирования-дросселирования тяги.
2. ИТЕРАЦИОННАЯ СХЕМА РЕШЕНИЯ ДИНАМИЧЕСКОЙ ЗАДАЧИ
Задачу оптимизации будем решать итерационным методом, используя модели различной степени точности. На начальном этапе КА с ЭРДУ представляется точкой переменной массы. Предполагается, что КА всё время движется по околокруговой орбите с "медленным" изменением среднего радиуса и наклонения, что позволяет получить "асимптотическую" модель движения с помощью процедуры усреднения (Модель А).
На следующем этапе проводится анализ условий реализации выбранных программных траекторий с точки зрения возможности управления угловым движением КА (Модель Б). Оценивается уровень потребных управляющих моментов и угловых скоростей, дополнительных затрат рабочего тела.
На заключительном этапе уточняются дополнительные затраты энергетики на реализацию перелёта и формируется закон управления в форме обратной связи (Модель В). Здесь учитываются основные возмущения, ошибки реализации тяги и ошибки исполнения программы управления.
2.1. Расчёт номинальной траектории и проектно-баллистических параметров перелёта
Примем, что радиальная составляющая ускорения равна нулю. Тогда направление тяги задаётся углом щ между трансверсалью и вектором тяги (рис. 4, а), тогда проекции реактивного ускорения на оси орбитальной системы координат равны:
PP — cos^, as = 0, aw =— sin^. MM
В работе [7] с помощью метода усреднения получена модель движения для околокруговой орбиты
dr
— = 2 • a dt
1
di
du dt \
М3 dt
Мз dVX
= a,
w .
cos u,
Мз
dt
=a.
(10)
Здесь г - средний радиус орбиты, I - наклонение орбиты, Ух - текущая характеристическая скорость, и - аргумент широты, а - полное реактивное ускорение, / - гравитационный параметр Земли.
При помощи составляющей реактивного ускорения ацг, нормальной к оскулирующей плоскости орбиты, можно изменять наклонение орбиты (рис. 4, б). Из второго уравнения системы (10) следует, что если в точках, аргумент широты
принимает значения
ж 2 :
Зж 2
изменить на-
правление бинормальной составляющей aW на противоположное, то наклонение будет изменяться монотонно.
Таким образом, управление вектором тяги при перелётах между некомпланарными орбитами требует изменения знака бинормальной составляющей реактивного ускорения aW дважды за виток (рис. 5).
Программа управления в рамках данной модели задаётся следующим образом [7]:
VVx , u ) = Vm (Vx ) sign (cos u), (11)
где
Vm (Vx ) = arctg-
2cos^
rK
Рис. 5. Управление вектором тяги на одном витке: - релейная программа; 2 - программа с учётом
продолжительности «переключения» бинормальной составляющей вектора тяги
Рис. 6. Измененение траекторных параметров и параметров управления в процессе полёта (гп = 7171 км, г„ = 42164 км, 1п = 51,6 0, = 0 0)
4 0 К 0 к /
амплитуда периодических колебаний угла ш (рис. 6).
В.Н. Лебедевым в рамках модели усреднённого движения получены аналитические решения, описывающие "универсальную" траекторию перелёта между некомпланарными орбитами, не зависящую от проектных параметров КА. Если конечный радиус гк больше начального радиуса г0 , то радиус г сначала возрастает до максимального значения, а затем убывает (участок торможения), при этом угол у больше 900 (см. рис. 6).
Пример оптимальной фазовой траектории перелёта в координатах г -1 на геостационарную орбиту приведён на рисунке 7.
Если переориентация вектора тяги осуществляется за счёт разворотов корпуса КА, то возникает необходимость совместной оптимизации траекторного и углового движений. Получить решение задачи при такой постановке достаточно сложно. В [6] для усреднённой модели получено оптимальное управление с учётом динамики углового движения. Дополнительные затраты характеристической скорости на управление вращательным движением составляют 4-5% [5]. Однако для реализации управления в этом случае необходимо создавать достаточно большие моменты.
Тогда оптимальная программа управления может быть модифицирована с учётом ограничений на максимальный радиус переходной орбиты, так, чтобы он не превышал радиуса геостационарной орбиты, а также с участком начального разгона (Мупр = 0 ).
Пример расчёта проектно-баллистических параметров для оптимальной программы управления и программы управления, учитывающей ограничения на максимальный радиус переходной орбиты, приведён в табл. 1.
з
r
r
3
r
w - г
w - г
к '0
K
2
2
2
(гк - г)/гк
1 0,8 0,6 0,4
0,2
Í- 1 0
-0,2
Vx = 0.37 км 'с^
V: с= 2,54 км
л V = d -И iri !—2
Vv = 5 ■ i Г _i
Vx=7,6 1 КМ.'с
"Я Т 0 ,7S км. с— 4 0 6 0 8 :
t = О
¡л0
Рис. 7. Фазовые траектории прямого и обратного перелётов: ( гп = 7171 км, г„ = 42164 км, и = 51,6 о, ^ = 0 о):
4 0 K 0 K /
1 - невозмущённая траектория 2 - возмущённая траектория
Таблица 1. Проектно-баллистические параметры перелёта на ГСО с возвращением (М = 30 т)
Ид, КМ И Р, Н с, км/с Ухь км/с Vx., км /с Т, сутки Т2, сутки МТ ] , кг МТ2, кг
с «идеальным» управлением
800 0,302 12,88 25 7,610 7,610 176,9 77,0 7873 3430
1000 0,313 12,60 25 7,516 7,516 178,9 76,5 7790 3332
1500 0,330 12,32 25 7,300 7,300 178,4 74,4 7597 3169
с ограничением на управление
800 0,301 12,88 25 7,627 7,627 177,2 77,2 7887 3436
1000 0,312 12,60 25 7,534 7,534 179,3 76,7 7806 3339
1500 0,329 12,32 25 7,321 7,321 178,9 74,6 7615 3177
2.2. Анализ условий реализации номинальных траекторий
Проводилось численное моделирование тра-екторного и углового движения МТА.
Для моделирования изменения траекторных параметров использовались уравнения в равноденственных элементах с учётом действия возмущающих ускорений от несферичности гравитационного поля Земли, действия гравитационных полей Луны, Солнца и светового давления, а также учитывались исполнительные ошибки в величине и направлении вектора тяги:
dh dt
hi й
"(aT + fr ),
^Т = hi(as +fs )ísinF+OcoF+ey +f )-erVew +fw)], dt й
der = h [- (S + fs ^cosF + ((1 + ü)sinF + er + fr) + ехп( + fw)], dt й
di
X
h
dt di r
(aw + fw )•1 V cos F ,
dt
h í
w
fw )• 2V sm F
(12)
dF Ü 2 ( - )f
— = TI— + \aw + fw )h V
h 3 м
где й=1+eX cosF+er sinF, v =iX sinF-ir cosF
V = 1 + iX + i2, h =
Мз
eX = e cos
(Q + ю),
er = esin(Q + ю), iX = tg2cosQ,ir = tg2sinQ,
F = v + Q + ю , p = Л(1 - e2) - фокальный параметр, A - большая полуось, ю - аргумент перицентра, Q - долгота восходящего узла, n - истинная аномалия.
Для имитации исполнительных ошибок принималось, что двигатели ЭРДУ реализуют малую тягу с ошибкой
P = PH0M (1 + 8пр ),
где P - номинальное значение тяги двигателя, 8пр - предельная относительная ошибка тяги.
Ошибки ориентации вектора тяги будут существенно зависеть от выбора схемы управления МТА с ЭРДУ.
й
Рассматривались схемы управления с поворотами отдельных блоков ЭРД и схема с использованием дополнительных ЭРД для создания требуемых составляющих вектора тяги.
Для схемы управления, при которой ориентация вектора тяги осуществляется за счёт разворотов только двигательных блоков, ошибки ориентации вектора тяги возникают за счёт неточной ориентации блоков ЭРД. Закон распределения ошибок можно принять нормальным:
Э^о»-(1 + ),
¥ = ¥„
•(1 + ),
где $ном , Уном - номинальные значения угла $ (угол между связанной осью ОХ и трансверса-
пр •
лью) и угла ориентации вектора тяги у, 8 §Пр - предельные относительные ошибки соответствующих углов, случайное число, распределенное по нормальному закону в интервале [-1,1].
Для схемы управления с использованием дополнительных ЭРД ошибки ориентации вектора тяги возникают за счёт разброса тяги относительно номинального значения, а также за счёт дискретности тяги в трансверсальном и бинормальном направлениях:
PT
cosy =
■yjl'r + PW
-, siny
PW
W
л1рг + PW
где Pr PW- трансверсальная и бинормальная составляющие вектора тяги:
PT = nT • P0 • (1 + §пр ) ,
PW = nW • P0 •(1 + §ПР ) •
Здесь P0 - тяга одного ЭРД; nT
количе-
ство работающих двигателей в трансверсальном и бинормальном направлениях соответственно, которые определяются следующим образом:
120 100 so 60 40 20 о -20 -40 -60 -80 -100
О
ъ д
¿5 Л* Лс
С о С 1 о
AiK, град
• от возмущений гр.п. Земли, лунно-солнечных возмущений Л от исполнительных ошибок для схемы с поворотом блоков ЭРД О от исполнительных ошибок для схемы с дополнительными ЭРД
nT = ENT
( PH Л
PT
P0
(-п H Л
W
= ENT
PW
W
P0
pw = phom ¡whom-
где PT = Phom COS&hom cosy номинальные значения трансверсальной, бинормальной составляющих вектора тяги.
Отклонения конечных траекторных параметров при действии на МТА возмущающих ускорений и возникновении исполнительных ошибок в реализации вектора тяги представлены в фазовых плоскостях (рис. 8, а) ААК — А/к и ААК — АеК (рис. 8, б).
Здесь ААк , А/'к , Аек - конечные отклонения большой полуоси, наклонения, эксцентриситета соответственно.
Для моделирования углового движения использовались динамические уравнения с учётом действия гравитационного момента и момента от тяги ЭРДУ, возникающий из-за смещения центра масс, и кинематические уравнения, записанные в кватернионной форме.
Примем, что МТА имеет две плоскости симметрии, в этом случае центробежные моменты инерции равны нулю. Тогда динамические уравнения углового движения имеют вид
da
X
M
X
Iz — Iy
dt Ix Ix
daY = M Ix — Iz
dt Iy Iy
daZ = Mz Iy — Ix
dt
I
Z
I
ШY ШZ
ш x rnz,
ш x шy
(13)
Z
где со х, СО у, СО Z - проекции угловой скорости на связанные оси координат ОХ, ОУ, О2; Мх, Му, М2 - проекции моментов всех действующих на
120 100 80 60 40
. 20 5 л
^ -20 -40
-
J
- - — - - - -
Л -1 A Л д с
- О
0 P С
-80 -100
/ье^Ю3
• от возмущений гр.п. Земли, лунно-солнечных возмущений Д от исполнительных ошибок для схемы с поворотом блоков ЭРД О от исполнительных ошибок для схемы с дополнительными ЭРД
а) б)
Рис. 8. Отклонения конечных траекторных параметров
МТА сил на связанные оси координат ОХ, ОУ, О2; 1Х, 1Г, 1г - моменты инерции МТА относительно осей ОХ, ОУ, О2.
Кватернионные кинематические уравнения углового движения, которые в скалярной форме имеют вид
¿4,
0 1 / \ dT = " 2 +ЮгЧ2 +
dq1 1 /
—= -(xqo dt 2
ю
dq2 dt
dq3
=-(■
2
coYqo + юх
qs);
(14)
юxq2
■a>Yqi),
1 (
¿г = 2(
где ц = 40 + 41/ + 42 3 + - кватернион, составленный из параметров Эйлера (параметров Родрига-Гамильтона), задающий ориентацию аппарата в пространстве относительно орбитальной системы координат в текущий момент времени.
Моменты относительно связанных осей от тяги маршевых двигателей, возникающие в результате смещения центра масс аппарата:
м р = 0,
Мр = ±Р Ах
М Р = 0.
Проекции гравитационного момента на свя-
занные оси в кватернионном представлении имеют вид
H
=4 (( - iy )-(q0- q2+q2- q2 - qoqi),
X 3 г
12
ну =— (( - 1х У(42 + 4о4з Мм - Ч2Ч3), г
нг =4(1¥ - 1х-42 + 42 -4))(2 + 4о4з). г
Для оценки управляющих моментов рассчитывались в первом приближении массово-инерционные характеристики КА, который принимается состоящим из набора типовых элементов с равномерным распределением масс. Моменты инерции 1Г и I имеют значения порядка 106 кг X м2, а момент инерции 1Х ~105 кг X м2. Смещение центра масс составляет порядка 0,15 - 0,35 м.
Существенное значение имеет гравитационный момент Н относительно оси О2, его величина составляет порядка 3,5 - 5 Н X м. Наибольший вклад в возмущающий момент относительно оси ОУ вносит момент от тяги ЭРДУ Мур, возникающий в результате смещения центра масс (рис. 9). Его величина составляет порядка 1,5 -5,5 Н X м. Величина возмущающего момента относительно ОХ (основной вклад вносит гравитационный момент) на два порядка ниже и составляет 0,04 - 0,06 Н X м .
На рис. 10 приведены угловые скорости относительно связанных осей аппарата, возникающие в результате действия возмущающих моментов.
N
~ 2
0
¡ 1 МГу /
VÍ Hz ,
20
<10
00
120
НО
160
80 t00 / f сугпки
Рис. 9. Амплитудные значения возмущающих моментов относительно связанных осей МТА:
1 - Нп = 800 км, i = 51,6о, Р = 12,88 Н, А х = 0,3 м;
и 7 и 7 7 7 7 max 7 7
2 - Нп = 1500 км, i = 51,6о, Р = 8,12 Н, Ах = 0,2 м
и 7 и 7 7 7 7 max 7
сутки а)
Рис. 10 - Угловые скорости относительно связанных осей МТА (г0 = 7171 км, Гк = 42164 км, 10 = 51,60, Iк = 00, Р = 12,88 Н, с = 25 км/с, Т1 = 176,9 сут, М0 = 30 т)
2.3. Алгоритм коррекции программы управления на участке дальнего наведения
Когда траектория КА мало отличается от расчётной, можно использовать алгоритмы, основанные на отслеживании номинальной траектории или программы управления.
Коррекция программы управления может производиться за счёт обновления модулирующей функции У т(Ух):
¥щ(Ух) =
К -
'к -
\г_к
2со$^-
К -/
у/ ГК
Здесь Т/, /, V/ - средние значения параметров движения, соответствуют моменту Ь, к которому приведены текущие навигационные измерения, уточняющие вектор состояния.
Если отклонения велики или не удаётся заранее построить точную модель действующих возмущений, существенно изменяющихся в процессе полёта, целесообразно применять многошаговые алгоритмы управления с прогнозом конечного состояния и идентификацией возмущающих факторов. Здесь в качестве возмущающих факторов выступают исполнительные ошибки в величине и направлении вектора тяги.
Номинальная траектория рассчитывается по заданным начальным условиям и принятой априорной модели возмущений.
Поскольку на КА действуют возмущения (возмущающие ускорения, возмущения вектора реактивного ускорения), отличающиеся от принятой априорной модели, математическую модель движения следует периодически приводить в соответствие с измеренными ошибками вектора фазовых координат. Так как в настоящее время разработаны достаточно точные модели возмущений гравитационного поля Земли, лунно-солнечных возмущений, то уточнению подлежит величина реактивного ускорения.
В качестве уточняемых параметров выступают предельные относительные ошибки величины тяги.
Опыт использования ЭРД в космосе показывает, что отклонение величины тяги (систематическая ошибка) относительно номинального значения (формуляра) составляет от - 2,5% до +6%.
После этого строится прогнозируемая траектория. Если прогнозируемый вектор хПр ^к )
не принадлежит области допустимых хк, то следует построить новую программу управления
у = уп (Ь). Настройка алгоритма коррекции производится из условия минимума прогнозируемого конечного промаха, но так, чтобы
(к)е ХК
хпр\1 К)^ ЛК .
Для сужения области отклонения конечных параметров траектории на порядок достаточно 2 - 3-х коррекций (табл. 2). Если каждая кор-
Таблица 2. Результаты коррекции программы управления
2
2
2
-
1-
Количество Отклонение среднего радиуса конечной Отклонение наклонения конечной орбиты
коррекций орбиты Дгк , км Д'к , град
0 -61,45 0,079
1 -4,11 0,066
2 -0,49 0,032
3 0,14 0,030
рекция производится после участка навигационных измерений продолжительностью ~ 1 сутки, то суммарное время перелёта увеличивается на 2 - 3 суток.
3. ФОРМИРОВАНИЕ АЛГОРИТМОВ ТЕРМИНАЛЬНОГО УПРАВЛЕНИЯ ПРИ ПЕРЕЛЁТАХ НА ГСО
3.1 Постановка задачи управления
Данная задача формулируется как минимаксная задача оптимизации с функционалом
* T
I = Т +Ахк ААхк ^ min (11)
с последовательным изменением элементов: большой полуоси, наклонения, эксцентриситета орбиты и долготы, определяющей положение КА на орбите.
T
Здесь Г = (ТМ + ТНАВ)/Тгсп, где Т
М
= jSdt
о
имеет вид А =
Лц Л
12
Л21 Л22
В этом случае терминальная часть функционала (9) принимает вид функции
^ (АЛ, АТ) = Л2АЛ2 + ЛпАТ2 + 2АТАЛ ■ Л,
(12)
Задача решается с помощью многошагового алгоритма с поэтапным уточнением параметров управления. В данном случае под шагом понимается сочетание пассивного и активного участков,
следующих друг за другом. На активном участке действует только трансверсальная тяга. Количество шагов соответствует числу коррекций параметров орбиты.
Последовательность продолжительностей работы двигателей принимается убывающей, и определяться выражением [8]
Tt = а
1 -
i -1 n -1
b
моторное время, ТНАВ - время навигационных измерений.
Эта задача решается в три этапа: на первом этапе КА приводится в экваториальную плоскость; на втором этапе эксцентриситет уменьшается до допустимых значений, на третьем этапе ведётся управление периодом вращения и долготой стояния. Такая схема является приближённо-оптимальной и по существу реализует локально-оптимальные управления.
Для коррекции наклонения на 0,5о - 1,0о в окрестности ГСО потребуется время работы ЭРДУ в пределах 0,8 - 1,6 суток, затраты Ух при этом будут составлять 0,042 - 0,085 км/с. Для коррекции эксцентриситета с 0,005 до 0,0005 потребуется моторное временя в пределах 0,5 суток, затраты Ух составят 0,020 км/с.
Решается задача синтеза закона управления, который минимизировал бы функционал (9), при доведении большой полуоси (периода орбиты Т) и долготы Л до заданных значений (ТгСО, Л Р).
Матрица постоянных коэффициентов Л
12
где I, п - номер и число коррекций соответственно; а, Ь - параметры, характеризующие закон убывания продолжительностей работы двигателей.
Тогда задача определения оптимального закона управления сводится к трёхпараметричес-кой задаче оптимизации, которая формулируется следующим образом: для заданных начальных значений элементов орбиты, трансверсального ускорения аТ, продолжительности пассивных участков найти параметры п, а и Ь, доставляющие минимум функционалу (11).
Минимизация продолжительности манёвра обеспечивается за счёт уменьшения числа коррекций. Параметры а и Ь выбираются из условия минимума функции (12) методом покоординатного спуска.
3.2. Алгоритм решения задачи управления
Общий алгоритм выбор оптимального закона управления будет выглядеть следующим образом:
1) для заданных значений А Т0и А Л 0 определяется минимальное количество коррекций п;
2) для найденного значения п и заданных А Т0 и АЛ 0 определяются начальные значения а0 и Ь0;
3) определяется закон управления т(а,Ь) с одновременным моделированием процесса перевода КА с учётом возмущений;
4) на каждом шаге после имитации навигационных измерений на пассивном участке определяется А Тп и А Л п, для которых повторяются пункты 1-3.
Процесс продолжается до тех пор, пока не будет достигнута заданная точность.
3.3. Определение оптимальных значений параметров управления
Для определения минимального количества коррекций п и начального значения параметра Ьд используется зависимость суммарного смещения долготы АЛ2 за время перелёта от параметра Ь и числа коррекций п (рис. 9). Суммарное смещение по долготе за время всего перелёта можно найти по формуле
А Л = АЛП + АЛа ,
Рис. 9. Зависимость AXS от параметра b и n ( A T0 = 1000 с, tn = 86400 с, aT = 0, 180 мм/с2)
где ДЛП - смещение долготы на пассивных участках, ДАд - смещение долготы на активных участках.
Начальное значение параметра а определяется выражением
а 0 =
A T0
3 • T0
If "TZ
i = 1
1 -
i - 1 n -1
Используя зависимость, представленную на рис. 2, для заданных значений Д А 0, и Д Т0 можно определить минимальное количество коррекций п. Для найденного значения п и заданных Д А 0, Д Т0 можно определить начальное приближение Ь0 и соответственно для параметра а0.
3.4. Определение оптимальных значений параметров управления
Таким образом, остаётся найти оптимальные значения параметров а и Ь, которые обеспечива-
Ш
Рис. 10. Зависимость Fот а при b = const: 1 - b = 0,3; 2 - b = 0,4; 3 - b = 0,5; 4 - b = 0,6 (A T0 =1000 c, X = 75,1о)
ли бы минимум функции F. Двухпараметричес-кую задачу оптимизации будем решать методом покоординатного спуска.
Функция Fв зависимости от параметра а при b = const хорошо аппроксимируется кривой второго порядка d1ai2+d2a+d3 (рис. 10).
Аппроксимирующие коэффициенты d1, d2, d3 определяются методом наименьших квадратов. Тогда оптимальное значение для параметра а определяется по формуле
d
2
aopt 2d1 .
Зависимость функции F от параметра b при а = const имеет более сложный характер (рис. 11). Аппроксимировать такую зависимость с достаточной точностью довольно сложно. Как видно из рис. 11, есть явный минимум функции F. Поэтому при оптимизации по параметру b можно использовать методы одномерной оптимизации. В данном случае удобно использовать метод золотого сечения.
F(b)
Рис. 11. Зависимость Fот b при а = const: 1 - а = 52552 с; 2 - а = 50552 с; 3 - а = 48552 с; 4 - а = 44552 с ( A T0 =1000 c, X = 75,1о)
b
о
3.5. Результаты решения задачи точного формирования геостационарной орбиты с приведением космического аппарата в заданную точку стояния
Пример расчёта программы управления в заданную точку ГСО для двух случаев, с поэтапным уточнением параметров управления и без уточнения, приведён в табл. 3 и 4. На рис. 13 приведены области конечных отклонений периода и долготы для двух случаев. Фазовая траектория
перевода приведена на рис. 12.
Таким образом, программа управления движением КА с ЭРДУ в заданную точку стояния будет состоять из последовательности пассивных и активных участков. Для случаев, представленных в табл. 3 и 4, программа управления будет состоять из 4 пассивных участков продолжительностью 86400 с и 4 активных участков Ti (выделено жирным шрифтом). Характеристическая скорость для Д Т0 = 300 - 1000 с составляет порядка 0,01 км/с.
Таблица. 3. Параметры управления при переводе КА в заданную точку ГСО ( A T = 700 c, a = 0,0009 м/с2, t = 86400 с, А Л 0 = 0,12рад)
№ S, тг, с Т>, с т3, с т4, с
без уточнения параметров управления
1 2657,9 2616,4 2325,7 1536,6
с уточнением параметров управления
1 2657,9 2616,4 2325,7 1536,6
2 2739,0 2359,9 1418,9 -
3 2484,2 1325,2 - -
4 1213,7 81,3 - -
Таблица. 4. Длительности активных участков ( A T0= 700 c, aT = 0,0009 м/с2, t = 86400 с, А Л 0 = 0,12рад)
№ S а b n VXK, км/с aTK, с AXk, рад F
без уточнения параметров управления
1 2657,9 3,000 5 0,008 7,5 -0,003 0,106
с уточнением параметров управления
1 2657,9 3,000 5
2 2739,0 1,800 4 0,008 6 -0,0004 0,020
3 2484,2 1,100 3
4 1213,7 0,1 3
41)1 0 1,01 0,02 0,0! 0.» 0,Ю 0,1» 0,07 0 0! 0.04 0,1 0,11 0,11
М„ рад
Рис. 12. Фазовая траектория перевода в заданную точку ГСО
Рис. 13. Область конечных Д Ти Д Л : 1 - без уточнения параметров управления, 2 - с поэтапным уточнением параметров управления
СПИСОК ЛИТЕРАТУРЫ
1. Механика космического полета (проблемы оптимизации) / Г.Л. Гродзовский, Ю.Н. Иванов, В.В. Токарев. М. Наука, 1975. 704 с.
2. Лопота В.А. Космическая миссия поколений XXI века // Полет. 2010. №7. С. 3-12.
3. Коротеев А.С. Ядерный космос России // Новости космонавтики. 2010. Том 20. №2(325). с. 44-47.
4. Перспективы и эффективность применения космических ядерно-энергетических установок и ядерных электроракетных двигательных установок / В.П. Легостаев, В.А. Лопота, В.В. Синявский // Космичес-
кая техника и технологии. 2013. №1. С. 4 - 15.
5. Салмин В.В., Ишков С.А. Оптимизация траекторий и параметров межорбитальных транспортных аппаратов с двигателями малой тяги // Космические исследования. 1989. Т.27. №1. С. 42-53.
6. Салмин В.В. Оптимизация космических перелетов с малой тягой: Проблемы совместного управления траекторным и угловым движением. М.: Машиностроение, 1987. 208 с.
7. Лебедев В.Н. Расчет движения космического аппарата с малой тягой. М: ВЦ АН СССР, 1968. 108 с.
8. Управление орбитой стационарного спутника / Г.М. Чернявский, В.А. Бартенев, В.А. Малышев. М.: Машиностроение, 1984. 144 с.
SELECTION OF CONTROL LAWS FOR TRAJECTORY AND ANGULAR MOTION OF SPACECRAFT WITH A NUCLEAR ELECTRIC PROPULSION ENGINE DURING NON-COPLANAR INTERORBITAL FLIGHTS
© 2013 V.V. Salmin, A.S. Chetverikov
Samara State Aerospace University named after Academician S. P. Korolyov (National Research University)
The problem of choice of control laws for the trajectory and the angular motion of a spacecraft and the electric propulsion and nuclear energy source for non-coplanar flights. The technique of solving a dynamic problem ( the selection of paths and laws of motion control ), based on an iterative scheme with a consistent pattern of movement complexity . In the initial iteration spacecraft represented by a point of variable mass . In subsequent iterations are used more sophisticated models that take into account the angular motion of the spacecraft , the restrictions on the trajectory and control modes , the perturbing forces and moments . At the final iteration of the terminal control law formed an orbital period of standing and longitude in geostationary orbit. Taking into account the main disturbance implementation errors and execution errors traction control program. An algorithm for solving the control problem , which is formulated as a minimax problem with integro - functional terminal . The results of mathematical modeling of the process of transfer of the spacecraft at a given point of the geostationary orbit , taking into account the disturbing factors. Key words: transfer vehicle, of electro-propulsion, nuclear power plant, geostationary orbit, angular motion, control law
Vadim Salmin, Doctor of Technics, Professor, Head at the Space Vehicles Department. E-mail: [email protected] Alexey Chetverikov, Engineer at the Space Vehicles Department. E-mail: [email protected]