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

УДК 629.785

ОПТИМИЗАЦИЯ БАЛЛИСТИЧЕСКИХ СХЕМ ПЕРЕЛЕТОВ МЕЖДУ НЕКОМПЛАНАРНЫМИ ОРБИТАМИ С ПОМОЩЬЮ КОМБИНАЦИИ ДВИГАТЕЛЕЙ БОЛЬШОЙ И МАЛОЙ ТЯГИ

© 2010 К.В. Петрухина, В.В. Салмин Самарский государственный аэрокосмический университет Поступила в редакцию 29.03.2010

Исследуется задача выбора оптимальных баллистических схем перелета КА на геостационарную орбиту с использованием двигателей большой и малой тяги. Предлагается схема решения многокритериальной задачи оптимизации для комбинированных перелетов на геостационарную орбиту: в качестве главных критериев принимаются масса полезной нагрузки на целевой орбите и продолжительность перелета, а остальные критерии (общее время пребывания КА в радиационных поясах Земли и длительность теневых участков) переводятся в разряд ограничений. На первом этапе решается двухкритериальная задача оптимизации посредством выделения множества Парето-оптимальных решений, на втором этапе проводится проверка выполнения ограничений. На третьем этапе производится синтез проектно-баллистических характеристик перелета в целом, с учетом ограничений. Метод решения динамической части задачи состоит в использовании принципа расширения допустимых состояний и управлений при переходе от задачи оптимизации в строгой постановке (в соответствии с принципом максимума Л.С. Понтрягина) к задаче локальной оптимизации. Приводятся результаты численного моделирования комбинированных схем перелетов, оптимизации их проектно-баллисти-ческих характеристик, расчета проектных параметров межорбитального транспортного аппарата. Ключевые слова: комбинированная двигательная установка, электрореактивный двигатель, химический разгонный блок, многокритериальная задача оптимизации, множество Парето, принцип расширения, принцип максимума Л.С. Понтрягина, локально-оптимальное управление

В задачах оптимизации космических маневров основными критериями оптимальности являются масса полезной нагрузки и продолжительность перелета. При выборе баллистических схем перелетов приходится искать компромисс между этими критериями. Электрореактивные двигатели малой тяги (ЭРД), использующие принцип ускорения заряженных частиц рабочего тела в электростатических или электромагнитных полях, обладают высокой скоростью истечения рабочего тела (15-30 км/с), что обеспечивает существенно меньший расход рабочего тела по сравнению с двигателем большой тяговоору-женности на химическом топливе. Однако для космических систем малой тяги реактивное ускорение, создаваемое двигателем, составляет 0.110 мм /с2, поэтому маневры с малой тягой в "сильных" гравитационных полях достаточно продолжительны и занимают от нескольких недель до нескольких месяцев [1-4].

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

Петрухина Ксения Вячеславовна, аспирант кафедры летательных аппаратов. E-mail: sputnik@ssau.ru. Салмин Вадим Викторович, доктор технических наук, профессор, заведующий кафедрой летательных аппаратов. E-mail: sputnik@ssau.ru

сти параметров маневра, где оба двигателя по отдельности доставляют примерно одинаковую полезную нагрузку за одинаковое время, применение комбинированных схем может быть более выгодным. В настоящее время близок к завершению проект Российского космического агентства "Фобос-грунт" [7], в котором предполагается достижение спутника Марса - Фобоса с применением комбинированной двигательной установки: химического разгонного блока "Фрегат" и электрореактивного энергодвигательного модуля с солнечной энергоустановкой.

Проблемы оптимизации комбинированных схем перелетов рассматриваются в монографии [3]. Строгая математическая постановка и формулировка условий оптимального сочетания двигателей дана в работе [5]. В работе [6] и монографии [8] при определенной идеализации модели перелета было показано, что оптимальная структура баллистической схемы маневра предполагает использование двигателя большой тяги лишь на начальном участке, а весь последующий маневр осуществляется двигателем малой тяги.

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

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

1. ФОРМУЛИРОВКА МНОГОКРИТЕРИАЛЬНОЙ ЗАДАЧИ ОПТИМИЗАЦИИ

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

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

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

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

Как правило, во многих задачах оптимизации множество допустимых состояний и управлений Б задается посредством некоторых условий, выделяющих его из более широкого множества Е. Задачу оптимизации в этом случае можно решать на указанном более широком множестве Е при некоторых дополнительных условиях. Идея принципа расширения [9, 10] состоит в том, что функционал доопределяется на более широком множестве Е так, что наименьшее значение он принимает в Б. Практически это означает, что путем отбрасывания некоторых связей, ограничений и граничных условий упрощается математическая постановка и вычислительная схема решения задачи оптимизации. Полученное решение "нулевого" приближения проверяется на принадлежность множеству Б. Если условие принадлежности выполняется, то задача решена, а отброшенные связи и ограничения носили надуманный характер, необоснованно усложняя постановку задачи. Если условие принадлежности не выполняется, то полученное решение служит предельной оценкой функционала на расширенном множестве допустимых траекторий и управлений. Такая оценка может оказаться полезной и содержательной на ранних стадиях исследования. Получение подобных ре-

