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

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

Вместе с тем группировки КА ДЗЗ, развёртываемые на орбитах малого наклонения (от ~0 до 20°), позволяют существенно повысить периодичность покрытия экваториальных областей. Даже в случае использования одного КА в составе космической системы периодичность полного покрытия экватора может достигать примерно 100 мин. Перспективы развёртывания КА на орбитах малого наклонения ранее рассматривались, например, в работе [1]. Однако, в силу относительно северного расположения отечественных космодромов, требуются существенные энергетические затраты средств выведения для доставки полезных нагрузок на подобные орбиты.

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

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

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

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

Для проектирования комбинированных схем выведения рассматриваемого класса полезных нагрузок предлагается следующая программа полёта:

• •    групповое выведение КА в составе головного блока ракеты-носителя (РН) на опорную круговую орбиту;

• •    перевод группы КА с помощью маршевой двигательной установки (МДУ)

РБ на некоторую промежуточную эллиптическую орбиту по схеме оптимального трёхимпульского апсидального перелёта;

• •    переход отдельных КА с помощью собственных ЭРДУ на рабочие орбиты.

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

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

При фиксированных параметрах ПО вычисляется масса КА, которая доставляется на эту ПО в рамках оптимального трёхимпульсного апсидального перелёта.

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

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

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

Далее на основе полученных приближённых решений рассматривается задача быстродействия в рамках формализма принципа максимума Понтрягина [3]. Задача нахождения оптимального управления сводится к численному решению краевой задачи. При этом для получения невязок на правом конце траекторий межорбитального перехода внутри процедуры решения краевой задачи решается задача Коши.

Оптимальный трёхимпульсный перелёт

• •    второй импульс обеспечивает одновременное увеличение радиуса перицентра и поворот плоскости орбиты на угол ∆ i ₂ ;

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

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

На этапе выведения КА с помощью РБ используется биэллипти-ческая (трёхимпульсная) апсидальная схема перелёта (рис. 1). Положение и размер опорной орбиты фиксированы, а в качестве выбираемых параметров

∆ V ₁ ( R _α, ∆ i ₁ ) =

у/ V ⁺ V 2 - 2 VV ^ cos A г 1 ;

А V 2 ( R a , A i 2 ) = V ^V/^V F² V ^{V cos} A ⁱ ? ;

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

Использование трёхимпульсной схе-

∆ V 3 ( R α , ∆ i 1 , ∆ i 2 ) =

+ V ² - 2 V V cos(A i - A i_t - A i₂ ) , “ к              ^JL 2 “ к

где

V π _к

2ц

r x ^{+ r}& ^;

V 0 =

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

V π ₁( R α ) =

2ц

2ц

r ₁

R 1 + ^R J

последовательность:

• •    первый импульс переводит КА на переходный эллипс и поворачивает плоскость орбиты на угол ∆ i ₁ ;

\ 2ц

V α ₁( R α ) =    _R

а

2ц

^R а + Г 1’

Рис. 1. Схема трёхимпульсного перелёта (рисунок создан авторами)

V _π2( R _α) =

2ц

r

л

2ц

Гk + R а ’

\ 2ц     2ц

V ( R ) =       –        ;

оЛ ^a/ R R _а    R _а + г _к

V0 — круговая скорость на опорной орбите; μ — гравитационный параметр небесного тела; r1   — радиус опорной орбиты; Vπ1 и Vα1 — скорости

в перицентре и апоцентре первой переходной орбиты; Vπ2 и Vα2 — скорости в перицентре и апоцентре второй переходной орбиты; Vπк — скорость в перицентре промежуточной орбиты, с которой КА будет совершать «раскрутку» до целевой траектории; rπ и rα — радиусы перицентра и апоцентра промежуточной орбиты; ∆i — разность между наклонением целевой траектории и наклонением ПО; Rα — радиус апогея переходных орбит.

Тогда суммарное приращение скорости РБ можно определить по формуле (1):

∆ V ∑ ( R α , ∆ i 1 , ∆ i 2 ) = ∆ V 1 ( R α , ∆ i 1 ) +

+∆ V ₂ ( R _α, ∆ i ₂ ) + ∆ V ₃ ( R _α, ∆ i ₁ , ∆ i ₂ ).         (1)

Далее для определения минимума суммарного приращения скорости РБ (∆ V _∑ → min) как функции трёх переменных используется алгоритм SolvOpt ( Solver for local optimization problems ) [4]. В качестве метода оптимизации выступает алгоритм минимизации Шора [5].

После определения оптимальных параметров трёхимпульсного перелёта вычисляется масса КА на ПО m _ПО по соотношениям (2) и (3):

