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

УДК 629.78

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

© 2011 А. А. Храмов

Самарский государственный аэрокосмический университет имени академика С.П. Королёва (национальный исследовательский университет)

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

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

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

Кроме того, в настоящее время разрабатываются ДУ, которые по уровню тяги занимают промежуточное положение между ДМТ и ЖРД. К примеру, тяга разрабатываемого в настоящее время магнитоплазменного двигателя с переменным удельным импульсом УАБГМК [1] составляет величину порядка нескольких Ньютон. Его лётные испытания планируются на международной космической станции (МКС) в 2013 году. Вследствие большой потребляемой мощности (порядка сотен киловатт) ДУ будет включаться эпизодически на несколько минут. При этом будет использоваться энергия, накопленная в

специальных аккумуляторах, которые, в свою очередь, во время пауз в работе плазменного ускорителя будут подзаряжаться от солнечных батарей МКС. Предполагается, что использование УАБГМК понизит массу топлива для поддержания орбиты МКС в год с 7500 кг до 350 кг.

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

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

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

В настоящее время при решении задач межорбитальных манёвров используются три возможные постановки: с двигателем большой, малой и ограниченной тяги. Параметром, определяющим постановку задачи, является отношение ускорения от тяги ДУ к гравитационному ускорению [3]. В первом случае характерной особенностью является малая продолжительность активных участков, что оправдывает применение импульсной аппроксимации и сводит задачу оптимизации массы полезной нагрузки к задаче минимизации характеристической скорости манёвра. Механика полёта с малой тягой отличается значительной протяжённостью активных участков, многовитковыми траекториями перелёта. Основным методом решения является усреднение модели движения по виткам вследствие малости реактивного ускорения и оптимизацию манёвра на значительных интервалах времени. Оптимизация межорбиталь-ных манёвров с двигателем ограниченной тяги в общей постановке без упрощающих допущений сводится к решению нелинейных двухточечных или многоточечных краевых задач для систем обыкновенных дифференциальных уравнений, для которых не существует надёжных методов решения, гарантирующих сходимость задачи.

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

поле Земли. Под коррекцией понимается манёвр КА с незначительным изменением большой полуоси орбиты по сравнению с её начальным значением. Предполагается, что эксцентриситет начальной и конечной, а также переходных орбит лежит в диапазоне: 0 < e < 0,1. Управление движением центра масс КА осуществляется с помощью последовательности чередующихся активных и пассивных участков полёта. На активных участках реализуется включение ДУ, создающей постоянную по величине тягу. Суммарная продолжительность пассивных участков на витке считается заданной и постоянной в течение манёвра. Необходимо определить оптимальную программу управления ДУ, энергетические и временные затраты на манёвр. В качестве критерия оптимальности рассматривается время межорби-тального перехода.

Для описания орбитального движения КА в плоскости орбиты в центральном гравитационном поле используется математическая модель в оскулирующих элементах [4]:

dA „ — = 2. dt 1

-wT

R

m

V1 - q2 - k2

dq = dt

dk = dt

du

dt

dVx

I—íwS smu + wT

+ cosu

q + cos u

R ,,

V yj

/ k + sin u Y|

+ sinu >

(1)

V

R

4

w¡ + w'T

dt Здесь

2 2

p = A (1- q - k ) - фокальный параметр орбиты,

A - большая полуось, q и k - компоненты вектора Лапласа, u - аргумент широты,

Vx - характеристическая скорость перелёта,

wT = w cos a и wS = w sin a - компоненты вектора ускорения тяги ДУ вдоль транс-версального и радиального направлений, соответственно,

3

w =P/m - ускорение от силы тяги,

P - тяга ДУ, m - текущая масса КА, a - угол между вектором тяги ДУ и трансверсалью,

/т - гравитационный параметр,

R=1 +qcosu+ksinu.

Компоненты вектора Лапласа связаны с эксцентриситетом e и аргументом перигея w орбиты следующими соотношениями:

q = e cos w, k = e sin w .