Рис. 1. Комбинированная схема выведения КА на геостационарную орбиту (ГСО)

шении связано с последовательной редукциеи математической модели задачи оптимизации, и соответствующий метод носит название метода последовательных расширений [11]. Затем реализуется обратный процесс: полученное решение уточняется в ходе итерационной процедуры, множество допустимых решений сужается до тех пор, пока не будут выполнены все условия и ограничения, описанные в первоначальной постановке задачи. По существу, здесь реализуется схема последовательных приближений, а полученные предварительно оценки позволяют охарактеризовать степень близости решения к абсолютно оптимальному на любой итерации.

В общем виде задача отыскания оптимальных схем выведения на ГСО является многокритериальной. В качестве основных критериев принимаются:

1. Масса полезной нагрузки, выводимой на целевую орбиту: M ПН ^ max .

2. Продолжительность перелёта: Т ^min.

3. Продолжительность пребывания КА в радиационных поясах Земли: Трад ^ min .

4. Продолжительность пребывания КА в тени

Земли: Тшеш ^ min .

Поставленная многокритериальная задача не имеет строгого решения. Предлагается следующий метод ее поэтапного решения. Проведем декомпозицию задачи на ряд частных задач. Введем критерии времени пребывания КА в тени и в радиационных поясах Земли в ограничения:

доп

T < T

1 рад рад

T < T

доп тени '

2. ПРОЕКТНАЯ МОДЕЛЬ КА С КОМБИНИРОВАННОЙ ДВИГАТЕЛЬНОЙ УСТАНОВКОЙ

Представим стартовую массу КА как сумму масс отдельных систем:

M 0 = M fE + Ы%Б + M ЭЭРДУ + M ЭДРДУ + + МДрДУ + МДРДУ + MДРДУ + M ПН

. (1)

J< доп т<доп „

— рад и 1тени - допустимые продолжительности пребывания КА в радиационных поясах и в тени Земли соответственно.

На первом этапе решается двухкритериаль-ная задача оптимизации:

[MПН ^ max,

T ^ min.

Практически реализуется схема максимизации полезной нагрузки при фиксированной (варьируемой) продолжительности перелёта, т.е. используется классический метод рабочих характеристик [12].

На втором этапе проводится проверка выполнения ограничений. Для этого в ходе моделирования движения КА рассчитываются время пребывания КА с ЭРДУ в тени Земли и время пребывания КА в радиационных поясах. Из множества полученных решений отбрасываются решения, не удовлетворяющие ограничениям.

Третий этап - этап синтеза проектно-баллис-тических характеристик перелета в целом, на основе решения частных задач с учетом ограничений. где 2 =

где М0 - начальная масса КА; М% - масса конструкции ХРБ (сухая масса); МТрБ - масса топлива ХРБ; МЭ,РДУ - масса источника и преобразователя энергии ЭРДУ; М:ДрРДУ - масса

двигательной установки ЭРДУ; Мд1ррду - масса

системы подачи и хранения рабочего тела ЭРДУ;

МДРДУ - масса рабочего тела ЭРДУ; МДРДУ -

масса корпуса КА, прочих элементов и систем КА с ЭРДУ; МПН - масса полезной нагрузки.

Масса конструкции ХРБ считается постоянной и равной сухой массе разгонного блока.

Масса топлива ХРБ складывается из массы топлива, предназначенного для выполнения маневра перехода на промежуточную орбиту, и массы топлива для увода отработавшего разгонного блока на орбиту с низким перигеем:

мХРБ = МХРБ + МХРБ

1У1ТЬ 1У1Т Г™Т увода ■

Массы отдельных компонент КА с ЭРДУ зависят от проектных параметров. Обычно применяются следующие зависимости [3, 13]:

М ЭЭРДУ = аЭУ N,

M^y =Гду (P + kPynp ) ,

m ЭПХУ = ГСПХ M ЭРДУ,

MКЭРДУ =YK P + Y"K N,

P

M ЭРДУ =_ T

1У1Т JM ,

c

(2)

(3)

(4)

(5)

(6)

где Р - тяга маршевых двигателей; РУПР - тяга управляющих двигателей; N - мощность энергоустановки; ТМ - моторное время (полное время работы ЭРДУ); аЭУ, Уду, уСПХ, у' К, у"К - соответствующие удельные массовые характеристики.

Мощность энергоустановки зависит от тяги двигателей и скорости истечения рабочего тела

Рс 1 + ^

(7)

N =

2 пПп

P

характеризует относительный

расход массы управляющих двигателеи, пТ -тяговый коэффициент полезного деИствия, пПЭ - КПД преобразователя энергии.

Если из уравнения баланса масс на начальной орбите (1) выразить массу полезной нагрузки и поделить ее на стартовую массу аппарата, получим следующее выражение для главного критерия оптимальности - относительной массы полезной нагрузки:

U = 1-

МГ - МТ - MZda *,-М+х)

