ДИНАМИКА, БАЛЛИСТИКА, УПРАВЛЕНИЕ ДВИЖЕНИЕМ ЛЕТАТЕЛЬНЫХ АППАРАТОВ
DOI: 10.18698/0236-3941-2015-6-4-21 УДК 629.782
ОПТИМАЛЬНОЕ УПРАВЛЕНИЕ КОСМИЧЕСКИМ АППАРАТОМ НА УЧАСТКЕ ПРЕДВАРИТЕЛЬНОГО АЭРОДИНАМИЧЕСКОГО ТОРМОЖЕНИЯ ПРИ ВЫВЕДЕНИИ НА ОРБИТУ ИСКУССТВЕННОГО СПУТНИКА МАРСА
Н.Л. Соколов
ЦУП ЦНИИмаш, г. Королев, Московская область, Российская Федерация e-mail: [email protected]
Исследована задача оптимального управления движением космического аппарата при применении комбинированной схемы его выведения на орбиту искусственного спутника Марса. В качестве основных критериев оптимальности использованы максимум скорости вылета космического аппарата из атмосферы и максимум ширины коридора входа аппарата в атмосферу. Разработаны аналитический метод и алгоритм ускоренного расчета квазиоптимальных траекторий выведения аппаратов на спутниковые орбиты. Полученные данные использованы в качестве первого приближения при решении задач оптимального управления космическим аппаратом в общей постановке. Анализ численных материалов показал принципиальную возможность реализации предложенной комбинированной схемы выведения при использовании космического аппарата с аэродинамическим качеством более 0,3. Полученные результаты имеют практическое значение и могут быть использованы при исследовании конкретных миссий дальнего космоса.
Ключевые слова: космический аппарат, выведение на орбиту, комбинированная схема, оптимальное управление, метод исследований, аналитический алгоритм, эффективность применения.
SPACECRAFT OPTIMAL CONTROL IN DRAG
BRAKING FLIGHT SEGMENT DURING INSERTION INTO ORBIT
OF MARS ARTIFICIAL SATELLITE
N.L. Sokolov
Central Research Institute for Machine Building, Korolev, Moscow Region, Russian Federation e-mail: [email protected]
The article analyzes the spacecraft optimal control problem while using a combined profile of the spacecraft insertion into the orbit of the Mars artificial satellite. The authors use both the maximum atmosphere escape velocity and the maximum width of the re-entry corridor of the spacecraft as the main criteria of optimization. Both an analytical method and algorithm are used for a high-speed calculation of the quasi-optimal trajectories of the spacecraft insertion into the satellite orbits. The obtained results are regarded as an initial approximation for solving the spacecraft optimal control problem in a general statement. The numerical data obtained during the research showed the possibility ofimplementing the proposed combined profile of the spacecraft insertion into the orbit, if a spacecraft with a lift-to-drag ratio of over 0.3 is used. The results are of practical importance and can be used while analyzing certain deep space missions.
Keywords: spacecraft, insertion into orbit, combined insertion profile, optimal control,
research technique, analytical algorithm, implementation efficiency.
Одной из важнейших проблем организации космических миссий для изучения Марса и других планет Солнечной системы является выработка и обоснование требований к энергетическим и массово-габаритным параметрам создаваемых перспективных космических аппаратов (КА). Это обусловливается необходимостью учета различных факторов, влияющих на проектный облик и массовые характеристики автоматических и пилотируемых КА. К ним, в первую очередь, относится разработка мер защиты от радиационных воздействий, средств и способов создания искусственных гравитационного и магнитного полей, технологий производства топлива из местных ресурсов и др. [1,2].
Все это в сочетании с необходимостью решения широкого спектра научно-исследовательских целевых задач предопределяет исключительную важность проблемы поиска путей увеличения доли научной аппаратуры в общем весовом балансе КА. Существенным резервом в решении этой проблемы является организация движения КА по энергетически оптимальным траекториям.
Практика показывает, что при планировании космических экспедиций к планетам Солнечной системы, как правило, предполагается использование ракетодинамических схем формирования спутниковых орбит (например, программа "ExoMars" [3]). Вместе с тем, представляется интересным оценить энергетическую эффективность применения комбинированной схемы, предусматривающей предварительное аэродинамическое торможение КА в атмосфере и последующий его разгон в апоцентре промежуточной орбиты.
Проблеме поиска оптимального управления КА при движении на припланетных участках и в атмосфере посвящен ряд работ отечественных и зарубежных авторов [4-15]. Проводимые исследования основывались на использовании фундаментальных методов оптимального управления КА [4-10]. Кроме того, известен ряд работ, посвященных решению отдельных задач оптимизации траекторий полета КА в атмосфере при его выведении на орбиту искусственного спутника Марса (ИСМ). Однако в них рассматривается, как правило, однопараметри-ческое управление КА углом крена на основе упрощенных математических моделей движения КА в атмосфере [11-15]. В результате, найденные решения носят частный характер, а рассчитанные траектории полета КА содержат большие вычислительные погрешности. В настоящей работе исследуется совместное оптимальное двухпараме-трическое управление углами крена y и атаки а при формировании орбит ИСМ для системы дифференциальных уравнений, описывающих пространственное движение КА.
Следует отметить, что решение задач оптимального управления космическими аппаратами (КА) связано с рядом существенных трудностей. Система дифференциальных уравнений движения КА с учетом влияния действующих на него сил не может быть преобразована к аналитическому виду, так как эти уравнения не имеют конечного решения в результате их интегрирования. Поэтому возникает необходимость в использовании численных методов оптимизации нелинейных систем уравнений, предусматривающих проведение многопараметрических итерационных процессов решения краевых задач [5-9]. Практика показывает, что на точность и быстродействие вычислительных процессов решения таких задач значительное влияние оказывает выбор первого приближения краевых значений, вектора состояния системы и вектора сопряженных переменных. При этом если количественная оценка параметров вектора состояния может быть дана на основе анализа динамики движения КА, то предварительно оценить значения всех сопряженных переменных, не имеющих явного физического обоснования, практически невозможно. Вместе с тем, неудачный выбор первого приближения может привести не только к большой продолжительности решения краевых задач, но и в ряде случаев к несходимости вычислительного процесса в принципе.
Одно из перспективных направлений поиска состоит в разработке аналитических методов определения структуры оптимального управления и расчета параметров движения КА. Кроме сокращения затрат расчетного времени это дает и другие преимущества. Форма решения получается более наглядной, что облегчает проведение сравнительного анализа различных вариантов. Оптимальные законы управления могут иметь вид явных зависимостей от начальных условий и параметров системы, что позволяет оценить степень влияния той или иной характеристики на управляющие параметры.
Постановка задачи. Рассмотрим основные критерии оптимального управления КА, определяющие принципиальную возможность и эффективность реализации предлагаемой схемы:
— минимум потребных энергозатрат на формирование орбит ИСМ;
— максимум коридора входа КА в атмосферу.
Движение КА в атмосфере по аналогии с работами [5, 6, 15] описывается системой дифференциальных уравнений в скоростной системе координат с учетом влияния гравитационных, аэродинамических, центробежных и кориолисовых сил в предположении центральности поля тяготения:
dV_ ~dt
pV 2Cx(a)S 2m
— g sin 9 — ш2т cos ф (sin ^ sin e cos 9 — cos ф sin 9);
d9 pVCy (a)S g V
cos y — — cos 9 +--cos 0 + cos ^ cos
dt 2m V r
w2r
+ cos w (sin w sin e sin 0 + cos w cos ;
de pVCy (a)S sin 7 V 2ш
— =-----cos 0 cos e tg w--- (cos 0 sin w—
dt 2m cos 0 r cos 0
л , w2r . cos e
— sin e sin 0 cos w)--— sin w cos w--; (1)
V cos 0
dh dA V cos 0 cos e dw V
— = V sin 0, — =--, — = —cos 0 sin e;
dt dt r cos w dt r
7 ^ Cy (a) m
r = R + h, g = ^, K = , Px = n ( w;
r2 Cx (a) Cx(a)S
здесь V — скорость КА, 0 — угол наклона вектора скорости к местному горизонту, e — курсовой угол, r — радиус-вектор, соединяющий центр планеты и положение КА, A и w — долгота и широта подспутниковых точек КА, m — масса КА, t — время, р — плотность атмосферы, Cx и Cy — аэродинамические коэффициенты лобового сопротивления и подъемной силы соответственно, R — радиус планеты, h — высота полета, g — ускорение силы тяжести, ^ — произведение постоянной притяжения на массу планеты, S — площадь миделева сечения.
Значения управляющих параметров a и y могут изменяться в пределах
0 < a < amax, —п < y < п. (2)
Для различных моделей атмосферы Марса (минимальная, номинальная, максимальная) плотность р в зависимости от высоты полета КА определяется в соответствии с методикой, изложенной в работах [16, 17].
Значения коэффициентов Cx и Cy зависят от форм КА. В качестве примеров рассматривались аппараты сегментно-конической формы с максимальным аэродинамическим качеством Kmax = 0,34; типа "несущий корпус" с Kmax = 1,5; самолетной формы с Kmax = 2,4. Для таких форм зависимости Cx, Cy и K от угла атаки a приведены на рис. 1 [5, 15].
Будем считать, что начальная точка траектории t = to соответствует моменту входа КА в атмосферу Марса. При этом все значения начальных параметров КА известны:
V (to) = Vo, 0 (to) = 0о, e (to)= eo, h (to) = ho, A (to) = Ao, w (to) = Wo.
Рис. 1. Зависимости коэффициентов лобового сопротивления Сх, подъемной силы Су и аэродинамического качества К от угла атаки а (М > 4):
--космический самолет с Ктах = 2,4;----КА типа "несущий корпус" с
Ктах = 1,5;------КА скользящего типа с Ктах = 0,34
В конечной точке траектории £ = £к (вылет КА из атмосферы) известно значение высоты полета Н (£к) = Натм или радиуса-вектора
КА г (£к):
Г (£к) = Гк = Я + Натм; (4)
здесь Натм = 100 км — высота условной границы атмосферы Марса.
При выведении на орбиту ИСМ с заданным радиусом апоцентра га должно выполняться соотношение, связывающее значения конечных параметров Ук, 9к и гк:
К =
2^Га (Га - Гк) Гк (Г2а - Г2 COS2 вк) '
(5)
В качестве критерия оптимальности использовался максимум скорости КА в конечной точке траекторий 3 = Ук = тах, что обеспечивает минимум потребных энергозатрат при формировании орбиты ИСМ [5, 15].
Задача максимизации коридора входа КА в атмосферу сводится к решению двух независимых вариационных задач о нахождении минимума и максимума высот условного перицентра, характеризующих верхнюю и нижнюю границы коридора входа:
или
О
min h^ = min
7 B
min hB = max
Ч2 V2 cos2 e0H_ R
к(е + 1)
Г2 Y02 cos2 вB
Ke + 1)
R
где е — эксцентриситет подлетной орбиты.
Сформулируем задачу оптимального управления КА в общем виде: для процессов, описываемых системой дифференциальных уравнений (1), требуется определить программу управления углами и 7(£), обеспечивающую экстремум функционала 7 при ограничениях (2) и краевых условиях (3)-(5).
Верхняя и нижняя границы коридора входа определяются путем итерационного поиска значений траекторных углов #0, при которых при оптимальном управлении КА обеспечиваются заданные ограничения и краевые условия.
Указанные задачи решались с помощью принципа максимума Пон-трягина [18].
Запишем гамильтониан:
h = £ № =
i=1
pV2Cx (a) S 2m
pVCy (a) S
+
Ф1 +
cos yФ2 +
pVCy (a) S
sin YФ3 + Ф,
2m ' * ' 2m cos 9
где Ф — функция, не зависящая в явном виде от управляющих параметров а и y •
Сопряженные переменные имеют вид:
dФ4
dH
dH dФ3
dV' dH
dt
dФ5
dH
dt
dФ6
öV dt
dV dt
dH
de '
dH
(6)
Из условия максимума гамильтониана H получим формулы для определения законов оптимального управления углами крена и атаки:
Фэ dCy /da V Ф1 cos 9
Y — arctg_ _^_—_
Ф2 cos 9' дСх/да Ф2 cos 9 cos y + Фэ sin y
С учетом условия трансверсальности определим значения сопряженных переменных и гамильтониана в конечной точке траектории [19]:
Ф1 (¿к) — 1, Фэ (¿к) — Ф5 (¿к) — Фб (¿к) — H (¿к) — 0. (7)
Учитывая, что в правых частях дифференциальных уравнений (1) не содержится в явном виде переменная Л, получаем соотношение
^ф5 — _ дн — 0 d¿ дЛ .
Сопоставляя это соотношение с условием равенства нулю сопряженной переменной Ф5 (¿к), приходим к выводу о том, что Ф5 (¿) = 0 на всем участке полета КА, включая граничные точки траектории
Ф5 (¿о) — Ф5 (¿к) — 0. (8)
Далее, из условия H (¿к) — 0 составим дополнительное уравнение, связывающее неизвестные параметры в конечной точке траектории,
+ + _о (9)
d¿ d¿ d¿
В результате задача определения оптимального управления, обеспечивающего максимум скорости при вылете КА из атмосферы, сводится к решению пятипараметрической краевой задачи для дифференциальных уравнений (1), (6) и краевых условий (3)-(5), (7)-(9).
Как было отмечено ранее, решение такого типа задач классическими методами сопряжено со значительными трудностями. Для упрощения поиска оптимальной структуры управления КА разработан следующий аналитический метод.
Метод расчета оптимальных траекторий. При разработке аналитического метода использовались общеизвестные допущения, обоснованные в ряде работ [5, 7, 12, 14, 15]:
h < R, р — ро exp (-eh), F + F„ < ^гр < Fa,
где Рк, Ргр, Fa — кориолисова, центробежная, гравитационная и аэродинамическая силы соответственно; р0 — плотность атмосферы на поверхности планеты; в — логарифмический коэффициент изменения плотности атмосферы от высоты.
В результате система (1) перепишется в виде
dv — _pv2gs d9 — pvCys cos Y - pVM;
d¿ 2m ' d¿ 2m
de pVCy S sin Y ydh у ■ q.
d¿ 2m cos 9 ' d¿ '
dЛ V cos 9 cos e dw V л . (10)
— — —-, —- — — cos 9 sin e;
d¿ R cos w d¿ R
m1 — Г^ _ Л —, M2 — cos 9 cos e tg w
V2 / рГ 2 рД
Следуя [5, 15], будем считать М! и М2 — кусочно-постоянными функциями.
Введем замены переменных
2т
^ = -¿V , г = - 1п V).
рV 2СЖ£
Это позволит без введения дополнительных допущений упростить анализ уравнений сопряженных переменных и законов оптимального управления.
В результате получим систему, не содержащую в явном виде аргумент г:
d0 _ Cy
dz Cx Y
2mMi de Cy sin y 2mM2
CxS ' dz Cx cos 0 CxS '
dp 2тв sin 0 dA 2m cos 0 cos e d^ 2m cos 0 sin e
dz
(11)
CxS dz pRCxS cos <p dz pRCxS Отметим, что при движении КА в атмосфере аргумент z возрастает.
Для определения оптимальных законов управления параметрами а и y воспользуемся принципом максимума Понтрягина [18]. При zo > z > zK гамильтониан и система уравнений сопряженных переменных запишутся следующим образом:
__ _ Cy _ 2mMi H = Фс + -f cos yФ1
Cx CxS
Ф1 +
Cy sin y Cx cos 0
Ф2 -
2mM2 CxS
1ф2-
2тв sin в 2m cos 0 cos £ 2m cos 0 sin eT
-Фа + ^ „-+-^тт^ФБ; (12)
CxS
pRCxS cos <p
pRCx S
d^i
dz
dH
"d0
Cy sin y sin 0 T 2me cos 0 T
У -Ф2 + —„ n Фа+
Cx cos2 0
+
CxS
2m sin 0 cos e pRCxS cos <p
2m sin 0 sin e T
Ф4 +-^г^—Ф5;
pRCxS
d^
dz
^фэ
dz
dH 2m cos 0 sin £ 2m cos 0 cos e
-Ф4--^TT,—Ф5;
de pRCxS cos <p
pRCxS
(13)
dH 2m cos 0 cos e T 2m cos 0 sin e T -Ф4 +--n Ф5
dp p2RCxS cos <p
dФ4
dz
dH
dA
= 0,
d^
dz
dH
d^
p2RCxS
2m cos 0 cos e sin <p pRCxS cos2 <p
Ф4.
При использовании в качестве аргумента параметра z, согласно [19], в систему (11) вводится дополнительное дифференциальное уравнение dz/dz = 1. В связи с тем, что правые части этой системы не содержат в явном виде аргумент z, соответствующее уравнение для сопряженной переменной Ф0 определяется формулой dФ0/dz = 0.
Согласно сделанному предположению, в уравнения (11)-(13) входят кусочно-постоянные разрывные функции MbM2. Однако в силу теоремы Вейерштрасса - Эрдмана [19] наличие разрывов в правых частях уравнений не нарушает непрерывности гамильтониана и сопряженных переменных:
Ф* [zj + O(z)] = Ф* [zj - O(z)], H [zj + O(z)] = H [zj - O(z)],
где Zj — значения аргумента, соответствующее j-му моменту разрыва функций Mi(z) или M1 (z), O(z) — величина меньшего порядка, чем z.
Законы изменения а и 7 при оптимальном управлении определяются в результате решения системы дН/да = 0, дН/д7 = 0, и их можно записать в виде
дСУ т дСУ sin 7т дСЖт Ф2
—^ cos 7Ф1 + —^ —7Ф2 + Фс = 0, tg 7 = 4. (14)
да да cos в да Ф1 cos в
Граничные условия для сопряженных переменных Ф^ (i = 0,1,... ..., 5) при z = zc и z = zK получим из условия трансверсальности [19]:
I - Hiz + Фс^ + Ф1^в + Ф2$е + Ф3£р + Ф4^Л + Ф5^ = 0. (15)
Таким образом, для определения оптимальных законов изменения управляющих параметров а и 7 необходимо решить уравнения (14) с учетом дифференциальных связей (11)—(13) и краевых условий (15).
В рамках предложенного метода применительно к широкому классу задач оптимального управления КА в атмосфере можно записать общие формулы для определения сопряженных переменных Фс и Ф4:
Фс = ас, Ф4 = а4. (16)
В связи с тем, что гамильтониан H в явном виде не зависит от аргумента z, справедливо соотношение Н = а, что позволяет записать дополнительное уравнение связи между неизвестными параметрами движения КА и сопряженными переменными:
de т de т dp т dw т dA
Ф1 + —Ф2 + -^Фз + -7^5 = а - ас - — а4. (17) dz dz dz dz dz
Другие неизвестные параметры, в том числе сопряженные переменные Ф1 и Ф2, в явном виде влияющие на законы оптимального управления КА, определяются в зависимости от условий поставленных вариационных задач.
Применительно к рассматриваемой задаче максимизации скорости КА при вылете из атмосферы, что соответствует критерию оптимальности J = zк = min, из условия трансверсальности (15) для конечной точки траектории получим:
Ф0к = -1, Ф2к = Ф4к = Ф5к = Нк = 0. (18)
Учитывая соотношения (16)-(18), из уравнений (13) найдем решения для сопряженных переменных:
Фс = -1, Ф2 = Ф4 = Ф5 = Н = 0, Ф3 = а3,
z
Ф1 = Ф1С + cos edz.
z 0
Оптимальный закон управления углом y, обеспечивающий максимум гамильтониана H, имеет вид (см. (12)):
cos y = sign Фь (19)
т.е. либо y = 0 при Ф1 > 0, либо y = п при Ф1 < 0.
Постоянную а3, знаком которой определяется знак производной сопряженной переменной Фь найдем из условия равенства нулю гамильтониана в конечной точке
- (— _ ф — \
a Up 1 ФУ z=zK'
Анализ динамики движения КА на участке его вылета из атмосферы показывает, что приращение угла 9 близко к нулю, тогда как интенсивность изменения плотности атмосферы с увеличением высоты достигает значительных величин. Это позволяет пренебречь вторым слагаемым последнего уравнения.
Поскольку при z = гк dz > 0, dp < 0, то а3 < 0 и dФ1/dz < 0.
Следовательно, функция Ф1(^) является монотонно убывающей и может менять знак с плюса на минус не более одного раза. Итак, в общем случае структура оптимального управления углом крена y представляет собой одноразовое переключение y с 0 на п. При этом принципиально возможны случаи, когда оптимальным является движение КА с постоянными значениями угла крена y: либо с y = 0, либо с y = п.
Для нахождения оптимального закона управления параметром а воспользуемся условием (14). Учитывая, что Ф2 = 0, Ф0 = — 1, получаем:
дСж /дСу
я = cos Y Фь
да да
В частности, для зависимостей
Cx (а) = Cx0 + A sin2 (та — n),
Cy (а) = Cy0 + В sin(ma — n) cos (та — n),
приведенных в работе [5], закон изменения а при оптимальном управлении принимает вид
1 /В cos y ФА n
а = 2m arct4—at") + m (20)
Анализ данного уравнения показывает, что зависимость угла атаки а от аргумента z имеет ярко выраженный минимум атщ = n/m, достигаемый в момент переключения угла крена y. В случаях, если Yopt = 0 или Yopt = п, минимальное значение угла а достигается в конечной точке траектории.
Таким образом, с помощью уравнений (19), (20) определяется структура двухпараметрического оптимального управления углами крена и атаки при движении КА на участке предварительного аэродинамического торможения при выведении на орбиту ИСМ.
Алгоритм ускоренного расчета траекторий движения КА. Для определения траекторий движения КА, соответствующих найденной структуре оптимального управления, разработан алгоритм ускоренного расчета, базирующийся на прогнозировании оставшихся участков полета на основе данных о текущих координатах КА. В качестве независимого аргумента при расчете траекторий использовалась переменная z», изменяющаяся в диапазоне от zmin = zc до zmax = zк с интервалом Az. Значения скорости КА V, траекторного угла в», высоты полета h и курсового угла e вычисляются в зависимости от z».
Скорость Vi определяется в соответствии с введенной заменой переменных для преобразования системы уравнений (10) по формуле
V = e-Zi. (21)
Зависимость траекторного угла в от аргумента z определяется на основе интегрирования первого уравнения системы (11) и может быть записана в виде
в»+1 = вг + (K cos 7 - 2PXM1) (zi+1 - z»). (22)
Из анализа функции M1 следует, что она имеет монотонно возрастающий характер. Сразу после входа КА в атмосферу функция M1 является отрицательной, затем обращается в ноль и принимает положительные значения. Из зависимости (22) следует, что траекторный угол в в процессе движения КА имеет аналогичный характер: меняется от отрицательных значений до положительных.
Поделив первое уравнение системы (11) на третье, получим дифференциальное уравнение с разделяющимися переменными
¿в _ Cy S cos 7 - 2mM1 dp
После его интегрирования с учетом экспоненциальной зависимости плотности атмосферы от высоты получаем соотношение между текущими значениями высоты полета h и траекторного угла в
hi+i = hi - 1 ln
ß cos 0i+1 — cos 6i
(23)
рс M1 - Kcos7/2Px. '
Для определения значений курсового угла e проинтегрируем дифференциальное уравнение, полученное путем деления второго уравнения системы (11) на первое. В результате запишем формулу
Cy sin y
— 2mM2
P = P + cos ^_ (24)
£i+i = £i + ^-ö—vy. (24)
CxS cos y — 2mM1
Таким образом, с использованием соотношений (21)-(24) рассчитываются траектории движения до выхода КА из атмосферы (Л,к = 100 км) с параметрами V^, 9к и ек. При этом расчеты могут проводиться при различных программах управления углами крена и атаки.
Для определения скорости, траекторного и курсового углов при выходе КА из атмосферы в инерциальной системе координат используются следующие соотношения [4]:
Ки = \J V2 + Кр + 2ККр cos 9к cos £к, Кр = ШГк cos ^ки,
д . ( . д Кк \ . ( . Кк cos 9к
9ки = arcsin sin 9к—- , еки = arcsin sin ек
W V VtM cos 0
Поскольку при использовании структуры оптимального управления (19), (20) плоскость движения КА не меняется, широту подспутниковой точки при выходе КА из атмосферы можно рассчитать по формуле в зависимости от наклонения орбиты i и курсового угла еки:
cos i
^ки = arccos-.
cos &ки
Определив параметры V^, 0ки, &ки и вычислив значения кеплеров-ских интегралов энергии Ci и площадей C2
Ci = ^ - К2И, C2 = гк2К2и cos2 0ки, Гк = R + 100 км,
Гк
получим формулы для расчета высоты и скорости КА в апогее переходной орбиты
К = " ^- ^ - Л, К, = ./к2 - 2"(Га- Гк). (25)
у ГаГк
Разработанные зависимости (21)-(25) положены в основу алгоритма определения управляющих параметров, при которых обеспечивается максимум скорости вылета КА из атмосферы и, соответственно, скорости аппарата в апогее переходной орбиты.
Как было отмечено, вход КА в атмосферу при использовании оптимальной структуры управления осуществляется с нулевым углом крена и углом атаки а*, соответствующим максимальному значению аэродинамического качества. Далее с интервалом Дг = — г по соотношениям (21)-(24) пересчитываются текущие значения скорости V, углов £ и высоты полета В результате определяются значения скорости Кк, траекторного #к и курсового £к углов при вылете КА из атмосферы для двух различных режимов полета: с 7 = 0, а = а* и 7 = п, а = а*. В соответствии с уравнением (25) вычисляются значения высот апогея переходных орбит Л.а1(7 = 0) и Л,а2 (7 = п).
Нетрудно видеть, что сразу после входа КА в атмосферу планеты будут справедливы соотношения: ha1 > h^ и ha2 < h^, где h^ — заданная высота апогея формируемой орбиты. В процессе полета КА с нулевым значением угла крена и более позднего его переключения на 7 = п высота ha2 будет возрастать, достигая в определенный момент времени t* заданную высоту Л.азад. Начиная с этого момента в соответствии с установленной структурой оптимального управления (19), (20) угол атаки а уменьшается, обеспечивая тем самым снижение действующей на КА подъемной силы и его полет в более плотных слоях атмосферы. Такое управление углом атаки а сначала обеспечивает снижение интенсивности роста высоты ha2, а затем по достижении некоторого угла а' прекращение роста ha2. После этого происходит разворот КА по углу атаки, что соответствует возрастанию подъемной силы и аэродинамического качества. В результате высота ha2 начинает уменьшаться, приближаясь к заданной высоте Л.азад. В момент достижения ha2 заданной высоты апогея Л,азад угол атаки а устанавливается равным а*, а угол крена 7 принимает значение, равное п. При таком режиме полета обеспечиваются необходимые условия вылета КА из атмосферы и достижение заданной высоты апогея ha2.
Таким образом, с применением описанного вычислительного алгоритма рассчитываются траектории движения аппарата на этапе предварительного аэродинамического торможения при обеспечении максимальной скорости выхода КА из атмосферы. Полученные результаты в дальнейшем использовались в качестве первого приближения при решении задач оптимального управления КА при его выведении на орбиту искусственного спутника Марса, что позволяет существенно сократить продолжительность вычислительного процесса решения вариационных задач.
Анализ численных результатов. Решение краевых задач, выполненное в широком диапазоне условий входа КА в атмосферу, высот формируемых орбит и проектных характеристик аппаратов, позволило определить оптимальные траектории движения КА при двухпараме-трическом оптимальном управлении углами крена и атаки. Показано качественное совпадение численного и аналитического решений. Так, в процессе движения КА в атмосфере угол крена 7 изменяется от 7с « 3 ... 8° до 7 ~ 172 ... 177° (при аналитическом решении угол 7 меняется от 0 до 180°). Угол атаки а при входе КА в атмосферу принимает значение а*, соответствующее максимальной величине аэродинамического качества Kmax. Далее происходит уменьшение а до значения, равного ~10... 12°, а затем а вновь увеличивается до а* (при аналитическом решении также установлен ярко выраженный минимум угла атаки в процессе полета КА).
Рис. 2. Зависимость скорости V, высоты Н, углов атаки а и крена 7 от времени £ при движении в атмосфере Марса (Ктах = 0,34, Нп = 10 км, Рх = т/СхВ = 300 кг/м2, На = 500 км)
В качестве примера на рис. 2 представлены зависимости, показывающие изменение углов атаки и крена, скорости и высоты полета от времени движения КА в атмосфере Марса.
Интенсивность изменения углов 7 и а существенно зависит от высоты условного перицентра траектории входа КА в атмосферу НП. Так, при увеличении НП возрастание угла 7 до значений 172... 177° осуществляется на более раннем участке полета КА в атмосфере, а интенсивность изменения угла атаки снижается. При входе КА в атмосферу по верхней границе коридора НП схема управления углами 7 и а вырождается в движение КА с постоянными значениями этих углов 7 « 175° и а = а*. При входе КА в атмосферу по нижней границе коридора НП угол 7 постоянен и близок к нулю, а изменения угла а от а* до атп и снова до а* осуществляются с максимальной интенсивностью. Отметим, что НП определяется из условия захвата КА атмосферой, а высота НП — некоторыми физическими ограничениями, например максимальной величиной перегрузки или температуры. Другим условием является обеспечение выхода КА на орбиту с заданной высотой апоцентра На.
Анализ результатов показал, что при формировании круговой орбиты ИСМ высотой Н = 500 км для КА с аэродинамическим качеством К = 0,34 и приведенной нагрузкой на лобовую поверхность Рх = 300 кг/м2, учитывая возможный разброс параметров атмосферы, максимальный коридор входа составляет ±25 км. При этом верхняя граница НП тах = 30 км соответствует наименее плотной (минимальной) модели атмосферы, а нижняя граница НП т|П = -20 км — наиболее плотной (максимальной) модели. Полученные значения коридора входа превосходят значение навигационного коридора (при использовании автономных систем навигации КА дНПав = ±10 ... 20 км [9, 15]), что позволяет сделать вывод о принципиальной возможности реализации предлагаемой схемы выведения.
Показано, что потребные энергозатраты д V на формирование заданной орбиты существенно зависят от высоты условного перицентра. Уменьшение высот НП от сначала приводит к незначительному увеличению энергозатрат д V, а затем с приближением кП к нижней границе коридора происходит интенсивное возрастание д V (рис. 3).
Так, при изменении высоты НП от 30 до и -10км энергозатраты д V практически не меняются и составляют д V^ 145 м/с, дальнейшее снижение НП до -20 км приводит к росту д V до 280 м/с (Н = 500 км). Зависимости д V ) имеют аналогичный характер и для случаев формирования более высоких орбит. При этом значение дV увеличивается с ростом высоты Н. Так, при Н = 2000 км энергозатраты дV в зависимости от значений кП составляют 310... 480 м/с, а для Н = 5000 км — дV^ 570 ... 610 м/с. Представленные данные показывают, что энергетически оптимальным является осуществление входа КА в атмосферу вблизи верхней границы коридора.
При использовании форм КА типа "несущий корпус" и самолетного типа, располагающих значительно большими значениями аэродинамического качества Ктах, максимально возможный коридор входа в атмосферу д^п существенно расширяется, в основном за счет уменьшения нижней его границы. Потребные энергозатраты при движении КА по верхней границе коридора д V(^П) практически не меняются, а при полете по нижней границе д V (^П) значительно увеличиваются с ростом Ктах. На рис. 4 приведены зависимости ширины коридора входа д^п и энергозатрат дV(^П) и дV(^П) от максимального значения аэродинамического качества КА.
Видно, что при Ктах = 1,5 (несущий корпус) значения д^п составляют ^200 км, энергозатраты д V(^П) ~ 140 м/с, д V(^П) ~ 340 м/с.
Рис. 3. Зависимости потребных энергетических затрат АУ от высоты условного перицентра траектории входа КА в атмосферу Нп при выведении на круговые орбиты ИСМ с высотами Н (Ктах = 0,34, Рх = ш/Ох Б = 300 кг/м2)
Рис. 4. Зависимости ширины коридора входа КА в атмосферу АКП и энергозатрат АУ(Ь^П) и АУ(ЬП) от максимального значения аэродинамического качества (Н = 500 км, Рх = ш/Ох Б = 300 кг/м2)
При Kmax = 2,4 (самолетная форма КА) Ahn « 320 км, AV(h^) « « 138 м/с, AV(hn) « 430 м/с. Менее существенное влияние на ширину коридора входа и энергозатраты оказывает нагрузка на лобовую поверхность аппарата Px. Изменения Px от 300 до 500 кг/м2 приводит к расширению коридора на -15 км, к росту энергозатрат AV(h^) на -1... 2 м/с, а AV(hn) на -3 ... 5 м/с (Kmax = 0,34, H = 500 км).
В целом для широкого диапазона исходных данных и проектно-баллистических характеристик потребные энергозатраты AV не превышают 650 м/с, а при входе КА в атмосферу вблизи верхних границ коридора hn и формировании орбит ИСМ с высотами H не более 500 км энергозатраты составляют -140 ... 150 м/с.
Для сравнения, в случае применения традиционной ракетодинами-ческой схемы формирования спутниковых орбит, описанной в работах [9, 15], потребные энергозатраты достигают 2,5.. .4 км/с, что в 6-10 раз больше, чем при использовании рассмотренной комбинированной схемы управления.
Заключение. Рассмотрена комбинированная схема выведения КА на орбиту ИСМ, предусматривающая предварительное аэродинамическое торможение КА в атмосфере, перевод на переходную эллиптическую орбиту и подачу разгонного импульса в ее апоцентре.
Разработан метод решения задач оптимального управления КА на основе введения ряда допущений и преобразования систем дифференциальных уравнений движения и сопряженных переменных.
Решены вариационные задачи минимизации потребных энергозатрат при формировании орбит ИСМ и максимизации коридора входа КА в атмосферу при использовании двухпараметрического управления углами крена и атаки.
Показано, что для КА, располагающих аэродинамическим качеством Kmax > 0,34 физически реализуемый коридор входа КА в атмосферу превосходит навигационный коридор, что обеспечивает принципиальную возможность осуществления предлагаемой схемы управления.
Для случаев входа КА в атмосферу вблизи верхней границы коридора потребные энергетические затраты на формирование орбиты ИСМ на порядок меньше, чем при реализации ракетодинамической схемы перевода КА с подлетной гиперболической траектории на заданную орбиту.
Полученные результаты имеют практическую значимость и могут быть использованы при исследовании конкретных миссий дальнего космоса.
ЛИТЕРАТУРА
1. Горшков Л.А. Полет человека на Марс // Наука и жизнь. 2007. №7. С. 4-12.
2. Jonathan Amos. Europe's Mars plans move forward. BBC News. 12.10.2009.
3. Peter B. de Selding. ESA Halts Work on Exo Mars Orbiter and Rover. Space News. 20.04.2011.
4. Челноков Ю.Н. Применение кватернионов в задачах оптимального управления движением центра масс космического аппарата в ньютоновском гравитационном поле // Космические исследования. 2003. № 4. Т. 32. С. 83-91.
5. Иванов Н.М., Мартынов А.И. Движение космических летательных аппаратов в атмосферах планет. М.: Наука, 1985. 384 с.
6. Синявская Ю.А., Корнилов В.А. Иерархическая оптимизация в задачах проектирования систем автоматического управления // Труды МАИ. 2011. № 44.
7. Сапрыкин О.А., Соболевский В.Г. Баллистический анализ вариантов посадки космических аппаратов в заданном районе // Космонавтика и ракетостроение. 2013. № 3 (72). С. 78-86.
8. Okhotsimsky D.E., Golubiev Y.F., Sikharulidze Y.G. Mars orbiter insertion by use of atmospheric deceleration // Acta Astronautica. 1978. Vol. 5. No. 9/10.
9. Матюшин М.М., Соколов Н.Л., Овечко В.М.Совершенствование методологии оптимального управления космическими аппаратами дальнего космоса // Актуальные проблемы российской космонавтики. Материалы XXXIX академических чтений по космонавтике. Москва, январь 2015. С. 297-298.
10. Иванов В.М., Лобачев В.И., Соколов Н.Л.Управление космическим аппаратом баллистического типа при спуске с орбиты в заданную область поверхности Земли // Фундаментальные исследования. 2014. № 8. Ч. 3. С. 577-582.
11. Hiltz A.A., Florense D.E., Low D.L. Selection, development and characterization of a thermal protection system for a Mars entry vehicle // AIAA Paper. 1968. No. 304.
12. Ярошевский В.А. Приближенный расчет траектории входа в атмосферу. Ч. I, II // Космические исследования. Изд. АН СССР. 1964. Т. 2. Вып. 4, 5. С. 15-21.
13. Griffin J.W., Vinh N.X. Three-dimensional optimal maneuvers of hyper velocity vehicles//AIAA. 1971.
14. Чепмен Д.Р. Приближенный аналитический метод исследования входа тел в атмосферы планет. М.: ИЛ, 1962. 298 с.
15. Иванов Н.М., Мартынов А.И.Управление движением космического аппарата в атмосфере Марса. М.: Наука, 1977. 134 с.
16. Данченко О.М. Математическая модель плотности атмосферы Марса // Труды МАИ. 2012. № 50. С. 10.
17. Мороз В.И. Рабочая модель атмосферы и поверхности Марса // Препринт ИКИ. 1975. Вып. 240. 241 с.
18. Понтрягин Л.С., Болтянский В.Г., Гамкрелидзе Р.В., Мищенко Е.Ф. Математическая теория оптимальных процессов. М.: Наука, 1969. 384 с.
19. Летов А.М. Динамика полета и управления. М.: Наука, 1969. 360 с.
REFERENCES
[1] Gorshkov L.A. Human flight to Mars. Nauka i zhizn' [Science and Life], 2007, no. 7, pp. 4-12 (in Russ.).
[2] Jonathan Amos. Europe's Mars plans move forward. BBC News, 12.10.2009.
[3] Peter B. de Selding. ESA Halts Work on Exo Mars Orbiter and Rover. Space News, 20.04.2011.
[4] Chelnokov Y.N. The use of quaternions in the optimal control problems of motion of the centre of mass of a spacecraft in the Newtonian gravitational field. Kosmicheskie issledovaniya [Space research], 2003, no. 4, vol. 32, pp. 83-91 (in Russ.).
[5] Ivanov N.M., Martynov A.I. Dvizhenie kosmicheskikh letatel'nykh apparatov v atmosferakh planet [Movements of spacecraft in atmospheres of planets]. Moscow, Nauka Publ., 1985. 384 p.
[6] Sinyavskaya Y.A., Kornilov V.A. Hierarchical optimization in the problems of automated control system design. Tr. MAI [Proc. MAI], 2011, no. 44 (in Russ.).
[7] Saprykin O.A., Sobolevsky V.G. Ballistic analysis of the variants of spacecraft landing in the specified area. Kosmonavtika i raketostroenie [Space navigation and rocket engineering], 2013, no. 3 (72), pp. 78-86 (in Russ.).
[8] Okhotsimsky D.E., Golubiev Y.F., Sikharulidze Y.G. Mars orbiter insertion by use of atmospheric deceleration. Acta Astronnautica, 1978, vol. 5, no. 9/10.
[9] Matyushin M.M., Sokolov N.L., Ovechko V.M. The development of optimal control methodology of deep space spacecraft. Aktual'nye problemy rossiyskoy kosmonavtiki. Materialy XXXIX akademicheskikh chteniy po kosmonavtike [Actual problems of Russian cosmonautics. Proc. of XXXIX academic readings on cosmonautics]. Moscow, January 2015, pp. 297-298 (in Russ.).
[10] Ivanov V.M., Lobachev V.I., Sokolov N.L. Control of a ballistic type spacecraft during the descent in the specified area of the Earth's surface. Fundamental'nye issledovaniya [Basic researches], 2014, no. 8, part 3, pp. 577-582 (in Russ.).
[11] Hiltz A.A., Florense D.E., Low D.L. Selection, development and characterization of a thermal protection system for a Mars entry vehicle. AIAA Paper, 1968, no. 304.
[12] Yaroshevsky V.A. Approximate calculation of the atmospheric reentry trajectory. Part I, II. Kosmicheskie issledovaniya. Izd. ANSSSR [Cosmic research], 1964, vol. 2, iss. 4, 5, pp. 15-21 (in Russ.).
[13] Griffin J.W., Vinh N.X. Three-dimensional optimal maneuvers of hyper velocity vehicles. AIAA, 1971.
[14] Chapman D.R. An approximate analytical method for studying entry into planetary atmospheres. NASA-TR-R-11[R]. Washington.
[15] Ivanov N.M., Martynov A.I. Upravlenie dvizheniem kosmicheskogo apparata v atmosfere Marsa [Spacecraft motion control in the atmosphere of Mars]. Moscow, Nauka Publ., 1977. 134 p.
[16] Danchenko O.M. Mathematical model of the Martian atmospheric density. Tr. MAI [Proc. MAI], 2012, no. 50, p. 10 (in Russ.).
[17] Moroz V.I. The working model of Mars atmosphere and surface. Preprint of IKI (Space Research Institute), 1975, iss. 240, 241 p. (in Russ.).
[18] Pontryagin L.S., Boltyanskiy V.P., Gamkrelidze R.V., Mishchenko E.F. Matematicheskaya teoriya optimal'nykh protsessov [Mathematical theory of optimal processes]. Moscow, Nauka Publ., 1969. 384 p.
[19] Letov A.M. Dinamika poleta i upravleniya [Flight dynamics and control]. Moscow, Nauka Publ., 1969. 360 p.
Статья поступила в редакцию 23.06.2014 Соколов Николай Леонидович — канд. техн. наук, старший научный сотрудник, заместитель начальника ЦУП ЦНИИмаш.
ФГУП ЦНИИмаш, Российская Федерация, 141070, Московская обл., г. Королeв, ул. Пионерская, д. 4.
Sokolov N.L. — Ph.D. (Eng.), Senior Researcher, Deputy Director of the Mission Control Centre of the Federal State Unitary Enterprise The Central Research Institute for Machine Building.
Federal Unitary State Enterprise Central Research Institute for Machine Building, ul. Pionerskaya 4, Korolev, Moscow Region, 141070 Russian Federation.
Просьба ссылаться на эту статью следующим образом:
Соколов Н.Л. Оптимальное управление космическим аппаратом на участке предварительного аэродинамического торможения при выведении на орбиту искусственного спутника Марса // Вестник МГТУ им. Н.Э. Баумана. Сер. Машиностроение. 2015. № 6. C. 4-21.
Please cite this article in English as:
Sokolov N.L. Spacecraft optimal control in drag braking flight segment during insertion into orbit of Mars artificial satellite. Vestn. Mosk. Gos. Tekh. Univ. im. N.E. Baumana, Mashinostr. [Herald of the Bauman Moscow State Tech. Univ., Mech. Eng.], 2015, no. 6, pp. 4-21.