CC BY
21
УДК 519.6
БОТ: 10.18698/2309-3684-2021-2-5467
Моделирование и оптимизация перелета спутников малой массы с Земной орбиты на орбиту Марса с помощью ионных двигателей
© Т.Ю. Мозжорина, Л.О. Чуванова
МГТУ им. Н.Э. Баумана, Москва, 105005, Россия
В данной работе рассматривается оптимизация перелета спутника малой массы с Земной орбиты на орбиту Марса с использованием ионных двигателей. Ионный двигатель позволяет минимизировать расход топлива и разогнать космический аппарат до довольно высоких скоростей вдали от планет солнечной системы. Рассмотрению подлежит гелиоцентрический участок полета. Ставится задача минимизации времени перелета. В работе приняты следующие допущения: орбиты Земли и Марса являются круговыми и лежащими в одной плоскости. В качестве управления выбирается угол между тангенциальной скоростью космического аппарата в гелиоцентрической системе и направлением действия тяги. При составлении алгоритма оптимизации использован принцип максимума Понтрягина, который приводит задачу оптимизации функционала к краевой задаче для системы обыкновенных дифференциальных уравнений. Решение краевой задачи найдено одним из численных методов — методом пристрелки, дающим наиболее точные результаты. Проведен анализ полученных результатов и проведено сравнение с данными, полученными ранее в подобных расчетах зарубежными авторами другим численным методом решения. Делается вывод о работоспособности метода пристрелки при решении подобных задач.
Ключевые слова: ионные двигатели, метод пристрелки, краевые задачи обыкновенных дифференциальных уравнений, оптимальное управление, принцип максимума Понтрягина, перелет между орбитами Земли и Марса
Введение. Рассматривается перелет между орбитами Земли и Марса с помощью двигателя с малой тягой, а именно, ионного. Этот двигатель относится к электрическим ракетным двигателям, принцип работы которых основан на преобразовании электрической энергии в направленную кинетическую энергию частиц и которые позволяют существенно увеличить мощность двигателя на единицу расхода рабочего тела, а, следовательно, уменьшить запасы топлива и увеличить полезную нагрузку [1]. Принцип работы ионного двигателя описан в [2]. К достоинствам данного вида двигателей относят малый расход топлива и длительное время активной работы, а к недостаткам — малую величину тяги. Ионные двигатели невозможно использовать для старта с планеты, однако, в условиях невесомости при достаточно долгой работе они позволяют разогнать космический аппарат до скоростей, которых невозможно достичь при использовании других двигателей.
Марс является наиболее перспективной кандидатурой в качестве
цели для исследуемых перелетов, примеры этого приведены в [3]. Начиная с 1960 года, было предпринято немало попыток исследовать Марс, и ещё несколько миссий по его исследованию запланировано в наше время [4].
Математическая модель задачи оптимизации перелета с Земной орбиты на орбиту Марса с использованием принципа максимума Понтрягина представляет собой краевую задачу для системы обыкновенных дифференциальных уравнений (состоящих из уравнений движения и кинематических соотношений), которая решалась методом пристрелки. Этот метод заключается в сведении краевой задачи к задаче Коши, для решения которой существует много численных методов решения, позволяющих получать результат с гарантированной точностью [5]. В качестве управления был выбран угол между касательной к окружности с центром в Солнце и направлением действия тяги. Следует отметить, что принцип максимума Понтрягина представляет собой необходимое условие существования экстремума рассматриваемого функционала, которое позволяет среди всех возможных допустимых процессов отобрать те, которые могут претендовать на роль оптимальных [6], [7].
Аналогичная задача была решена в работе [8], где она решалась с помощью метода функций нагружения для рассмотрения связей, наложенных на конечные значения. Метод пристрелки позволяет получить наиболее точное решение при правильно выбранных начальных значениях вектора независимых переменных (параметров пристрелки), что представляет собой нетривиальную задачу. Сложности, связанные с выбором начального приближения, к сожалению, приводят часто к несходимости метода, что ограничивает возможность его применения при решении задач оптимизации в реальных инженерных задачах. Оптимальность полученного решения была проверена сравнением с данными работы [8] и сравнением полной и частной производных функции Понтрягина, которые должны быть равны между собой при выполнении необходимого условия существования экстремума [6].
Постановка задачи математическое моделирование перелета с Земной орбиты на орбиту Марса. Уравнения движения для расчетов и некоторые технические характеристики двигателя были взяты аналогичными тем, что использовались в работе [8]:
- = -" А \ЯТ .^
* Я I Я ) (1 " ботн * *) Л
* Я (1 - ботн • *)
йЯ/Л — и, йф/ * — у/Я,
<
где t — время, с; u — скорость движения ракеты вдоль радиуса, м/с; v — скорость движения ракеты по касательной к окружности с центром в Солнце, м/с; R — расстояние до Солнца, м; р — полярный угол положения ракеты, рад; R — начальное положение радиуса
Земной орбиты, м; в — угол между касательной к окружности с центром в Солнце и направлением действия тяги (управление), рад;
A — F - G—c, м/с2; F — сила гравитационного притяжения в
mo R0
начальный момент времени, Н; m0 — масса ракеты, кг; G — гравитационная постоянная, м3/с2 • кг; —с — масса Солнца, кг; QomH = Qlmo — отношение расхода топлива двигателя к исходной массе, с"1; ^ - T/m0 — отношение тяги двигателя к исходной массе, м/с2. Минимизируемый функционал
t
J — J dt — tx ^ min.
о
Функция Понтрягина в условиях данной задачи принимает следующий вид:
H = ¥и
— - A 1 +, T°mH—ч sine
R A { R ) (1 - Q0mH • t)
J
С _ T ^
UV + --m—- cos в
R (1 - QomH • t)
J
Здесь щи, y/v, y/R , у — сопряженные переменные.
Традиционно при решении задачи быстродействия методами проекции градиентов или условного градиента последнее слагаемое (подынтегральную 1, умноженную на =— 1) не учитывают, так как
она не влияет на нахождение экстремума. Для метода пристрелки ее наличие принципиально, так как будет использована невязка по равенству нулю функции Понтрягина в конце процесса, что следует из принципа максимума при свободном времени на правом конце.
Вычислив производную по управлению и приравняв ее к нулю дИ/дв = 0, получим для оптимального управления tg#* =yM/yv, откуда, используя тригонометрическое соотношение
tg2** — 1, cos в
выводятся следующие выражения:
8Ш в -
¥и
№ +¥V
008 в -
№ +¥1
Знаки у выражений выбирались при условии выполнения следующего неравенства, обеспечивающего наличие максимума функции Понтрягина:
д2 н
- (~Уи 8Ш в-у 008 в) < 0.
дв2 (1 - 0-отн ■ г) Отсюда получается для управления:
sign(8inв) = sign(уu ), sign(c08в) = Sign(Щv ). Сопряженная система имеет вид:
йг V - Я у - у Я ,
йг и уи - - у Я '
йг N у >2 v(2 - 2 Ау Я2 ич V Я3Я2
. йг - 0.
Краевые условия — известные значения для Земли у0 , щ, Я0, % (начальные условия) и для Марса V, щ , Я (конечные условия). Из условий трансверсальности получим, что у (г1) - 0. Так как эта сопряженная функция постоянна по времени, у- 0 У г е [0, гх ] и, соответственно, можно исключить из сопряженной системы последнее уравнение, а в остальных убрать все слагаемые, зависящие от у .
В качестве параметров пристрелки выбираются недостающие начальные условия, а, именно, у, у, у. Дополнительного
параметр пристрелки — время перелета г , начальное приближение
которого выбирается, как известное примерное время перелета с Земли на Марс, равное 210 дням [9]. Этот параметр пристрелки необходим при переходе к системе ДУ по другому аргументу: готн - гД е [0,1], который обеспечит более аккуратный выход из
<
модуля численного интегрирования. Выберем в качестве первого приближения начальный угол в, равный 45° => у/и0 =у/х,0-
Особенности численного алгоритма решения задачи. Полученная задача решается с помощью метода пристрелки, описанного в [10] и [11], схема которого приведена на рис. 1.
Рис. 1. Схема алгоритма численного решения простейшей задачи оптимального управления методом пристрелки
Интегрирование во внутреннем цикле проводилось методом Рунге-Кутта четвертого порядка [12]. Метод Рунге-Кутта обладает значительной точностью и широко используется при численном решении дифференциальных уравнений [13]. Внешний цикл — модифицированный метод Ньютона [14] вычисляет последующие значения параметров пристрелки, линеаризуя систему алгебраических нелинейных уравнений, уравнений невязок. Для оптимизации шага в модифицированном методе Ньютона вычислялась локальная норма [15], улучшающая сходимость в ряде отдельных случаев.
За невязки приняты:
дх- uf ^ 0,
vf -vn S2-0 ^ 0,
V0
¿3 = ^ 0, 3 R
- н ^ 0.
В нашем случае переменные, из которых состоят дифференциальные уравнения, очень разных порядков, что увеличивает вычислительную погрешность. Чтобы избежать этого, все параметры пристрелок и уравнения невязок используются в относительном виде. Программа расчета написана на языке программирования С++ с переменными, описанными двойной степени точности.
Для решения системы линейных алгебраических уравнений используется метод LUP-разложения, алгоритм которого приведен в [16].
Для проверки оптимальности полученного решения сравниваются полная и частная производные функции Понтрягина по времени. Частная производная в данной задаче имеет следующий вид:
ЛН T Q
- - ошн • QornH sin 0 + ¥ cos Q)
^ (1 - О^отн ■ t fWU ^ )
Полная производная может быть численно получена, как:
dH _H1-H1_1 dt a t
Результаты расчетов. Для проведения численных расчетов были использованы следующие исходные данные:
U0 - 0; V0 - 29,8; R - 149,6-106; % - 0; uf - 0; vf - 24,1; R - 227,9 -106.
Начальные значения для неизвестных сопряженных переменных выбирались в качестве параметров пристрелки и были выбраны равными
¥и(t0) = ¥v(t0) « 842; Wr (О - -10-5.
Время перелета с Земли на Марс — 200 суток. Также 0(to) = 45°.
Остальные переменные, участвующие в системе и выбранные в качестве константных, принимают следующие значения:
аотн - 1,29-103 сутки1; Тотн - 8,299-104 ^;
в - 6,6743 -10"
м
с2 - кг
М.
- 1,9891 -1030 кг; А - 0,00593 м
м
Солнца
с
Точность решения, при которой выполнялся выход из метода Ньютона, составляла 106. Для метода Рунге-Кутта было выбрано 1000 шагов.
На рис. 2—5 показаны графики фазовых переменных, полученные в результате решения поставленной задачи.
и, м/с 10000
8000
6000
4000
2000
V, м/с
30000
28000
26000
24000
0 25 50 75 100 125 150 175 (, дни Рис. 2. График скорости ракеты и вдоль радиуса Л
0 25 50 75 100 125 150 175 дни Рис. 3. График касательной скорости ракеты V
с
0
R, -1011 м 2,3 2,2 2,1 2,0 1,9 1,8 1,7 1,6 1,5
р, рад 140 120 100 80 60 40 20 0
0 25 50 75 100 125 150 175 t, дни Рис. 4. График изменения радиуса от Солнца до ракеты R
0
25 50 75 100 125 150 175 /, дни Рис. 5. График изменения полярного угла р
На рис. 6—8 показаны графики для сопряженных переменных.
V 400 200 0
-200 -400 -600
0 25 50 75 100 125 150 175 I, дни Рис. 6. График изменения сопряжённой переменной \уи
1000
800
600
400
200
0 25 50 75 100 125 150 175 t, дни Рис. 7. График изменения сопряжённой переменной
¥я 0,00018
0,00016
0,00014
0,00012
0,00010
0 25 50 75 100 125 150 175 I, дни Рис. 8. График изменения сопряжённой переменной ук
На рис. 9 показан график изменения функции Понтрягина, а на рис. 10 — график изменения управления в .
Н 0,00 -0,02 -0,04 -0,06 -0,08 -0,10 -0,12
0 25 50 75 100 125 150 175 дни Рис. 9. График изменения функции Понтрягина
0
Моделирование и оптимизация перелета спутников малой массы... в, рад 300
250
200
150
100
50
0 25 50 75 1 00 125 150 175 г, дни Рис. 10. График изменения управления
Используя радиус и полярный угол, можно построить траекторию полета космического аппарата, показанную на рис. 11. Для большей наглядности на графике также изображено Солнце в соответствующем масштабе.
у, -1011 м 2,00
1,75
1,50
1,25
1,00
0,75
0,50
0,25
0,00
-1,5 -1,0 -0,5 0,0 0,5 1,0 1,5 х, -1011 м Рис. 11. Траектория полета
Далее проводится проверка корректности решения. Для этого строится совместный график полной и частной производных функции Понтрягина, показанный на рис. 12. Совпадение данных по производным составляет 6 значащих цифр.
Решение было получено за 11 итераций метода Ньютона. Время полета получилось равным около 193 дней.
дН 0,00175
дt '
0,00150
<Н
И 0,00125 0,00100 0,00075 0,00050 0,00025 0,00000
0 25 50 75 100 125 150 175 г, дни
Рис. 12. Совместный график полной и частной производных функции Понтрягина
Заключение. Анализ результатов расчета показал возможность применения метода пристрелки для решения задачи перелета между орбитами планет в солнечной системе, однако, возникла необходимость применения модифицированного метода Ньютона вместо классического, поскольку классический метод не давал сходимости. Метод пристрелки применим и при решении задач оптимального управления с переключением [15], [17]. Кроме этого, применение локальной нормы улучшило сходимость метода. Полученное решение прошло проверку корректности. Результаты решения оказались близкими полученным в [8]. Точность интегрирования системы ДУ была проверена просчетом от конца к началу для сошедшегося варианта. На всех шагах сохранилось совпадение данных по всем переменным с точностью до 6-го знака.
ЛИТЕРАТУРА
[1] Политех. Электрический ракетный двигатель [Электронный ресурс]. Режим доступа: https://polytech.bm.digital/ontology/368214095364644864/el-ektricheskij-raketnyij-dvigatel (дата обращения: 05.03.2021)
[2] Хайтек. Ионная тяга: как человечество использует электрические двигатели для полётов в космос [Электронный ресурс]. Режим доступа: https://hightech.fm/2019/09/24/ion-space (дата обращения: 20.03.2021)
[3] Жарков В.Н., Мороз В.И. Почему Марс? Природа, 2000, № 6, с. 58-67.
[4] Информационное телеграфное агентство России (ИТА-ТАСС). Исследование Марса космическими аппаратами. Досье [Электронный ресурс]. Режим доступа: https://tass.ru/info/5178916 (дата обращения: 21.03.2021).
[5] Моршнева И.В., Овчинникова С.Н. Численное решение краевых задач для обыкновенных дифференциальных уравнений. Метод стрельбы. Методические указания для студентов 3 и 4 курсов мехмата. Ростов-на-Дону, УПЛ РГУ, 2003, 29 с.
[6] Ногин В.Д. Введение в оптимальное управление. Учебно-методическое пособие. Санкт-Петербург, ЮТАС, 2008, 92 с.
[7] Григорьев К.Г., Григорьев И.С, Заплетин М.П. Практикум по численным методам в задачах оптимального управления (дополнение 1). Москва, Изд-во МГУ, 2007, 184 с.
[8] Лейтман Дж. Методы оптимизации с приложениями к механике космического полета. Москва, Наука, 1965, 538 с.
[9] Equity.today. Портал о финансовых рынках. Сколько лететь с Земли до Марса [Электронный ресурс]. Pежим доступа: https://equity.today/polet-na-mars.html (дата обращения: 21.03.2021).
[10] Крайнов А.Ю., Моисеева К.М. Численные методы решения краевых задач для обыкновенных дифференциальных уравнений. Учебное пособие. Томск, STT, 2016, 44 с.
[11] Ахмеров P.P., Cадовский Б.Н. Очерки по теории обыкновенных дифференциальных уравнений [Электронный ресурс]. Pежим доступа: http://www.nsc.ru/rus/textbooks/akhmerov/ode_unicode/index.html (дата обращения: 21.03.2021).
[12] Мудров А.Е. Численные методы для ПЭВМ на языках Бейсик, Фортран и Паскаль. Томск, МП «^CKO», 1991, 272 с.
[13] Демидович Б.П., Марон И.А., Шувалова Е.З. Численные методы анализа. Приближение функций, дифференциальные и интегральные уравнения. Москва, Наука, 1967, 368 с.
[14] Асламова B.C., Колмогоров А.Г., Огупакова Н.Н. Вычислительная математика. Часть первая. Учебное пособие для студентов дневного и заочного обучения технических и химико-технологических специальностей. Ангарск, АГТА, 2003, 82 с.
[15] Федоренко P.^ Приближенное решение задач оптимального управления. Москва, Наука, 1978, 486 с.
[16] Куксенко СП., Газизов ТР. Итерационные методы решения системы линейных алгебраических уравнений с плотной матрицей. Томск, Изд-во ТГУ, 2007, 208 с.
[ 17] Мозжорина Т.Ю. Численное решение задач оптимального управления с переключением методом пристрелки. Математическое моделирование и численные методы, 2017, № 2 (14), с. 94-106.
Cтатья поступила в редакцию 12.05.2021
^втку на эту статью просим оформлять следующим образом:
Мозжорина Т.Ю., Чуванова Л.О. Моделирование и оптимизация перелета спутников малой массы с Земной орбиты на орбиту Марса с помощью ионных двигателей. Математическое моделирование и численные методы, 2021, № 2, с. 54-67.
Мозжорина Татьяна Юрьевна — канд. техн. наук, доцент кафедры «Вычислительная математика и математическая физика» МГТУ им. Н.Э. Баумана. e-mail: mozzhorina@mail.ru.
Чуванова Людмила Олеговна — студент кафедры «Вычислительная математика и математическая физика» МГТУ им. Н.Э. Баумана. e-mail: lusie1031@gmail.com
Optimization of low-mass satellites flight from Earth orbit to Mars orbit using ion engines
© T.Yu. Mozzhorina, L.O. Chuvanova Bauman Moscow State Technical University, Moscow, 105005, Russia
In this paper, optimization of the transfer of a low-mass satellite from Earth orbit to Mars orbit using ion thrusters is considered. The ion engine allows you to minimize fuel consumption and accelerate the spacecraft to fairly high speeds far from the planets of the solar system. The heliocentric section of the flight is subject to consideration. The task is to minimize the flight time. The following assumptions are made in the work: the orbits of the Earth and Mars are circular and lying in the same plane. The angle between the tangential velocity of the spacecraft in the heliocentric system and the direction of thrust action is selected as a control. When compiling the optimization algorithm, the Pontryagin maximum principle was used, which leads the optimization problem of a functional to a boundary value problem for a system of ordinary differential equations. The solution to the boundary value problem was found by one of the numerical methods — the false position method, which gives the most accurate results. The analysis of the results obtained is carried out and a comparison with the data obtained earlier in similar calculations by foreign authors by another numerical solution method is carried out. The conclusion is made about the efficiency of the false position method when solving such problems.
Keywords: ion thrusters, the false position method, boundary value problems of ordinary differential equations, optimal control, Pontryagin's maximum principle, flight between the orbits of the Earth and Mars
REFERENCES
[1] Polytech. Electric rocket engine [Politekh. Elektricheskij raketnyj dvigatel'] [Electronic resource]. Access mode: https://polytech.bm.digital/ontology/ 368214095364644864/elektricheskij -raketnyij-dvigatel (accessed: 05.03.2021).
[2] Hajtek. Ionnaya tyaga: kak chelovechestvo ispol'zuet elektricheskie dvigateli dlya polyotov v kosmos [Hi-tech. Ion thrust: How humanity uses electric motors to fly into space] [Electronic resource]. Access mode: https://hightech.fm/2019/09/24/ion-space (accessed: 20.03.2021).
[3] Zharkov V.N., Moroz V.I. Pochemu Mars? [Why Mars?] Priroda [Nature], 2000, no. 6, pp. 58-67.
[4] Russian news agency. IssledovanieMarsa kosmicheskimi apparatami. Dos'e [Exploration of Mars by spacecraft. Dossier] [Electronic resource]. Access mode: https://tass.ru/info/5178916 (accessed: 21.03.2021).
[5] Morshneva I.V., Ovchinnikova S.N. Chislennoe reshenie kraevyh zadach dlya obyknovennyh differencial'nyh uravnenij. Metod strel'by. Metodicheskie ukaza-niya dlya studentov 3 i 4 kursov mekhmata [Numerical solution of boundary value problems for ordinary differential equations. The method of shooting. Methodological guidelines for students of the 3rd and 4th courses of mehmat]. Rostov-on-Don, UPL RSU Publ., 2003, 29 p.
[6] Nogin V.D. Vvedenie v optimal'noe upravlenie. Uchebno-metodicheskoe posobie [Introduction to optimal control. Educational and methodical manual]. Saint Petersburg, UTAS Publ., 2008, 92 p.
[7] Grigoriev K.G., Grigoriev I.S., Zapletin M.P. Praktikumpo chislennym metodam v zadachah optimal'nogo upravleniya (dopolnenie 1) [Workshop on numerical methods in optimal control problems (Supplement 1)]. Moscow, MSU Publ., 2007, 184 p.
[8] Lietmann G. Optimization techniques: with applications to aerospace systems. Academic Press, 1962, 453 p.
[9] Equity.today. Portal o finansovyh rynkah. Skol'ko letet's Zemli do Marsa [Eq-uity.today. Portal about financial markets. How long to fly from Earth to Mars] [Electronic resource]. Access mode: https://equity.today/polet-na-mars.html (accessed: 21.03.2021).
[10] Krainov A.Yu., Moiseeva K.M. Chislennye metody resheniya kraevyh zadach dlya obyknovennyh differencial'nyh uravnenij. Uchebnoe posobie [Numerical methods for solving boundary value problems for ordinary differential equations. Study guide]. Tomsk, STT, 2016, 44 p.
[11] Akhmerov R.R., Sadovsky B.N. Ocherkipo teorii obyknovennyh differencial'nyh uravnenij [Essays on the theory of ordinary differential equations] [Electronic resource]. Access mode: http://www.nsc.ru/rus/textbooks/akhmerov/ode _unicode/index.html (accessed: 21.03.2021).
[12] Mudrov A.E. Chislennye metody dlya PEVMna yazykah Bejsik, Fortran i Paskal' [Numerical methods for PCs in Basic, Fortran and Pascal languages]. Tomsk, MP "RASKO", 1991, 272 p.
[13] Demidovich B.P., Maron I.A., Shuvalova E.Z. Chislennye metody analiza. Prib-lizhenie funkcij, differencial'nye i integral'nye uravneniya [Numerical methods of analysis. Approximation of functions, differential and integral equations]. Moscow, Nauka Publ., 1967, 368 p.
[14] Aslamova V.S., Kolmogorov A.G., Stupakova N.N. Vychislitel'naya matematika. Chast' pervaya. Uchebnoe posobie dlya studentov dnevnogo i zaochnogo obucheniya tekhnicheskih i himiko-tekhnologicheskih special'nostej [Computational mathematics. Part I. A textbook for full-time and part-time students of technical and chemical-technological specialties]. Angarsk, ASTA Publ., 2003, 82 p.
[15] Fedorenko R.P. Priblizhennoe reshenie zadach optimal'nogo upravleniya [Approximate solution of optimal control problems]. Moscow, Nauka Publ., 1978, 486 p.
[16] Kuksenko S.P., Gazizov T.R. Iteracionnye metody resheniya sistemy linejnyh al-gebraicheskih uravnenij splotnoj matricej [Iterative methods for solving a system of linear algebraic equations with a dense matrix]. Tomsk, TSU Publ., 2007, 208 p.
[ 17] Mozzhorina T.Yu. Numerical solution to problems of optimal control with switching by means of the shooting method. Mathematical modeling and Computational Methods, 2017, no. 2 (14), pp. 94-106.
Mozzhorina T.Yu., Cand. Sc. (Eng.), Assoc. Professor of Department of Computational Mathematics and Mathematical Physics, Bauman Moscow State Technical University. e-mail: mozzhorina@mail.ru.
Chuvanova L.O., Student of Department of Computational Mathematics and Mathematical Physics, Bauman Moscow State Technical University. e-mail: lusie1031@gmail.com.
