Научная статья на тему 'РАСЧЕТ ТРАЕКТОРИИ МЕЖПЛАНЕТНОГО ПЕРЕЛЕТА ЗЕМЛЯ-МАРС С МАЛОЙ ТЯГОЙ БЕЗ ИСПОЛЬЗОВАНИЯ МЕТОДА ГРАВИСФЕР'

РАСЧЕТ ТРАЕКТОРИИ МЕЖПЛАНЕТНОГО ПЕРЕЛЕТА ЗЕМЛЯ-МАРС С МАЛОЙ ТЯГОЙ БЕЗ ИСПОЛЬЗОВАНИЯ МЕТОДА ГРАВИСФЕР Текст научной статьи по специальности «Физика»

CC BY
172
11
i Надоели баннеры? Вы всегда можете отключить рекламу.
Журнал
Труды МАИ
ВАК
Область наук
Ключевые слова
ТРАЕКТОРИЯ МЕЖПЛАНЕТНОГО ПЕРЕЛЕТА / МЕТОД ГРАВИСФЕР / ОПТИМИЗАЦИЯ / ПРИНЦИП МАКСИМУМА

Аннотация научной статьи по физике, автор научной работы — Синицын Алексей Андреевич

Рассматривался перелет с низкой круговой околоземной орбиты на низкую круговую околомарсианскую орбиту с двигательной установкой малой тяги при постоянных по величине тягой и удельным импульсом тяги. Для расчета траектории перелета разработана методика прямого расчета траектории (без использования метода сфер действия) с учетом гравитационного воздействия (в рамках ньютоновской модели) на аппарат Земли, Солнца и Марса. Направление вектора тяги определялось в рамках вариационной постановки с функционалом продолжительность перелета (задача быстродействия) с использованием необходимых условий оптимальности в форме принципа максимума Понтрягина. Применение принципа максимума позволило заменить вариационную задачу на трехточечную краевую задачу, которая решалась численно. Приводится сопоставление продолжительности перелета и конечной массы, рассчитанных с использованием разработанной методики и с использованием метода сфер действия.

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

CALCULATION OF INTERPLANETARY EARTH-MARS LOW THRUST TRANSFER WITHOUT METHOD OF GRAVISPHERES

The article considers the spacecraft interplanetary transfer from circular low Earth orbit to circular low orbit near Mars by propulsion system with constant values of thrust and specific impulse. The main goal of the article consists in length gravispheres method inaccuracy estimation. Due to the nearness of planets’ phases, the simulation revealed substantial for heliocentric phase value of excessive hyperbolical velocity on the edge of gravispheres. while applying conventional technique, imposing the value of excessive hyperbolical velocity for the trajectory conjugate phases and the planets’ heliocentric phase. These simplifications affect the transfer effectiveness figures, namely, transfer duration and final mass. Another factor, defining the value of gravispheres method inaccuracy consists in neglecting the long duration of gravitational attraction of the Sun on planet phases, as well as of the planets on heliocentric phase of spacecraft trajectory. The article presents the developed technique for direct determination of spacecraft trajectory (without applying gravispheres method) including Earth, Sun and Mars gravitational attraction effect on a spacecraft (as Newton’s dynamic model). Thrust vector direction is defined from solving variation problem with transfer duration as performance index (minimum time problem) applying necessary optimality condition in the form of Pontryagin’s maximum principle. Maximum principle application allowed simplify variation problem to the three points boundary value problem, solved numerically. Comparison of transfer performance indexes, i. e. duration and final spacecraft mass, obtained by the developed technique implementation, and application of gravispheres method is presented. In the last case on planet phases variants of parabolic velocity and hyperbolic excess velocity with value similar to the direct solution (without gravispheres method applying) is considered. Presented computation results revealed worsening of transfer duration and final mass in the case of applying gravispheres method vs. direct trajectory calculation (4-11% for duration and 2-6% for mass).

Текст научной работы на тему «РАСЧЕТ ТРАЕКТОРИИ МЕЖПЛАНЕТНОГО ПЕРЕЛЕТА ЗЕМЛЯ-МАРС С МАЛОЙ ТЯГОЙ БЕЗ ИСПОЛЬЗОВАНИЯ МЕТОДА ГРАВИСФЕР»

Труды МАИ. Выпуск № 94

http://trudymai.ru/

УДК 629.788:523.43