M

Мо Мо МПтПпэ

уду{Р+kPyJ PTM(1 + Yçj

Mn

M0c

(8)

p y"k pcp+x

Здесь M0, MХ = fixe.

Моторное время перелета с двигателем малой тяги при заданных граничных условиях зависит от величины тяги ЭРД (или начального реактивного ускорения) и от скорости истечения и рассчитывается по формуле [8]:

T =■

1м

an

1 - expf - —

(9)

где ¥хк - конечная характеристическая скорость перелета.

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

чем Тм ).

Оптимальная скорость истечения для ЭРД в первом приближении определяется из (8):

opt

V

2ТмПэ

а

(10)

ЭУ

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

1) динамическую - нахождение оптимальных программ управления и траекторий;

2) параметрическую - нахождение оптимальных проектных и баллистических параметров.

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

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

3. ПОСТАНОВКА ДИНАМИЧЕСКОЙ ЗАДАЧИ ОПТИМИЗАЦИИ

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

В качестве математической модели движения КА с ЭРДУ возьмем систему дифференциальных уравнений в оскулирующих элементах:

dA 2p

dt p-e) de

P -[esin^-S+(1+ecos)))], U

dt \

di dt \

dm = 1 dt e\

. n „ ecos2 3+2cos)+e ^ sin)S +---f

1+ecos)

cosu

U 1+ecos)

-■W,

p U

esin u ■ ctgi 1+ecos)

- cos^ s + sin3(2+ecos)) f -

1+ecos)

■W

(11)

dD dt

sinu

U sini(1+ecos)) dU = #■[(+ecos))2 -

■W,

dt

(1+ecos))u

ctgi sinu ■ W

где р = А (1 - в2) - фокальный параметр;

3 = и -а - истинная аномалия; в - эксцентриситет;® - угловое расстояние перицентра от узла; О - долгота восходящего узла; / - наклонение орбиты; т - время прохождения через перицентр; ^ - время; 3 - истинная аномалия; и -аргумент широты; £, Т, Ж - проекции реактив-

c

2

ного ускорения на направление радиуса-вектора, на перпендикулярное к нему в плоскости орбиты и на перпендикулярное к плоскости орбиты; j = fM - произведение гравитационной константы на массу притягивающего центра.

Как известно, система уравнений (11) имеет особенности при e = 0 и i = 0. На практике наиболее распространенным приемом является переход к равноденственным элементам. В данном же случае при моделировании движения КА с ЭРДУ система (11) является более предпочтительной, поэтому перед интегрированием будем задавать фиксированные конечные значения эксцентриситета и наклонения, отличные от нуля и соответствующие требуемой точности решения задачи.

Введем две правые системы координат: орбитальную ( Onrb ) и связанную с КА ( OXYZ ). Вектор тяги p направлен вдоль оси OX.

Запишем выражения для компонент реактивного ускорения в орбитальной системе координат:

T = Sa cos X cos у, S = Sa sin Xcosy, (12)

W = Sa sin у.

Здесь a - модуль полного реактивного уско-

рения ( a =

a

1 - a0t / c

), S - функция включения-

ориентации вектора тяги в плоскости орбиты ( X е [о°;180°]); у- угол ориентации вектора тяги в плоскости местного горизонта ( у е [- 90°; 90° ]) (рис. 2). Очевидно, что минимальное время перелета достигается при постоянной работе двигателя (5 = 1) с переключениями компонент реактивного ускорения.

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

Запишем граничные условия:

t = to t = tK

A(to )= Ao A(tK )= A

e(t0 )= e0 ^ e(tK )= eK

i(t0 ) = i0 i(tK ) = iK

(13)

выключения двигателей (S = {0, 1} ); X- угол

Расширим множество допустимых режимов. Условия на переменные а>, О, u, не наклады-

da dО du

ваются, поэтому уравнения-, -, — мо-

dt dt dt

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

Мгнов*?нная

плоскость орбиты

Мгновенная плоскость местного горизонта

Рис. 2. К определению углов ориентации вектора тяги КА с ЭРДУ

тимальности динамическом задачи - моторного времени перелета (совпадающего с общей продолжительностью при 8 = 1).

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

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

AT = (1 + e cos¿)A +

(1 - e )