С целью упрощения математической модели с учётом специфики решаемой задачи вводятся следующие ограничения на управление и параметры орбиты:

- рассматриваются ДУ с нерегулируемой тягой и скоростью истечения;

- вектор ускорения ориентируется вдоль трансверсали wS = 0 ;

- вследствие малого расхода рабо-

чего тела в процессе коррекции ускорение от тяги считается постоянным:

wT =dw0, где 8е {-1,0,1} - функция включения тяги;

w 0 - ускорение от тяги в начальный момент времени;

- орбиты движения КА в процессе манёвра определяются эксцентриситетами от 0 до 0,1, что позволяет принять компоненты вектора Лапласа q и k в правых частях уравнений движения равными нулю;

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

После перехода к новой независимой переменной, аргументу широты u с учётом введённых ограничений, уравнения движения центра масс КА запишутся в виде

dA 0 A3 „

— = 2w0— 8 ,

du

m

dq _ A2 c

— = 2w0— о cos u,

du

m

(2)

dk „ A2

— = 2w0 — о sin u .

dK

du

A3

— о . m

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

u

0, A = A«, k = k

u = uk, A = Ak

о, q = q0,

*k> л~лк> k = kk’ q = qk-

Метод усреднения предполагает разделение движения КА на периодическую (быстро меняющуюся) и вековую (медленно меняющуюся) составляющие и отдельную их оптимизацию [5, 6]. Таким образом, поиск оптимального управления состоит из двух этапов. Первый этап - определение оптимального управления в пределах витка (локальная оптимизация); второй - построение усреднённой модели движения и определение оптимального управления вековым изменением параметров орбиты.

Для решения задачи локальной оптимизации воспользуемся принципом максимума Понтрягина. В качестве критерия оптимальности принимается минимум затрат характеристической скорости (Vx ® min)

при фиксированной продолжительности манёвра, равной одному витку. Гамильтониан для неавтономной системы (2) и принятого критерия оптимальности записывается в виде

H* = H + yu ® max,

dV dA dq dk

H =----- + — yA +— yq +— yk,

du du du q du

где

yA ,yq yk ,yu - сопряжённые множители. После преобразований получаем A2

{AyA +¥q cosu+yк

H = 2w0 — d{AyA +yq cosu+yk sinu}-m

A3

—О. m

Запишем соотношения для сопряжённых множителей:

du

м

m

3

2 aVa+

+y cosu+yk slnu

+

du

3 Alsl

+ö woJ~ 0 ’ 2 ]m

m

dy і _ _0

du du yq _ const, yk _ const ,

(3)

-2w,

A"

d {

yk cosu —yq sinu

}•

ёи ди т

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

включения и выключения ДУ определяются уравнением

Г (и ) = 2 ААуА + уа 0со$и + ук 08т Ы -1 = 0, (3)

1

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

8орг = (А У А + Уд0 С0§ Ы + У к 0 ^ Ы )| =

= ^п( I(ы))-

Рис. 1. Оптимальное управление на витке

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

двумя одинаковыми пассивными участками (рис. 2). В частном случае один из активных участков может отсутствовать.

Используя полученную структуру управления, перейдём к построению модели векового изменения параметров орбиты. В качестве параметров управления вводятся следующие величины: X - половина ширины

разгонного участка с аргументом широты его центра Ц (рис. 2). Ширина одного пассивного участка а задана и постоянна в течение манёвра. Размеры участков измеряются по аргументу широты.

йХ 1 г йХ ,

= — ------------йи.

9тг і

йи

2г 0 йи

йХ

йи

2г

йХ

йи

йи

2 г +ц-%-а

I

П+Х+а

йХ

йи

йи

йА ~йї ~ 1

А3 w,

а

)];

dq . А»0 . — = 4------------081П

йі

Ґ е (X л

Х+—

2

008

йк . м ^ гк0 .

— = 4)---------------08іп

йі

г е а л

Х+—

2

008

'а л

,2,

V /

/ал

2

оо8(п); (А) 8іп(п);

йі

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