Расчет траектории межпланетного перелета Земля-Марс с малой тягой без использования метода грависфер

Синицын А.А.

ГНЦ ФГУП «Центр Келдыша», ул. Онежская, 8, Москва, 125438, Россия e-mail: AlSinitsin@yandex.ru

Аннотация

Рассматривался перелет с низкой круговой околоземной орбиты на низкую круговую околомарсианскую орбиту с двигательной установкой малой тяги при постоянных по величине тягой и удельным импульсом тяги. Для расчета траектории перелета разработана методика прямого расчета траектории (без использования метода сфер действия) с учетом гравитационного воздействия (в рамках ньютоновской модели) на аппарат Земли, Солнца и Марса. Направление вектора тяги определялось в рамках вариационной постановки с функционалом продолжительность перелета (задача быстродействия) с использованием необходимых условий оптимальности в форме принципа максимума Понтрягина. Применение принципа максимума позволило заменить вариационную задачу на трехточечную краевую задачу, которая решалась численно. Приводится сопоставление продолжительности перелета и конечной массы, рассчитанных с использованием разработанной методики и с использованием метода сфер действия.

Ключевые слова: траектория межпланетного перелета, малая тяга, метод грависфер, оптимизация, принцип максимума, задача быстродействия.

Введение

Использование ядерных электроракетных двигательных установок (ЯЭРДУ) в межпланетных полетах является перспективным направлением развития ракетно-космической техники, так как это позволит доставку высокоэнерговооруженной полезной нагрузки большой массы к различным планетам Солнечной системы. ЯЭРДУ, как правило, относят к двигательным установкам малой тяги. Для аппаратов с ЯЭРДУ траектория межпланетного перелета состоит из трех последовательных фаз: раскрутка у планеты отлета, гелиоцентрический участок полета, скрутка у планеты назначения. Такое деление как правило связано с применяемой методикой поиска траектории межпланетного перелета, подразумевающей использование метода грависфер (или, по другому, сфер действия) нулевой протяженности.

Метод грависфер нулевой протяженности широко используется в расчетах траекторий межпланетных перелетов с малой тягой (смотри, например работу [1]). Часто используемое положение об использовании метода грависфер нулевой протяженности и нулевых гиперболических избытках скорости для стыковки участков, принимаемое в качестве упрощения, позволяет при использовании принципа максимума получить формулировку краевых задач перелета в рамках ограниченной задачи двух тел, в решении которых накоплен определенный опыт.

Однако, при всей распространенности, использование метода грависфер нулевой протяженности связано с разного рода неопределенностями. Одним из требующих прояснения вопросов является оценка методической погрешности от использования метода грависфер нулевой протяженности. Для получения корректной оценки методической погрешности от применения этих упрощений требуется сопоставление расчета траектории, сделанного с использованием упрощений, с расчетом, их не использующим. Методическая ошибка от использования грависфер нулевой протяженности (по сравнению с методикой расчета траектории без упрощений, связанных с разбиением траектории на участки, соответствующие грависферам) включает погрешность, связанную с влиянием возмущающих ускорений (на припланетных участках главным образом от гравитационного воздействия Солнца), а также погрешность, определяемую неоптимальностью условий стыковки участков межпланетной траектории (главным образом -использование нулевых или любых других заданных гиперболических избытков скорости на концах гелиоцентрического участка).

Последняя составляющая методической ошибки от использования метода грависфер может влиять на требуемые для организации межпланетного перелета энергозатраты, так как принятые упрощения условий стыковки участков перелета определяют краевые условия гелиоцентрического участка траектории.

1. Решение задачи разгона у Земли

В работе [2] приведены постановка и результаты задачи минимизации

продолжительности разгона до достижения параболической скорости. При замене

условия достижения параболической скорости на условие достижения заданного ненулевого гиперболического избытка скорости все остальные условия краевой задачи остаются неизменными. На рисунке 1 показаны результаты расчетов отлетных от Земли траекторий с низкой круговой орбиты высотой 300 км с различным значением требуемого гиперболического избытка скорости в виде зависимости характеристической скорости и радиуса в конце траектории от задаваемого значения гиперболического избытка скорости. Расчеты проведены при начальном ускорении аппарата а0 =3,72-10-4 м/с2 и удельном импульсе тяги I д =1500 с.

Характеристическая скорость раскрутки, м/с 9500

9000 8500 8000 7500 7000 6500

скор )СТЬ ч