е cos25 + 2cos5 + e1T( , +-

1 + e cos5

. cos u .T(

aw =-^

1 + e cos 3

С учетом (12) выражение для гамильтониана принимает вид:

H =

I = AxKaAAxK ^ min ,

(14)

V

[AS cosicosa

+

где AxK = [AA, Ае, Ai]

Здесь

AA = A(t)-Ак,

Ae = e(t)- ек,

Ai = i(t)- iK ■

+ AT sin X cos y + AW sin y]

В соответствии с принципом максимума, необходимое условие оптимальности имеет вид:

AS

a

dH дХ

+ AT cosXcos^] = 0

V

[- AS sin Xcos^ +

"a11 0 0 "

a = [a] ]= 0 a22 0 dH

0 0 a33 _ ' ду \

Bal

V

[- AS cos X sin у -

E

a=1 ,

где aA = a11, ae = a

22

a =a

33

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

Решение вариационной задачи будем проводить в соответствии с принципом максимума Л.С.Понтрягина [14].

Будем рассматривать случай, когда ЭРДУ работает без выключений (8 = 1), поскольку в этом случае обеспечивается минимум общей продолжительности перелета.

Для этого введем сопряженную вектор-функцию X, составим гамильтониан и найдем его максимум по управляющим функциям X и у:

х = (, X, х )г

Последнее выражение можно переписать в следующем виде:

- AT sin X sin у + AW cos у] = 0

Откуда получаем стационарные точки:

■At

sin X = ±

(15)

JA

T + AS

cos у = ±

sin у = ±

-\JAT + AS

(16)

■J AT + AS + AW

cos у = ±

,¡A2T + A

AT + AS + AW

H=

'Aden

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

V

[ass + att + aww ],

A =

X ii X i'

где

AS = 2 ( +)e) e smS-x¥A + smS-x¥e

(1 - e)

где

x =d2H = lA (1 - e2) X11 = dX2 V V - AT sin Xcos у]

a[- AS cos Xcos у -

S

W

д2 H

Хп Х21

дЛ ду \ - AT cos Л sin у]

A(l - е2)

—*-'-а

[.AS sin Л

cosy

д2 H

Х22=

Sr

[- AS cosЛcosy-V (17)

ду2 \

- AT

sin Л cos у - AW sin у]

Подставляя в (17) значения для углов ориентации вектора тяги (14), найдем среди них такое сочетание, которое обеспечивает выполнение условий А = x11 x22 - x12 x21 > 0 и x11 < 0 .

Таким образом, максимум гамильтониана достигается при

sin Л=

cos Л=

sin у =

*JaT + A¡

As

4 AT + A,

A

cos у =

4 AT + As + A

4aT+А

(18)

■J AT + As + A¡

При получении выражений (18) считалось, что двигатель работает без выключений ( § = 1), т.к. этот случай обеспечивает минимальную продолжительность перелета.

Сопряженная система уравнений имеет вид:

• дН

Ч A =--Г7 = fA А, вЛ u,k,w, a, ЧА , Те, Ч, ),

дА

• дН

Чв =-— = fA (А, вЛ иЛ¥, a, Та , Че, Ч ),(19)

дв

Ч, = - — = 0, Ч = const дг ' .

Условия трансверсальности в задаче Майе-ра принимают вид:

дА

дв

* a (t )=-дп,

¥ е (T )=--

(20)

* ' (T )=-f дг

4. МЕТОД ЛОКАЛЬНОЙ ОПТИМИЗАЦИИ

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

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

(I

— в каждый момент времени. Если подынтег-&

ральное выражение не меняет знака и представляет собой монотонную функцию, описанная постановка эквивалентна исходной [15].

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

В общем случае синтез локально-оптимальных управлений не гарантирует абсолютного оптимума в исходной постановке задачи, однако, существует класс задач, содержащих малый параметр (в частности, реактивное ускорение, создаваемое ЭРД), в которых локально-оптимальные управления тем ближе к оптимальным, чем меньше малый параметр, т.е. чем слабее корректирующее управление. Этот результат был доказан и опубликован в работе [4], где были найдены локально-оптимальные управления для широкого класса задач механики полёта с малой тягой.

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

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

T

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

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

Заменим этот функционал локальным критерием, обеспечивающим максимальную скорость изменения I:

dl_

dt

= 2а

. di

v Ao Ao J

1 dA de

■ + 2а e — +

Ao dt

dt

(21)

+ 2a 3i--> max

3 dt

Введем обозначения:

f

A* — AS =а1

A A,

Л

v Ao Ao J

-2A(l±el .e Slni+a2 e -АЛ Ao (1-e) 2

Л*

A =а

T 1

AA

Л

v Ao Ao J

2 A(1 + e) A (1 - e)

(l + e cos5)+

+ a2e •

ecos 5 + 2cos5 + e

A*

A = ai •■

1 + e cos5

cos u

1 + ecos 5 '

Тогда функционал запишем в виде:

di „

— = 2 dt \

А* • S + A* • T + A* • W

(22)

sin Л =

V(t))+(* )2

A*

cosy =

V(t))+A))