А =А -

А0

т =

-і

ух =

А^г

х

т

где

£0 =-

т

гравитационное ускорение на

начальной орбите,

г - текущее безразмерное время.

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

^=2л/ач(1 - а ^2Х - (ж - а)] •

ж

йт

Рис. 2. Структура управления на витке

Определение структуры управления на витке позволяет провести процедуру усреднения системы уравнений движения (2):

2ж

^

йт

^=4ҐАц:

йт

\ -а'

'і - “'

г

81П

Х +

-а

о8(П)

(5)

і-а

. 81П

Г

Х +

1п(п)

г-а

йГх

йт

где X и X - соответственно исходный и

усредненный вектор состояния.

С учётом введённых параметров управления можно записать

После преобразований и перехода к независимой переменной I получаем усреднённую систему:

В соответствии с методом усреднения [6] точность полученной этим методом модели движения имеет порядок малого параметра £ на интервале движения 1/е В системе (5) в качестве малого параметра выделяется произведение безразмерного ускорения и относительную величину активных участков на витке:

£ =

V

Это позволяет использовать полученную модель без потери точности для расчётов с более высоким уровнем тяги при соответствующем уменьшении продолжительности активных участков на витке. На рис. 3 представлены изолинии равных точностей в плоскости уровня ускорения от тяги и суммарной длительности активных участков для средней высоты орбиты Н0 = 500 км. С помощью этого графика можно оценить максимально допустимую продолжительность активных участков на витке при из-

2

2

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

АЛ = е Л0, Ак = Ад = е.

Оценим ошибки расчёта эксцентриситета и аргумента перигея орбиты, связанные с компонентами вектора Лапласа следующим образом:

е = у] д2 + к2 , о = аг^ — .

д

Используя известное соотношение для полного дифференциала функции

Эг Эг

двух переменных: аг = — ах +--------ау, можно

Эх Эу

записать:

Эе , Эе

Ае = —Ак +----------------Аа = е

Эк Эа

= е (бІП С0к + 008 СОк ),

Эю ., Эю .

Аю =----------Ак +--------Аа = е

Эк Эа

Эе Эе

--------1--------

Эк Эа

Эю Эю

--------1------

Эк Эа

=е

(008 Юк - 8ІП Юк )

где ек, о к - соответственно эксцентриситет и аргумент перигея конечной орбиты.

0 30 60 90 120 150 180 210 240 270 300 330 360

X, град

Рис. 3. Оценка ошибок расчёта по усреднённой модели

Задача определения оптимального управления на значительных интервалах времени математически формулируется следующим образом. Требуется определить закон изменения параметров управления X (О, Л (0, переводящих систему (4) из начального состояния Х0 = {А0, Ч0, к0 } в конечное X к = {Ак, ч к, к к } за минимальное время при заданной длительности пассивных участков на витке.

Предположение о малости изменений параметров орбиты за виток эквивалентно ограничению на прираще-

ние характеристической скорости за

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

А и = 2(р - а):

У

т-і =

>3 мот Аср

т ау:

АСР

где

tмот - моторное время на витке,

Лср = Лк +Л - среднее значение боль-

ср 2

шой полуоси орбиты.

Длительность одного пассивного участка будет равна 117

е

к

a = ж -

m av;

4 2wo

(6)

При этом параметр X ограничен диапазоном значений

^ аг;

0 < X < ж - a или о < x <

A 2wo

ср 0

(7)

Для решения задачи воспользуемся принципом максимума Понтрягина. Составим гамильтониан системы (4):

, йА йа йк Н =-' + ЛУа + йУа + * Ук =

: 1 + 2

Awo

m ж

ау a [2Х -(ж - a)]+

+ 2sin

cos

a

2

X

x\yq cosh +Vk sinh]

Уравнения для сопряжённых множителей записываются в виде:

У л [2£ -(р - а)]+

d¥ a dt

Awo

m ж

2 .

+— sin

A

X

[У q

Ґ к a Л

X + —

v 2

a

cos

cosh + ¥k sin

inh ]

2

X

dyq dy k

dt dt У q = У q0 = COnSt, Уk = yk0 = COnSt •

Оптимальные значения параметров управления определятся из необходимых условий максимума гамильтониана. Для положения центра активных участков h (линия переключения) можно записать: dH

-----= -Уq0 sin h + yk0 cos h = 0,

dhopt h opt

arctg

y k0 У q0

(8)

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

dk — = tgh, dq

tgh =

k к - k0

qk- qc

(9)

При этом количество неизвестных гамильтониана понижается на единицу, так как с учётом (8):

yq cosh +Уk sinh = yqo cosh + sin h ----ГЛ

+ yqo---- = ^yqo (1 + tg h j

cosh

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

dH = 4 A wo dX ]¡ m ж

АУа +

+ cos

є аЛ

X+—

, 2

v /

cos

^аЛ

~2 V /

Xopt =arccos

ХЛІУlo(l + tg 2h) АУA

= 0:

cos

ґ a Л

2

v /

-y/У q'o (1 + tg 'h )