радиус

Конечный радиус, км 6 000 000

5 000 000 4 000 000 3 000 000 2 000 000 1 000 000 о

0 500 1000 1500 2000 2500 Гиперболический избыток скорости, м/с

Рис. 1. Зависимость набора характеристической скорости при отлете от Земли и расстояния до центра Земли от гиперболического избытка скорости. Начальная

орбита - круговая высотой 300 км.

Из рисунка 1 видно, что увеличение гиперболического избытка скорости с нуля до 500 м/с приводит к увеличению характеристической скорости лишь на 1,6%. При этом радиус сферы действия Земли (~925 000 км) достигается при гиперболическом

избытке 700-800 м/с. Следует заметить, что возможность получения гиперболического избытка скорости при разгоне в сфере действия Земли хорошо известна (смотри, например, работу [3]). Из результатов, приведенных на рисунке 1, следует, что существует принципиальная возможность обеспечения больших величин гиперболического избытка скорости (порядка нескольких км/с и более) практически вне зависимости от характеристик аппарата и двигательной установки. При этом радиус в конце раскрутки определяется требуемой величиной гиперболического избытка скорости, а также характеристиками аппарата - его ускорением и удельным импульсом тяги. Таким образом, при использовании метода сфер действия выбор величины радиуса сферы действия планеты связан с оптимальными величинами гиперболического избытка скорости и методической погрешностью.

Следует отметить, что возможность выбора величины гиперболического избытка скорости в некоторых случаях будет существенно изменять концевые условия краевой задачи и, по-видимому, может существенным образом влиять на результаты решения краевых задач. Действительно, при ограниченной стартовой массе и нулевых гиперболических избытках часто требуется дополнительный виток, чтобы обеспечить выравнивание скоростей аппарата и планеты назначения. Если же можно сформировать краевую задачу с ненулевой скоростью аппарата относительно планеты назначения, то дополнительный виток может не потребоваться, что означает существенное сокращение продолжительности гелиоцентрического участка полета.

В процессе решения краевых задач с заданным значением гиперболического

избытка скорости выявилась особенность, заключающаяся в наличии

неединственности решения. На рисунке 2 приведены траектории разгона с низкой

околоземной орбиты высотой 300 км до гиперболического избытка скорости =500 м/с. Начальное ускорение аппарата и удельный импульс тяги приняты а0 =2.44-10-4 м/с2 и I уд =3000 с, соответственно.

Экстремали на рисунке 2 отличаются угловой дальностью перелета (количество витков для варианта а) - 2774,3; б) - 2775,5; в) - 2778,3), однако в данном случае многоэкстремальность не может быть обусловлена только различием угловой дальности перелета, так как например для семейства экстремалей варианта а) можно получить аналогичные траектории с заданной угловой дальностью, как для экстремалей вариантов б) и в). В этой связи возникает вопрос о классификации экстремалей. Характеристическая скорость, требуемая на осуществление перелета, составляет для варианта а) - 7391 м/с, для б) - 8392 м/с, для в) - 9844 м/с. Таким образом «лишние» экстремали для задачи быстродействия могут быть отбракованы по значению функционала. Однако, такая многоэкстремальность существенным образом влияет на сходимость при решении краевых задач - практика расчетов показывает, что область притяжения в пространстве начальных значений сопряженных переменных у экстремали с минимальной характеристической скоростью по-видимому меньше, чем у других экстремалей.

а)

б)

в)

Рис. 2. Варианты экстремалей траекторий отлета. Начальная орбита - круговая высотой 300 км. =3000 с, а0 =2.44-10-4 м/с2. Набор гиперболического избытка

скорости =500 м/с. 7

Неизвестно связано ли появление подобной многоэкстремальности на припланетном участке траектории с применением метода грависфер, либо варианты оптимальных траекторий с различной фазой отлета существуют вне зависимости от принятых упрощений.

Следует отметить, что аналогичная (симметричная) картина должна быть и при подлете к планете назначения.

2. Методика сквозного расчета траектории межпланетного перелета

Сквозная оптимизация траекторий межпланетного перелета (без использования метода грависфер) требует составления специальных постановки краевой задачи и методики расчета, так как практика расчетов показала, что только учет дополнительных гравитационных ускорений при краевых условиях слева - точка на орбите планеты отлета - и справа -точка на орбите планеты назначения - приводит к невозможности решить краевую задачу (не удается получить скрутки у планеты назначения). Ниже представлена методика сквозного расчета межпланетной траектории с малой тягой, в которой отсутствует указанная выше сложность.