A

sin y =

cosy =

v(A; ))+(* )2+(aw )) V(a; )+(A; )2

VA )2+(As )2+(, )2

(23)

Результатом поиска максимума выражения

(18) по двум переменным — = — ((f), y(t)) — max

di di

являются аналитические выражения для углов ориентации вектора тяги Л и У, где y - угол отклонения тяги от мгновенной плоскости орбиты, Л - угол между проекцией вектора тяги на плоскость орбиты и трансверсалью.

Полученный закон управления у(() х(() имеет достаточно простую структуру и позволяет провести расчет динамического манёвра без процедуры решения краевой задачи.

Как следует из (18), от значений весовых коэффициентов аА, ае, ai зависит скорость изменения большой полуоси, эксцентриситета и наклонения орбиты.

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

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

В первом приближении можно принять

аА =ае = а .

В результате применения данного алгоритма

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

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

Путем подбора весовых коэффициентов можно оптимизировать схему перелета по времени.

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

5. МЕТОД РАСЧЕТА ПРОДОЛЖИТЕЛЬНОСТИ ПРЕБЫВАНИЯ КА С ЭРДУ В ТЕНИ ЗЕМЛИ

При проектировании космического аппарата с ЭРДУ важным условием является определение условий его освещенности, поскольку от этих условий зависит функционирование КА, использующего в качестве источника тока бортовые солнечные преобразователи. Эффективность же последних существенным образом зависит от взаимного расположения КА, Солнца и Земли, и, в первую очередь, от расположения КА на освещенных и теневых участках орбиты. Кроме того, тепловые

к

S

T

д

Рис. 3. Зависимость большой полуоси (а), эксцентриситета (б), наклонения (в) орбиты, весовых коэффициентов (г) от времени перелета и траектория пространственного движения КА с ЭРДУ (д) для параметров начальной ( А0 = 14000 км, е0 = 0.3, /051.5°) и конечной ( Ак = 42163 км, ек < 0.001, гК < 0.01° )орбит

режимы конструкций и элементов КА, а также функционирование датчиков системы ориентации также зачастую решающим образом определяются астрономо-баллистическими условиями и, в первую очередь, условиями освещенности [16].

В монографии [8] описана методика расчета времени пребывания КА с ЭРДУ в тени Земли

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

Рис. 4. К определению тени планеты

Рис. 5. К положению КА относительно тени Земли

жуточной слабоэллиптической орбите ( е < 0.05 ). Однако, как следует из работы [17], наиболее предпочтительными для комбинированных схем являются орбиты с большим эксцентриситетом, что требует иных методов расчета продолжительности теневых участков.

Решение задачи определения моментов времени входа КА в тень Земли и выхода из нее предполагает известными на рассматриваемом интервале времени эфемериды Солнца и координаты КА.

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

Запишем так называемое "уравнение тени" [17], определяющие моменты входа-выхода КА в тень Земли (см. рис.5)

ФТ = a32 (1 + е cosS)2 + + p2(AcosS + BsinS)2 -p

2

(24)

где a3 - радиус Земли; rC - радиус-вектор Солнца в геоцентрической системе координат.

Условием входа в область тени или выхода из нее является равенство ФТ = 0, а физический смысл ( cos р < 0 ) имеют те решения уравнения (21), для которых выполняются условия

A ■ cosS + B ■ sin S< 0 .

При выполнении условий

A ■ cosS + B ■ sin S = 0. (25)

КА будет освещен прямыми лучами Солнца.

При анализе условий освещенности необходимо учитывать следующее правило: если КА входит в тень, то ФТ меняет знак с минуса на плюс, если КА выходит из тени - ФТ меняет знак с плюса на минус [17].

6. ОПТИМИЗАЦИЯ ДАТЫ СТАРТА

С УЧЕТОМ ВРЕМЕНИ ПРЕБЫВАНИЯ КА С ЭРДУ В ТЕНИ ЗЕМЛИ

В зависимости от даты старта Dcmapma, определяющей эфемериды Солнца, и начальной долготы восходящего узла Q0, задающей начальную ориентацию плоскости орбиты относительно Солнца, траектория космического аппарата с солнечной ЭРДУ будет характеризоваться различным временем затенения. Варьируя параметры Q0, Dcmapma, получаем различные значения времени пребывания КА в тени.

Поставим задачу построения изолиний времен пребывания КА, совершающего перелет по многовитковой пространственной траектории перехода на ГСО, в тени Земли. Для этого выполним серию расчетов для различных значений параметров О0, Остарта с шагом по долготе восходящего узла и по дате старта. Обработка массива расчетных значений проводилась с помощью пакета MathGL.

При этом равными считались значения времени пребывания КА в тени Земли, отстоящие от среднего значения на 0,1 суток.

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

На рис. 6 представлены результаты расчета продолжительностей теневых участков при варьировании значений даты старта и долготы восходящего узла при фиксированных значениях большой полуоси А = 14000 км и эксцентриситета е = 0.3 промежуточной орбиты.

Из рис. 6 видно, что оптимальные ( Ттени = 0 ) и неоптимальные ( Ттени = max ) даты старта повторяются с периодичностью 6 месяцев. Кроме того, существуют достаточно широкие окна старта, при которых Ттени = 0. Ранее этот факт был установлен для движения по околокруговым траекториям [8].

При фиксированной дате старта можно добиться уменьшения времени пребывания КА в тени за счет оптимального выбора начальной ориентации плоскости орбиты (угла Q0 ). Описанная приближённая методика позволяет существенно сузить границы области поиска оптимальных дат старта для космических аппаратов с солнечными электрореактивными двигателями.

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

Рис. 6. Области равной продолжительности пребывания КА с ЭРДУ в тени Земли: а - поверхность в пространстве; б - изолинии ( Ттени = fixe )

Первоначальное расширение множества допустимых решений оказалось достаточно эффективным, т.к. диапазон параметров, при которых Ттени = 0, оказался весьма широким.

Отметим, что для серии расчетов, выполненных с различными значениями большой полуоси переходного эллипса, время пребывания КА в тени составляет от 0 до 7 суток. Это существенно меньше, чем для траекторий перелетов на ГСО, соответствующих движению по квазикруговой орбите только с двигателем малой тяги -от 0 до 20 суток [8].

7. МОДЕЛЬ РАСЧЕТА ВРЕМЕНИ ПРЕБЫВАНИЯ КА С ЭРДУ В ТЕНИ ЗЕМЛИ

Движение космического аппарата в областях повышенной радиационной активности сопровождается деградацией солнечных батарей, а также приводит к увеличению требований к радиационной защите приборных отсеков КА, что значительно увеличивает их массу. Для расчета полученной дозы радиации, при геоцентрическом движении аппарата, используется модель распределения интенсивности радиации [18]. Радиационные пояса в этой модели считаются симметричными и на рис. 7 показано их сечение плоскостью, перпендикулярной экватору и проходящей через ось Земли (Ь - расстояние, измеренное в радиусах Земли). Эти линии уровня приближенно аппроксимируются уравнениями в полярной системе координат.

Внутренний радиационный пояс имеет максимальную интенсивность радиационного потока 1,9-103 Гр/сут , внутренний -2,8-102 Гр /сут.

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

В табл. 1 приведена сравнительная характеристика перелетов, отражающая время пребыва-

ния КА в радиационных поясах Земли в зависимости от параметров промежуточной орбиты -большой полуоси и эксцентриситета.

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

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

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

8. РЕШЕНИЕ МНОГОКРИТЕРИАЛЬНОЙ ЗАДАЧИ СИНТЕЗА ПРОЕКТНО-

БАЛЛИСТИЧЕСКИХ ХАРАКТЕРИСТИК

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

На первом этапе на основе перебора параметров баллистической схемы b (большой полуоси и эксцентриситета промежуточной орбиты) с заданным шагом производится моделирование каждого расчетного варианта. В плоскости параметров "время перелета - масса полезной нагрузки" строится множество решений п1 . Далее применяется метод рабочих характеристик (MПН = fixe, ГЕ = var ) и исходное множество решений п сужается до множества Парето п2 сп, которое является верхней границей множества п1.

На втором этапе производится проверка ог-

«-» т7 т'доп т ^ грдоп г,

раничений 1 рад < 1рад , lmgHU < 1тени . Шр^ТЬ^ не удовлетворяющие ограничениям, исключаются из множества п2. Таким образом, результатом

Таблица 1. Сравнительная характеристика пребывания КА в радиационных поясах Земли при варьировании параметров промежуточной орбиты

Рис. 7. Модель радиационных поясов Земли

№ Aif, е пр 1 ЭРДУ трРАД

км сут. сут.

1 17000 0,059 81,773 21,28

2 19000 0,158 71,578 14,47

3 21000 0,238 64,780 10,10

4 23000 0,304 59,955 7,36

5 25000 0,360 56,315 5,63

6 27000 0,407 53,570 4,48

7 29000 0,448 52,666 3,70

8 31000 0,484 49,613 3,15

решения задачи оптимизации является некоторое множество "неулучшаемых" решений п3 ^ п2, в котором можно выделить две области:

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

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

На рис. 8 представлены множества "неулучша-емых" решений для случаев последовательного и совместного управления элементами орбиты.

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

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

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

г . @

* г 17 £ Р* 5 6 1 ч И 13 г- ь 13 • 1 I ® •

& ¥ 3 ,ЗЕ

за 1за 5

1 43 а? •-1

<

I , Г.уТ.

б)

Рис. 8. Построение множества Парето для оптимальной даты старта (21.05.2011г): а - схема последовательного управления (А, е) ^ г; б - схема совместного управления (А, е, г)

дываются ограничения на угловое движение КА с ЭРДУ и возникает необходимость разделения каналов управления.

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

1) разворот КА с помощью ЭРД (в этом случае каждый ЭРД помещается в кардановый подвес);