При вычислении оптимального значения X необходимо учитывать ограничение (7), которое можно записать в виде

0 < arccos

АУ A

cos

ra Л

2

V у

¡У lo(l + tg 2h)

a.

Разрешая это неравенство относительно a , получаем

А\УA

< cos2

a

2

(11)

VУ lo I1 + tg h )

Если это условие не выполняется, то параметр X лежит на одной из границ допустимого интервала (7) в зависимости от знака

Ул-

1 + sign (ул)

?opt

-(ж - a).

(12)

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

В итоге применение формализма Пон-трягина сводит задачу оптимальной коррекции орбиты к двухточечной краевой задаче для системы четвёртого порядка:

X

йЛ . Л3 wn г_ „ .

— = 2-----------1 [2Х - (р - а)],

м \ т р

йу wn .

-1 = 4 —П81П

Ж

йк

йУЛ

м

81П

X +-

I т р

2

+---81П

Л

а' (а

— С08 —

2 ^ 2

/

\ / л

а а

— С08

2 2

/ V

а 1 1 )]+

/ .. \

а х+— 2 С08

V /

С08(^ )>

8ш(^ ) ,

2

V 7

(13)

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

мости фазовых параметров. Так как решение сопряжённой системы определяется с точностью до произвольного множителя, то можно принять: уЛп = 1.

Задача решалась модифицированным методом Ньютона. В качестве невязок для подбора параметров 1к и щ0 использовались ошибки по большой полуоси Л и компоненту вектора Лапласа ч, соответственно.

На рис. 4 и 5 представлены результаты решения краевой задачи для различного уровня ускорения от тяги при следующих фиксированных граничных условиях и ограничении :

во = 0.014, Оо = Пп;

ек = П.ПП2, ок = 36п;

Л0 = 6978.2 км. Лк =7038.2 км,

А Vе = 0.1 м/с.

Рис. 4. Управление при различных уровнях ускорения от тяги

Л

г а '

х

ш,

град.

t, сут.

Рис. 5. Зависимость параметров орбиты от времени

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