Рассматривалась задача перелета с заданной орбиты около планеты отлета на заданную орбиту около планеты назначения. При этом полагалось, что заданными в начале и в конце являются все элементы кроме аргумента широты и долготы восходящего узла, которые подлежат определению. Дата старта также была принята заданной. В качестве функционала рассматривалась продолжительность перелета (задача быстродействия), а неизвестными управления полагались законы ориентации

вектора тяги. При расчете траектории законы ориентации вектора тяги определялись согласно принципу максимума Понтрягина.

При постановке и в решении краевой задачи расчета траектории межпланетного перелета использовался метод многоточечной стрельбы. Его применение к рассматриваемому случаю в значительной мере схоже с подходами работ [4-5] и состояло в следующем.

Известно, что задача скрутки (спуск по многовитковой спирали с

параболической или гиперболической пролетной траектории на опорную круговую

орбиту) является плохо обусловленной задачей (то есть приводит к итерационному

процессу с плохой сходимостью). В то же время обратная задача раскрутки не

обладает этим негативным свойством. Типичный перелет с малой тягой между

круговыми орбитами около планет отлета и назначения представляет собой раскрутку

у планеты отлета, продолжительный гелиоцентрический участок полета и скрутку у

планеты назначения. Тогда фазовые координаты в начале и в конце траектории

известны (с точностью до аргументов широты и долгот восходящего узла).

Использованный в настоящем исследовании подход заключался в разбиении

траектории на два участка произвольно по времени полета так, чтобы момент гх,

разделяющий оба участка, приходился на гелиоцентрическую часть полета. Первый

участок траектории рассчитывался с использованием приведенных выше уравнений

движения в геоцентрической системе координат относительно планеты отлета.

Второй участок рассчитывался с применением уравнений движения в системе

координат относительно планеты назначения, причем интегрирование проводилось с

отрицательным шагом по независимой переменной. В рассматриваемом случае

9

неизвестными краевой задачи являются сопряженные переменные и часть фазовых на левом и правом концах траектории перелета. Невязку составляет разница фазовых и сопряженных переменных в момент ^, посчитанных интегрированием слева и справа. Большим достоинством рассматриваемого подхода является простота подбора начального значения сопряженных переменных, обеспечивающих скрутку у планеты назначения и раскрутку у планеты отлета. Отрицательным моментом является увеличение размерности краевой задачи вдвое.

Для расчета траектории межпланетного перелета использовались сразу две декартовы эклиптических системы координат с началом отсчета в центрах планет отлета и назначения. Использование декартовых эклиптических систем координат связано с тем, что уравнения движения и им сопряженные имеют в декартовой системе координат наиболее простой вид.

Уравнения движения относительно планеты отлета приведены ниже. Уравнения движения относительно планеты назначения и им сопряженные имеют аналогичный вид (с точностью до замены соответствующих индексов).

(1)

• К , К - Р С + Р 0 - К - Р С + Р к

у = -Р0 —Т-/лС \-77-М«\-7Т+а

К |К - Р С + Р о| |К - Р С + Р к| (2)

Принятые выше и далее обозначения: К и V - вектор положения и вектор скорости КА относительно планеты отлета, а р 0, р к, р с, - вектора положения планеты отлета, планеты назначения и Солнца относительно барицентра солнечной системы,

Мс Мк ^

соответственно; /ис =— ; /ик =— ; м с - гравитационная постоянная Солнца;

Мо Мо

¡и0 - гравитационная постоянная планеты отлета; /ик - гравитационная постоянная планеты назначения; а - вектор ускорения, обеспечиваемого ДУ

Эфемериды планет рассчитывались с использованием модели ЕРМ2008 [6]. Следует отметить, что при использовании модели эфемерид ЕРМ2008 координаты, скорости и ускорения являются явными функциями времени: р 0 = р 0 (г).

К уравнениям движения целесообразно добавить уравнение для изменения массы

т

и времени г. Таким образом уравнения движения являются неавтономными. Для приведения уравнений движений в автономный вид целесообразно введение новой независимой переменной т.

т=- т 8

(3)

йг йт

(4)

Здесь 8 - функция включения/выключения ДУ (принимает значения 0 когда ДУ выключена и 1 когда ДУ включена).