m _ПО = m ₀ – m _ТРБ – m _РБ ; (2)

m _ТРБ = m ₀

1 – exp

где m ₀ — начальная масса КА; m _ТРБ — масса топлива МДУ; m _РБ — «сухая» масса РБ; η — коэффициент потерь скорости, η = 1,025; w — удельный импульс тяги РБ с ЖРД.

Перелёт космического аппарата на рабочую орбиту с помощью электроракетной двигательной установки

На этапе перелёта КА с промежуточной эллиптической орбиты на рабочую траекторию с помощью ЭРДУ для вычисления значений конечной массы m _к и времени перелёта t _пер используется интерполяция

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

При известном значении массы КА на промежуточной орбите и требуемой характеристической скорости для совершения перелёта с ПО на конечную целевую орбиту с помощью ЭРДУ определяется конечная масса КА на целевой орбите и суммарное время перелёта по формулам (4) и (5).

m _к ( r _π, r _α, ∆ i ) = m _ПО exp

^V хар ^{( Г} . ’ ^Г a ’ ^{Д 1} )

^w ЭРДУ

^m ПО ^{- m} к ( ^R Я ^R ^ ^{Д 1} )

t _пер ( r _π, r _α, ∆ i ) =                        w _ЭРДУ , (5)

ЭРДУ где Vхар — интерполированная характеристическая скорость манёвра перехода с ПО на рабочую орбиту; mПО — масса единичного КА на ПО; wЭРДУ — удельный импульс тяги ЭРДУ; PЭРДУ — тяга ЭРДУ.

Оптимизации конечной массы с учётом ограничений

Для нахождения условного экстремума конечной массы m _к с учётом ограничений используется метод неопределённых множителей Лагранжа. Функция Лагранжа L будет иметь следующий вид:

L ( r _π, r _α, ∆ i , λ) = m _к ( r _π, r _α, ∆ i ) –

– λ · [ t пер ( r π , r α , ∆ i ) – t f ],                        (6)

где λ — множитель Лагранжа; t_f — ограничение на суммарное время перелёта.

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

д

L ( r _π, r _α, ∆ i ) = 0;

д г

Л

д

L ( r _π, r _α, ∆ i ) = 0;

< '                                               (7)

д

L ⁽r п , r а , ^А i ) = 0;

t _пер ( r _π, r _α, ∆ i ) – t_f = 0.

С помощью интерполяции B -сплайнами удаётся вычислить гладкие производные от функции Лагранжа по основным параметрам ПО, которые необходимы для оптимизации конечной массы КА за заданное время.

Для решения системы нелинейных уравнений (7) использовался гибридный метод Пауэлла (алгоритм hybrid [6]).

Формализм принципа максимума

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

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

tf

J = j dt ^ min.

t ₀

В качестве фазовых переменных, описывающих движение КА, выступают равноденственные элементы, которые не имеют особенностей в окрестностях нулевых значений наклонения и эксцентриситета [7]. Система дифференциальных уравнений движения КА с ЭРДУ в равноденственных элементах имеет вид (9).

dA dt

ξ f_T + ( e_x sin F – e_y cos F ) f_S

;

de

x

dt

(1 + ξ) cos F + e_x f_T + ξ f_S sin F + e_y η f_W ;

de y dt

(1 + ξ) sin F + e_y f_T – ξ f_S cos F – e_x η f_W ;

di

x

dt

фС08 F

f W ;

di y dt

dF dt

ф81П р

^^

f W ;

ξ ²

1    A Y

^ V Ц ⁽ i x sin F ^- i y cos F ) f w ,

где γ = 1 – e_x ² – e_y ² ;   ξ = 1 + e_x cos F + e_y sin F ;

η = ixsinF – iycosF; ϕ = 1 + ix2 + i2y; A — большая полуось; ex, ey, ix, iy — равноденственные элементы; e = ecos(ω + Ω); e = esin(ω + Ω); xy i = tg i cosΩ; i = tg i sinΩ; F — истинная x                     y долгота; t — время; fS, fT, fW — проекции реактивного ускорения на радиальное, трансверсальное и бинормальное направления орбитальной системы координат (рис. 2) соответственно; fS = fsinθ; fT = fcosθcosψ; fW = fcosθsinψ; θ — тангаж; ψ — рыскание (предполагается, что вектор тяги направлен по продольной оси КА); f =

^P ЭРДУ ^w ЭРДУ

^{– P} ЭРДУ ^{t + m} ПО ^w ЭРДУ

модуль

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