2) разворот КА с помощью двигателей малой тяги;

3) разворот КА с помощью маховиков;

4) без разворотов (двигатели располагаются неподвижно на КА с противоположных сторон, их переключение производится в соответствии с законом управления).

Приведенный массив расчетов, отраженный на графиках (см. рис. 8), показывает преимущества совместного управления элементами орбиты. Для случая "а" масса ПН при стартовой массе 6900 кг составляет порядка 1350 кг, а время перелета изменяется от 65 до 85 суток. Дальнейшее увеличение продолжительности перелета не дает выигрыша в полезной нагрузке (точки 8-3).

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

Первый участок продолжительностью от 50 до 60 суток: здесь масса полезной нагрузки не превышает 1300 кг, т.к. основные энергетические затраты берет на себя ХРБ при реализации переходного эллипса (короткие перелеты).

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

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

риситетом и наклонением орбиты при постоянных весовых коэффициентах.

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

На рис. 9 показаны зависимости углов ориентации вектора тяги КА с ЭРДУ. Видно, что управление углом ориентации вектора тяги в плоскости орбиты (случай "а") на начальном этапе (от 0 до 18 суток) происходит по закону, близкому к трансверсальному, после чего (от 18 до 75 суток) амплитуда колебаний значительно возрастает. В то же время изменение угла ориентации в местной горизонтальной плоскости (случай "б") на начальном этапе (от 0 до 30 суток) происходит за счет переключения бинормальной компоненты вектора тяги дважды за виток с постепенным возрастанием амплитуды от 40 до 75 градусов, после чего закон управления медленно приближается к бинормальному.

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

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

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