Вышеприведенные уравнения движения представлены в безразмерном виде для удобства использования метода интегрирования с автоматическим выбором шага (использовался метод Рунге-Кутта-Фельберга 7(8) порядка) [7,8]. Обезразмеривание переменных было проведено как в работе [2].

Координаты вектора ускорения а =

с \ а

V аг У

, могут быть выражены через тягу к д и

углы ориентации вектора тяги у и а (у - угол между проекцией вектора тяги на

а

у

опорную плоскость и выбранным направлением оси ОХ; а - угол между направлением вектора тяги и опорной плоскостью) следующим образом.

Я,

Я

Я,

а х =-д-Созу-Соза а =-д-Бту-Соза а 2 =-д-Бта

т

т

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

т

(5)

Углы ориентации вектора тяги у и а, а также функция включения/выключения д ДУ подлежат выбору как законы управления из условия минимизации принятого функционала (задача быстродействия).

Гамильтониан системы имеет следующий вид.

Н=рк-У+р.

.. Я , К р С +Р 0 - К"Р С +Р к ^ р 0 с--—-V к--—+а

Я

К р С +р 0

К р С +р к

+Рт ( - т д )-Рг

У

(6)

Здесь рк и ру- вектора сопряженных переменных координатам и скоростям. Последний имеет координаты рг , ру , рг .

Из принципа максимума можно получить следующие соотношения для определения законов управления.

Бту=

/2 2 ^ +Руу

Соз у=■

Р•

(7)

Бта=-

2 2 2 Ру +Ру + Р-

Соза=-

V у г

/

2,2,2 Ру +Ру +Ру

(8)

Для задачи быстродействия д=1 на всей траектории движения [9].

Уравнения для сопряженных переменных:

Р

V

у

у

у

2

Pv Pv ^ „ pv ( R-R С )