Рис. 2. Связь между геоцентрической экваториальной и орбитальной системами координат: r, n, b — рад_–иаль_–-ное, трансверсальное и бинормальное направления; T, S, W — проекции возмущающего ускорения на оси орбитальной системы координат; ^—r — радиус-вектор; u — аргумент широты; i — наклонение орбиты; Ω — долгота восходящего узла; ^— σ — вектор площадей (рисунок создан авторами)

Гамильтониан задачи быстродействия будет иметь вид (10).

H = –1 + λ A

dA dt +

de       d ^e

λ ^x + λ ^y e_x dt      e_y dt

+

8 H

8 H

8^

8 2 H

<

0;

0;

8^ 2

а 2 н	а 2 н
ае 2	аеа^	> 0
а 2 н	а 2 н
а^ае	а^ 2

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

k ₁ f_T = f _{| k|} |;

k ₂

f_S = f _|k| |;

f W

где

| k | =

di di x y

+ λ _{i x dt} + λ _{i y dt} + λ _{F dt}

k 1 =

где λ = [λ A λ e λ e λ i λ i λ F ] ^T xy x y

вектор

сопряжённых переменных к соответст-

k ₃

= f _|k| |,

2'' a A ^

+ ψ e_y

+ ψ _e x

sin F +

cos F + ^x

e_y + sin F

;

+ cos F

^

+

вующим равноденственным элементам.

Сопряжённую систему уравнений в силу её громоздкости запишем в общем виде, а именно:

dx = ah dt 8x ,

где ^– x = [ A e_x e_y i_x i_y F ] ^T — вектор фазовых переменных.

Закон управления вектором тяги определяется из условия максимума гамильтониана по управлению:

k„ = —A( e sin F - e cos F ) + v sin F — V cos F ;

2 у x y              ex            ey k 3   ^

n(V e e_v - V e e x - V f ) — ? (V i cos F + V i sin F ) . e_x y e_y x F ₂ i_x i_y

Интегрируя совместно системы уравнений (9) и (11), можно получить вектор невязок g на правом конце траектории:

g ( V o , t_f ) = [ A ( t f ) ^- A * e x ( t f ) e y ( t f ) i x ( t f )

где A *, i_x *, i_y * — требуемые значения параметров на правом конце траектории.

В качестве выбираемых параметров

— i* iy(tf) — i* VF(tf) H(tf)] T, рассматриваемой задачи выступает набор сопряжённых переменных на левом конце траектории и время перелёта:

V o = [V A ( t 0 ) V e x ( t 0 ) V e y ( t 0 ) V i x ( t o ) V i y ( t o') V f ( t o ^)] T ; t f .

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

g ( V o , t f ) = ⁰ .

Для решения системы нелинейных уравнений (13) использовался гибридный метод Пауэлла.

Анализ численных результатов

Были рассмотрены следующие виды космических транспортных систем:

• •  на  базе  РН  «Союз-2»  этапа  1а

и РБ «Фрегат»;

• •  на базе  РН  «Союз-2» этапа  1б

и РБ «Фрегат»;

• •    на базе РН «Ангара-А5» и РБ «Персей».

В качестве опорной орбиты искусственного спутника Земли (ОИСЗ) выступает круговая траектория высотой 2oo км и наклонением 51,8°. Целевая орбита также является круговой траекторией высотой 500 км и наклонением 2o°. Долгота восходящего узла, а также аргумент широты целевой орбиты принимаются равными 0°. В качестве ограничений на суммарное время перелёта был рассмотрен диапазон 50…100 сут полёта.

Предполагалось, что каждый единичный КА оборудован ЭРДУ на базе СПД-100Д.

Основные энергетические ристики РН, РБ и ЭРДУ [8–11]: характе- тяга МДУ РБ «Персей» 8 100 кгс; конечная масса РБ «Персей» 3618 кг; масса полезной нагрузки РН «Ангара-А5» на ОИСЗ 22 700 кг; удельный импульс тяги МДУ РБ «Персей» 353 c; тяга МДУ РБ «Фрегат» 2030 кгс; конечная масса РБ «Фрегат» 1 050 кг; удельный импульс тяги МДУ РБ «Фрегат» 333,2 с; масса полезной нагрузки РН «Союз-2» этапа 1а на ОИСЗ 7020 кг; масса полезной нагрузки

РН «Союз-2» этапа 1б на ОИСЗ 8 200 кг; тяга ЭРДУ 90 мН; удельный импульс тяги ЭРДУ 1 740 с.

На рис. 3–5 представлены зависимости от длительности выведения оптимальных параметров промежуточной орбиты, а также конечной массы КА, полученные в результате решения задачи условной максимизации. Для всех рассмотренных комбинаций средств выведения наблюдаются следующие характерные изменения: высота перицентра ПО монотонно убывает на рассматриваемом временном промежутке, а высота апоцентра, наклонение и конечная масса возрастают. Причём конечная масса монотонно возрастает почти линейно.

В сравнении с выведением только с использованием РБ были получены нижеследующие результаты.

Для комбинации РН «Союз-2.1а» с РБ «Фрегат» выведение только с помощью МДУ на целевую орбиту обеспечило доставку 975,6 кг полезной нагрузки. Таким образом, для комбинированной схемы выведения прирост по конечной массе составил 12, 17 и 22% для 50, 75 и 100 сут перелёта соответственно (рис. 3).

Для РН «Союз-2.1б» с РБ «Фрегат» выведение только с помощью МДУ на целевую орбиту обеспечило доставку 1 264,1 кг полезной нагрузки. Таким образом, для комбинированной схемы выведения прирост по конечной массе составил 12, 15 и 19% для 50, 75 и 100 сут перелёта соответственно (рис. 4).

Аналогично для РН «Ангара-А5» с РБ «Персей» выведение только с использованием МДУ РБ для группировки из шести аппаратов за раз составило 698,3 кг полезной нагрузки для каждого единичного КА (4 189,8 кг для космической системы в целом). Соответствующий прирост составил 11, 16 и 21% для 50, 75 и 100 сут перелёта соответственно (рис. 5).

в)

г)