Рис. 9. Углы ориентации вектора тяги: а - в плоскости орбиты, б - в плоскости местного горизонта

Таблица. 2. Результаты расчета проектно-баллистических параметров перелета на ГСО (точка №38, см. рис. 9 "б")

Исходные данные

Начальная масса КА, кг 6900

Масса рабочего запаса топлива ХРБ, кг 3164

Конечная (сухая) масса ХРБ, кг 980

Тяга двигателя ХРБ, Н 20000

Удельная тяга двигателя ХРБ, с 328

Тяга одного ЭРД (СПД-140), мН 280

Удельный импульс ЭРД, с 2600

Потребляемая мощность одного ЭРД, кВт 6

Масса одного ЭРД, кг 1.2

Количество рабочих ЭРД 8

Коэффициент резервирования ЭРД 1.5

Режим работы ЭРДУ в тени Земли не выключается

Удельная масса энергоустановки, кг/кВт 10

Удельная масса системы подачи и хранения рабочего тела 0.07

Относительная масса конструкции КА (после отделения МДУ РБ «Фрегат») 0.15

Высота начальной орбиты, км 350

Наклонение начальной орбиты, град. 51.5

Большая полуось промежуточной орбиты, км 16000

Эксцентриситет промежуточной орбиты 0.5

Долгота восходящего узла промежуточной орбиты, град. 60

Дата старта КА с ЭРДУ с промежуточной орбиты 01.01.2011

Баллистические параметры перелета

Большая полуось (ГСО), км 42163.992

Эксцентриситет конечной орбиты 0.0001

Наклонение конечной орбиты, град. 0.01

Время перелета, сут. 73.105

Моторное время, сут. 73.105

Время пребывания КА в радиационных поясах Земли, сут. 23.934

Время пребывания в тени Земли, сут. 1.296

Период обращения, ч. 23.934

Количество витков 121

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

Проектные параметры перелета

Масса КА на промежуточной орбите, кг 2756

Масса энергоустановки, кг 480

Масса ЭРДУ, кг 14.400

Масса СПХ, кг 38.830

Масса рабочего тела ЭРДУ, кг 554.709

Масса конструкции, кг 413.445

Масса полезной нагрузки, кг 1254.919

Полная тяга, мН 2240

Начальное ускорение, мм/с2 0.813

Площадь солнечных батарей, м2 173.913

Авторы статьи выражают благодарность за полезные замечания профессору О.Л. Стариновой

СПИСОК ЛИТЕРАТУРЫ

1. Tsien H.S. Take-Off from Satellite Orbit // J.Am.Rocket Soc. 1953. 23. 233.

2. Ирвинг Д. Полеты с малой тягой в гравитационных полях при переменной скорости истечения // Космическая техника / Под ред. Г. Сейферта. М.: Наука, 1964. С. 286-324