Рассмотрим зависимости параметра Х(0> определяющего структуру

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

АУХ. Структура управления при этом имеет один разгонный участок с центром Ц или один тормозной с центром Ц+р максимально возможной продолжительности. Наклонные участки отображают зависимость управления от времени внутри области допустимого управления (7), структура которого имеет два активных участка - разгонный и тормозной. Видно, что при уровне ускорения от тяги, равном 2 -10"4 м/с2 , управление в течение всего манёвра лежит внутри допустимого диапазона. При уровнях ускорения 4 -10-4 м/с2 и 6 -10-4 м/с2 управление на первом этапе состоит из одного разгонного участка, а на втором появляется тормозной участок. При уровне ускоре-

—3 2

ния 2 -10 м/с имеются два временных интервала с управлением по ограничению: первый с разгонным, второй с тормозным участками, между которыми находится интервал, где управление имеет и разгонный, и тормозной участки.

В пределе, при бесконечно большой тяге, наклонная прямая, соответствующая структуре управления с двумя активными участками, стремится к вертикальной линии, т.е. продолжительность этого этапа стремится к нулю. В итоге получаем решение задачи коррекции, при котором управление будет состоять на первом этапе только из разгонных участков, а на втором -только из тормозных максимально возможной продолжительности. Это соответствует решению задачи в импульсной постановке с разгонным и тормозным импульсами. И наоборот, при уменьшении тяги наклонная прямая X(t) становится всё более пологой и в пределе стремится к горизонтальной прямой, соответствующей X(t) = const; при этом длительности активных участков остаются постоянными.

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

40

35

30

25

15

10

5

витке. В случае ограниченной тяги используется численное решение краевой задачи; структура управления меняется в процессе манёвра. Для большой тяги применяется управление по ограничениям, соответствующее импульсному решению с разгонными и тормозными импульсами.

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

Таблица 1. Оценка методических ошибок

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

1. Ускоритель плазмы осилит дорогу до Марса за месяц [Электронный ресурс]: Портал membrana: люди, идеи, технологии. URL:

http://www.membrana.ru/articles/technic/ 2009/10/07/185700.html (дата обращения 12.03.2010).

2. Усталов, Ю. М. Анализ некоторых схем перевода КА на геостационарную орбиту с использованием солнечной тепловой двигательной установки ограниченной тяги [Текст] / Ю. М. Усталов, А. А. Храмов // Сборник трудов девятого Всероссийского научнотехнического семинара по управлению движением и навигации ЛА. - 1999. -Ч.1. - С. 141-144.

3. Гришин, С. Д. Проектирование космических аппаратов с двигателями малой тяги [Текст] / С. Д. Гришин. - М.:

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

Машиностроение. - 1990. - 224 с.

4. Эльясберг, П. Е. Введение в теорию полета искусственных спутников Земли [Текст] / П. Е. Эльясберг. - М.: Наука. -1965. - 338 с.

5. Салмин, В. В. Методы решения вариационных задач механики космического полета с малой тягой [Текст] / В. В. Салмин, С. А. Ишков, О. Л. Старинова. - Самара: Издательство Самарского научного центра РАН. - 2006. - 164 с.

6. Моисеев, Н. Н. Асимптотические методы нелинейной механики [Текст] / Н. Н. Моисеев. - М.: Наука. - 1981. - 400 с.

7. Храмов, А. А. Расчет маневров коррекции слабоэллиптических и круговых орбит с двигателем малой и конечной тяги [Текст] / А. А. Храмов, С. А. Ишков // Известия Самарского научного центра РАН. -2002. - Т.4. - №1. - С. 144 -152.

OPTIMUM PROGRAMS OF CORRECTION OF QUASIELLIPTICAL AND CIRCULAR ORBITS OF SPACECRAFT WITH THE LIMITED-THRUST ENGINE

© 2011 A. A. Khramov

Samara State Aerospace University named after academician S. P. Korolyov

(National Research University)

The paper deals with the problem of optimal correction of quasielliptical and circular orbits with limited thrust engines in the Earth central gravitational field with the use of spacecraft motion model for Laplace vector components. The control structure on the turn is justified the motion equations are averaged and their accuracy is assessed. The influence of the level of acceleration caused by thrust on the optimal control structure is investigated.

Spacecraft, limited thrust, optimal control structure, averaging method, estimation of accuracy.

Информация об авторах

Храмов Андрей Александрович, инженер-программист кафедры динамики полёта и систем управления. Самарский государственный аэрокосмический университет имени академика С.П. Королёва (национальный исследовательский университет). Область научных интересов: оптимальные законы управления космическими аппаратами. E-mail: hramovaa76@rambler.ru.

Khramov Andrey Alexandrovitch, programming engineer, department of flight dynamics and control systems, Samara State Aerospace University named after academician S. P. Korolyov (National Research University), hramovaa76@rambler.ru. Area of research: optimal laws of spacecraft control.