Р R = — -3 R-—+v -t--3 v c ( R-R С h-—

R3 R5 --Л3 --Л5

(lR-R с I)3 (I R-R с I)

Pv э - /О O ,Pv (R R к )

+// « Ti-Tv " к ( R-R к ^-лГ

(l R-R к 0 (| R-R к \)

PR

(10)

. R д P m = TI Р

m

(11)

Необходимо отметить, что полученная запись сопряженных уравнений содержит в правых частях составляющие, имеющие возведение расстояния от КА до притягивающего центра в пятую степень. В этом случае выбор характерного расстояния для обезразмеривания переменных является значимым для организации вычислительного процесса (неправильный выбор характерного расстояния вблизи планет отлета и назначения приводит к ситуации деления на малое число, близкое нулю, что вызывает останов выполнения программы). Например, при выборе характерного размера равным одной астрономической единице безразмерное

7171

расстояние на орбите высотой 800 км у Земли -«510 5. Возведение этой

150-10б

величины в пятую степень ведет к значению 3 10-22. Следует отметить, что точность представления вещественного числа (рассчитана согласно алгоритму, приведенному в работе [10]) составляет ~ 5 10 16 для стандартной точности машинного представления (восемь байт на число) и ~ 5 10 20 для двойной точности машинного представления (10 байт на число). Таким образом, использование декартовых уравнений движения относительно барицентра Солнечной системы и выбор в качестве характерного расстояния астрономической единицы будет приводить к

V

потере точности при расчете правых частей сопряженных уравнений вблизи планет назначения и отлета.

Вышеприведенные вычислительные трудности обусловили использование для расчета траектории межпланетного перелета декартовых эклиптических системы координат с началом отсчета в центрах планет отлета и назначения. Характерные расстояния было выбраны соответствующими радиусам начальной орбиты у планеты отлета и конечной орбиты у планеты назначения (хотя такой выбор характерных расстояний не вполне удачен, так как приводит к значительным изменениям правых частей уравнений движения и сопряженных на концах траектории).

Использование сопряженных переменных в декартовой системе координат даже в случае ограниченной задачи двух тел приводит к сложностям с построением продолжения решения по выбранному параметру исходных данных от найденного решения краевой задачи. Для ограниченной задачи двух тел в значительной степени эту вычислительную трудность преодолевает переход от использования декартовой системы координат к набору орбитальных элементов. В настоящее время большое распространение получили равноденственные элементы к , вх , е у , ¡х , 1у , ^ [11,12].

Выражения для уравнений движения и им сопряженных в равноденственных элементах при наличии трех гравитационных центров являются достаточно громоздкими, а также требуют большого объема отладочных и верификационных работ. Поэтому, в расчетах невязки и неизвестных краевой задачи использовались равноденственные элементы, а интегрирование дифференциальных уравнений для фазовых и сопряженных переменных проводилось в декартовых систем координат, как описано выше.

Необходимый перевод между фазовыми переменными и им сопряженными в различных системах координат (декартовой и равноденственных элементах) с различными характерными множителями выполнялись с использованием канонических преобразований [13,14].

Так как начальные значения долготы восходящего узла О0 и аргумента широты и0, а также соответствующие конечные Ок и ик являются свободными, то необходимо воспользоваться условиями трансверсальности для замыкания краевой задачи (как, например, в работе [12]). Следует отметить, что условия трансверсальности в случае круговых орбит несколько упрощаются, по сравнению с приведенными в работе [12]. Все фазовые переменные слева и справа заданы, кроме О0, и0,0к и ик.

Из условия нормировки: у 0=1. Из условия трансверсальности у 0 = 0 и у К = 0.

Неизвестные краевой задачи (14 штук): у0 , у0 , у0 , у0 , ^0, у К, ук , ук , ук ,

е х е у г х I у е х е у I х

ук , ^к , О0 ,0к , гк. Здесь под уj понимается сопряженная переменная равноденственному элементу j из набора { И, е х, еу, г х, г у, ^ }.

Так как нормирующие множители при интегрировании слева и справа могут различаться при расчете невязки в момент гх вводится дополнительный нормирующий множитель следующим образом:

у 1

у 1 =",

2 , 2 , 2, 2, 2, 2 I у . + у + у + у . + у . +у _

(12)

2 2 2 2 2 2 у/у2 +у ех +уеу +угх +у2 +у F

где j - любой из рассматриваемых равноденственных элементов. Компоненты вектора невязки £:

5=

к (г - )-к (г + ) е х (г 1- )-ех (г 1+ ) е у (г 1- )-е у (г 1+ )

¿X ( г 1- )-'* ( г 1+ )

1у (г - Ну (г 1+ )

Е (г - )-Е (г + ) ¥ к (г - )-¥ к (г 1+ ) ¥ ех ( г Г )-¥ ( г 1+ )

¥ е, (г Г )-¥ ^ (г 1+ )

¥и (г 1- )-¥,х (г 1+ )

(г 1- )-¥ ^ (г 1+ )

¥ е (г 1- )-¥ F (г 1+ ) 0-0 0-0 ¥ 1_ ¿0 ¥ **

¥ к 4К -¥ к -1К

¿х У Г ¿у х

(13)

Здесь знаки «-» и «+» обозначают принадлежность переменной к интегрируемому участку траектории («-» - слева, расчет от планеты отлета; «+» - справа, расчет от планеты назначения).

Величина вектора невязки минимизировалась с использованием метода

Ньютона и градиентного спуска.

3. Сопоставление результатов расчетов с использованием метода грависфер

нулевой протяженности и без

Рассматривался перелет Земля-Марс, соответствующий начальному этапу марсианской экспедиции. Расчет траектории перелета проводился без разбиения ее на участки, соответствующие сферам действия планет, согласно приведенной выше методики. При этом на всей траектории учитывалось притяжение Солнца, Земли и Марса как ньютоновское.

Исходные данные для расчетов приняты следующими. Начальное ускорение полагалось равным 9,56 -10-4м/с2 (что соответствует одному из вариантов марсианской пилотируемой экспедиции с начальной массой экспедиционного комплекса 450 т и мощностью ЯЭРДУ 24 МВт), а удельный импульс тяги 7000 с. Стартовая околоземная орбита полагалась круговой высотой 1000 км и наклонением 51,6°. Конечная околомарсианская орбита полагалась также круговой высотой 400 км и наклонением 45° (в качестве опорной рассматривалась плоскость эклиптики). Долгота восходящего узла начальной и конечной орбит полагалась свободной и определялась. Дата старта для определенности полагалась 01.06.2035 12:00:00.

Поставленная выше краевая задача оказалась плохо обусловленной и была

решена лишь при одном наборе исходных данных, приведенном выше. При этом

плохая сходимость (происходили отказ метода Ньютона и медленный спуск вдоль

направления антиградиента с малым рабочим шагом) наблюдалась даже в достаточно

близкой окрестности решения. Количество итераций при решении краевой задачи (и,

соответственно, количество численных расчетов матрицы Якоби и вектора градиента)

можно приблизительно оценить в несколько тысяч. Найденные неизвестные краевой

17

задачи: у0 =-0.004614324, у0 =-0.002137668, у0 =-0.002315907, у0 =1.43098Е-05,

е х е у ^ х г у

^0 =4.707837017, у * =-0.099632077, ук =-0.000601248, ук =0.000929648,

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

е х е у

ук =-1.6551Е-05, ук =5.6059Е-06, ^к =0.778730631, О0 =-0.006178844,

I х I у

Ок =-0.32657648, гк =300.5914423 суток.

Полученное решение может быть использовано для оценки методической погрешности от введения сфер действия. Однако для получения полной картины единичных расчетов недостаточно, а требуется решение ряда аналогичных краевых задач с вариацией исходных данных, так как разница в величинах показателей энергобаллистической эффективности для различных методик расчета может быть объяснена не только методической погрешностью от введения сфер действия. Другими объяснениями могут являться: многоэкстремальный характер решений краевой задачи перелета без использования сфер действия (в этом случае для различных семейств решений будет разница в энергобаллистических показателях); отличие от оптимальных гиперболических избытков в точках пересечения траектории со сферами действия для расчетов с использованием метода грависфер. Также, такая разница получается при непрерывно работающей ЭРДУ Введение пассивных участков также может изменить различие в результатах при сравнении методик.

Результаты расчетов основных показателей траектории перелета при использовании приведенной методики построения траектории с учетом гравитационного действия Солнца, Земли и Марса в сравнении с аналогичными показателями, полученными с использованием метода грависфер нулевой протяженности, показаны в таблице 1. При этом для метода грависфер нулевой

протяженности сделано два расчета. Расчет 1 проведен при условии нулевых гиперболических избытков скорости при стыковке участков. Параметры припланетных участков перелета определялись согласно аналитической звисимости из работы [2]. Расчет 2 проведен при тех же гиперболических избытках скорости, что и полученные при сквозном расчете значения. При этом характеристики припланетных участков рассчитывались из решения соответствующих краевых задач, а на гелиоцентрическом участке гиперболические избытки скорости полагались коллинеарными базис-вектору Лоудена.

Таблица 1. Основные количественные показатели траекторий, рассчитанных различными способами.

Сквозной расчет с Расчет 1 Расчет 2

учетом гравитации с использованием с использованием

трех небесных тел метода грависфер метода грависфер

на всей траектории нулевой нулевой

полета протяженности протяженности

Начальная масса, т 450 450 450

Стартовая дата 01.06.2035 01.06.2035 01.06.2035

Продолжительность раскрутки у Земли, сутки 81.81 76.9 81.8

Масса КА после раскрутки у Земли, т 405.51 408.2 405.5

Гиперболический избыток 10781

скорости после раскрутки у Земли, м/с 0 1078

Дата отлета из сферы действия Земли 22.08.2035 17.08.2035 22.08.2035

Продолжительность 189.62

гелиоцентрического 234.2 202.6

участка, сутки

Масса в конце

гелиоцентрического 302.52 281 295.5

участка, т

Гиперболический избыток 11903

скорости до скрутки у Марса, м/с 0 1190

Продолжительность скрутки у Марса, сутки 29.13 21.7 28.4

Суммарная

продолжительность 300.6 332.9 312.8

перелета, сутки

Конечная масса, т 286.8 269.2 280.1

Изменение

продолжительности перелета по сравнению со - +10.7% +4.1%

сквозным расчетом

Изменение конечной

массы по сравнению со - -6.1% -2.3%

сквозным расчетом

1 На границе грависферы на расстоянии 930 тыс. км от центра Земли.

2 Под гелиоцентрическим участком здесь понимался участок полета между границами грависфер Земли (930 тыс. км) и Марса (580 тыс. км).

3 На границе грависферы на расстоянии 580 тыс. км от центра Марса.

Следует заметить, что оптимизация величины гиперболических избытков скорости на границе грависфер планет при использовании метода грависфер нулевой протяженности представляет собой определенную трудность, так как для этого требуется реализация методики сквозного моделирования всех трех участков межпланетного перелета. По видимому сложность такой методики не меньше, чем у изложенной выше методики сквозного моделирования без использования метода сфер действия.

Как видно из таблицы 1, имеются расхождение результатов для различных способов расчета, которые могут быть существенны для дальнейшего анализа. Так, например, сравниваемые методики при использовании условия параболической стыковки участков в методе грависфер нулевой протяженности (расчет 1) дают расхождение более чем в 10% по продолжительности межпланетного перелета и почти 6% по конечной массе. Такое различие в показателях энергобаллистической эффективности обусловлено различием в траектории гелиоцентрического участка, что показано на рисунке 3.

В случае (расчет 2), когда условия параболической стыковки не используются, а гиперболические избытки скорости равны величинам, полученным при сквозном расчете на границах грависфер Земли и Марса, расхождение со сквозным расчетом значительно меньше. Так по продолжительности межпланетного перелета оно составляет 4%, а по конечной массе лишь 2%. Как видно из рисунка 3, при этом наблюдается большее совпадение с траекторией без использования метода грависфер.

Рис. 3. Проекция на плоскость эклиптики траекторий перелета Земля-Марс, рассчитанных по различным методикам.

Выводы

Таким образом, использование метода грависфер приводит к получению несколько ухудшенных величин показателей энергобаллистической эффективности перелета (продолжительности перелета и конечной массы) по сравнению с расчетом без использования метода грависфер.

Библиографический список

1. Константинов М.С., Леб Х., Петухов В.Г., Попов Г.А.. Проектно-баллистический анализ пилотируемой марсианской миссии с ядерной электроракетной двигательной установкой // Труды МАИ. 2011. № 42. URL: http: //www. mai.ru/science/trudy/pub-lished.php?ID=24274

2. Лебедев В.Н. Расчет движения космического аппарата с малой тягой. - Сер. Математические методы в динамике космических аппаратов. Вып. № 5. - М.: Вычислительный центр АН СССР, 1968. - 108 с.

3. Ефимов Г.Б. Оптимальный разгон в центральном поле до гиперболических скоростей // Космические исследования. 1970. Т. VIII. № 1. С. 26-47.

4. Vadali S.R., Nah R., Braden E., Jhonson Jr.I.L. Fuel-Optimal Planar Earth-Mars Trajectories Using Low-Thrust Exhaust-Modulated Propulsion // Journal of Guidance, Control and Dynamics, May-June 2000, vol. 23, no. 3. pp. 476-482.

5. Григорьев И.С., Заплетин М.П., Самохин А.С., Самохина М.А. Оценка возможного выигрыша по массе при использовании двигателей малой тяги в экспедиции к Марсу // XI Всероссийский съезд по фундаментальным проблемам теоретической и прикладной механики, Казань, 20-24 августа 2015: сборник докладов. С.1063-1065.

6. Питьева Е.В. и др. Эфемериды EPM2008. URL: ftp://quasar.ipa.nw.ru/incom-ing/EPM/EPM2008

7. Хайрер Э., Нерсетт С., Ваннер Г. Решение обыкновенных дифференциальных уравнений. Нежесткие задачи. - М.: Мир, 1990. - 512 c.

8. Fehlberg E. Classical Fifth-, Sixth-, Seventh-, and Eighth-Order Runge-Kutta Formulas with Stepsize // Control NASA TR R-287, October 1968, 88 p.

9. Petukhov V.G. Minimum-Thrust Problem and its Application to Trajectory Optimization with Thrust Switchings. IAC-13-C1.6.2. // 64th International Astronautical Congress, Beijing, China, 23-27 September 2013.

10. Форсайт Дж., Мальколм М., Моулер К. Машинные методы математических вычислений. - М.: Мир, 1980. - 280 с.

11. Гришин С.Д., Захаров Ю.А., Оделевский В.К. Проектирование космических аппаратов с двигателями малой тяги. - М.: Машиностроение, 1990. - 224 с.

12. Петухов В.Г. Оптимизация многовитковых перелетов между некомпланарными эллиптическими орбитами // Космические исследования. 2004. Т. 42. № 3. С. 260279.

13. Пауэре, Тэпли. Применение канонических преобразований при анализе оптимальных траекторий // Ракетная техника и космонавтика. 1969. № 3. С. 14-21.

14. Ивашкин В.В. Оптимизация космических маневров при ограничении на расстояния до планет. - М: - Наука, 1975. - 392 с.

i Надоели баннеры? Вы всегда можете отключить рекламу.