3. Гродзовский ГЛ., Иванов Ю.Н., Токарев В.В. Механи-

ка космического полета. Проблемы оптимизации. М.: Наука. 1975. 702 с.

4. Лебедев В.Н. Расчет движения космического аппарата с малой тягой, - М: ВЦ АН СССР, 1968. - 108 с.

5. Иванов Ю.Н. Оптимальное сочетание двигательных систем / Изв. АН СССР. Механика и машиностроение. 1964. №2. С. 3-14.

6. Салмин В .В. Оптимальное управление комбинированной системой, состоящей из двигателя ограниченной скорости истечения и двигателя ограниченной мощности // Космические исследования. 1970. Т.УШ. №4. С.535-541.

7. Попов Г.А., Обухов В.А., Константинов М. С., Федотов Г.Г., Мурашко В.М. Применение электроракетной дви-

гательной установки в проекте "Фобос-грунт"// Фундаментальные и прикладные проблемы космонавтики. 2001. №4. С. 26-31.Брусов В.С., Салман В.В. Комбинированная двигательная система, универсальная для диапазона маневров // Космические исследования. 1974. Т. XII. Вып. 3. С. 368-373.

8. Салмин В.В. Оптимизация космических перелётов с малой тягой. Проблемы совместного управления траекторным и угловым движением. М.: Машиностроение, 1978. 208 с.

9. Гурман В.И. Принцип расширения в задачах управления. М.: Наука, 1985

10. Кротов В.Ф., Гурман В.И. Методы и задачи оптимального управления. М.: Наука, 1973.-446 с.

11. Салмин В.В., Васильев В.В, Ишков С.А., Романенко В.А., Соколов В. О., Старинова О.Л., Юрин В.В. Приближенные методы расчета оптимальных перелетов космических аппаратов с двигателями малой тяги // Вестник СГАУ. 2007. №1 (11). С. 37-52.

12. Лебедев А.А. Введение в анализ и синтез систем. М.: Изд-во МАИ, 2001. 352 с.

13. Ишков С.А., Салмин В.В. Оптимизация траекторий и параметров межорбитальных транспортных аппаратов с двигателями малой тяги / / Космические исследования. 1989. T.XXVII. Вып. 1. С.42-53.

14. Математическая теория оптимальных процессов / Л.С.Понтрягин,А.Г.Болтянский [и др.]; [под ред. Л.С. Понтрягина]. М.: Наука, 1976. 392 с.

15. Моисеев Н.Н. Элементы теории оптимальных систем, М.: Наука, 1975, 528с.

16. Попович П.Р., Скребушевский Б.С. Баллистическое проектирование космических систем. М.: Машиностроение, 1987. 240 с.

17. Транспортные модули на базе комбинации современных химических двигателей и электроракетных двигательных установок для транспортных перевозок "орбита - орбита". Отчет о НИР. НИИПМЭ, 1998.

18. Вернов С.Н., Вакулов П.В., Горчаков Е.В. Радиационные пояса Земли и космические лучи. М.: Просвещение, 1970. 128 с.

OPTIMIZATION OF FLIGHTS BALLISTIC SCHEMES BETWEEN NON-COPLANAR ORBITS BY MEANS OF THE COMBINATION OF HIGH AND LOW THRUST ENGINES

© 2010 K.V. Petrukhina, V.V. Salmin Samara State Aerospace University

The problem of a choice of spacecraft flight optimum ballistic schemes into a geostationary orbit with use of high and low thrust engines is investigated. The scheme of optimization multicriterial problem decision for the combined flights into a geostationary orbit is offered: as the main criteria are accepted weight of payload in a target orbit and duration of flight, and other criteria (the general time of spacecraft stay in the Earth radiating belts and duration of shadow sites) are translated in the category of restrictions. At the first stage the two-criterial problem of optimization by means of generation of Pareto-optimum decisions, at the second stage check of performance of restrictions is carried out. At the third stage synthesis of design-ballistic characteristics of flight as a whole, taking into account restrictions is carried out. The method of a dynamic part problem decision consists in use of an expansion principle of admissible conditions and controls at transition from an optimization problem in strict statement (according to a Pontrjagin's maximum principle) to a problem of local optimization. Results of numerical modeling of the flights combined schemes, optimization of their design-ballistic characteristics, calculation of design parameters of the interorbital transport spacecraft are resulted. Key words: combined engine, electrojet engine, chemical booster, multicriterial problem of optimization, set Pareto, expansion principle, Pontrjagin's maximum principle, locally-optimum control

Ksenya Petrukhina, Graduate Student. E-mail: sputnik@ssau.ru. Vadim Salmin, Doctor of Technics, Professor, Head at the Aircraft Department. E-mail: sputnik@ssau.ru.