Рис. 3. Изменение оптимальных параметров промежуточной орбиты (ПО) и массы космического аппарата (КА) при выведении с помощью ракеты-носителя «Союз-2» этапа 1а и разгонного блока «Фрегат»: а — высота апоцентра ПО; б — высота перицентра ПО; в — наклонение ПО; г — конечная масса КА (рисунок создан авторами)

а)                                                       б)

в)

г)

Рис. 4. Изменение оптимальных параметров промежуточной орбиты (ПО) и массы космического аппарата (КА) при выведении с помощью ракеты-носителя «Союз-2» этапа 1б и разгонного блока «Фрегат»: а — высота апоцентра ПО; б — высота перицентра ПО; в — наклонение ПО; г — конечная масса КА (рисунок создан авторами)

а)                                                    б)

в)

г)

Рис. 5. Изменение оптимальных параметров переходной орбиты (ПО) и массы космического аппарата (КА) при выведении с помощью ракеты-носителя «Ангара-А5» и разгонного блока «Персей»: а — высота апоцентра ПО;

б — высота перицентра ПО; в — наклонение ПО; г — конечная масса КА (рисунок создан авторами)

Далее из полученного множества результатов было выбрано решение для выведения РН «Союз-2.1б» с РБ «Фрегат» за 65 сут. Этот результат использовался в качестве начального приближения для решения задачи быстродействия в рамках строгой постановки принципа максимума Понтрягина.

На рис. 6 представлен характер изменения углов тангажа и рыскания оптимального управления в течение полёта. Для удобства угол рыскания был нормирован в диапазоне 0^360°.

а)                                                    б)

Рис. 6. Изменение углов тангажа и рыскания для оптимального решения задачи перелёта за минимальное время: а — за всё время полёта; б — за первые три часа полёта;--рыскание;--тангаж (рисунок создан авторами)

На рис. 7 представлена эволюция параметров оскулирующей перелётной орбиты и изменения массы КА во времени для оптимального решения. Время перелёта на целевую орбиту оказалось меньше на 1,5 сут по сравнению с ранее полученным приближённым решением, что обеспечило увеличение массы полезной нагрузки КА на целевой орбите на 2 кг.

s и ч о с к

-а

а)

S =

3 о ч св

К

в)

Рис. 7. Изменение параметров орбиты и массы космического аппарата (КА) в течение оптимального перелёта с электроракетной двигательной установкой: а — большая полуось; б — эксцентриситет; в — наклонение; г — масса КА (рисунок создан авторами)

г)

Выводы

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

Анализируя результаты, можно отметить, что выведение КА посредством РН «Союз-2» с РБ «Фрегат» является наиболее уместным при запуске единичного КА, а для группового пуска — РН «Ангара-А5» с РБ «Персей».

При использовании комбинированных схем прирост массы полезной нагрузки для РН «Ангара-А5» составил до 21% в сравнении со схемами выведения, использующими только возможности РБ с ЖРД. Для РН «Союз-2» аналогичный показатель составляет до 22%.

Исследование выполнено за счёт гранта Российского научного фонда № 22-79-10206.