CC BY
98
Том XXXVIII
УЧЕНЫЕ ЗАПИСКИ ЦАГИ 20 0 7
№ 3 — 4
УДК 629.78.015.076.8:525.7
УПРАВЛЕНИЕ ВОЗВРАЩАЕМЫМ В АТМОСФЕРУ КОСМИЧЕСКИМ АППАРАТОМ НА НИЖНЕМ УЧАСТКЕ ТРАЕКТОРИИ
А. В. БОБЫЛЕВ, В. А. ЯРОШЕВСКИЙ
Рассматривается движение космического аппарата с умеренным аэродинамическим качеством на нижнем участке траектории входа в атмосферу, соответствующем малым скоростям полета (V < 1000 м/с). Оценивается точность аналитических формул,
описывающих изменение параметров траектории на этом участке. На основе анализа геометрических свойств траектории предлагается алгоритм определения момента выполнения последнего разворота аппарата по крену. Приводятся результаты статистического моделирования движения аппарата, управляемого только путем регулирования угла крена, по определению точности прибытия в расчетную точку раскрытия парашюта на высоте 6 км.
При синтезе законов управления возвращаемым аппаратом (ВА) необходимо учесть высокие требования по точности приведения ВА в точку, где предполагается раскрытие парашюта: ошибка не должна превышать нескольких километров при условии, что скорость ВА не должна быть больше 200 — 250 м/с. Исходя из этого, рассмотрим движение ВА на нижнем участке траектории. Будем считать, что управление траекторным движением ВА осуществляется по закону, близкому к использованному при управлении космическими летательными аппаратами «Шаттл» и «Буран» ([1, 2]), основанному на отслеживании программ изменения продольной перегрузки в зависимости от скорости.
Естественным требованием является выбор такого семейства программ пха (V), которое
обеспечивает полет ВА на заключительном этапе траектории с некоторым средним программным углом крена, близким к постоянному (например, 45°), имея в виду возможность компенсации ошибок по дальности непредсказуемого знака. Тогда и программное значение эффективного
суа ^ у
аэродинамического качества Ке =------------------------- следует выбрать близким к постоянному.
сха
Приближенные аналитические оценки продольной перегрузки пха и угла наклона траектории 0 для условий спуска в экспоненциальной атмосфере при плавном изменении коэффициентов подъемной силы и сопротивления определяются формулами (2.101) и (2.103) из [1]:
(Пха )а
1 + 2 К „и + и
1 + к2 + 2Ке (і + К2 )и + К2еи2 ’
(1)
( sin 0)
V fa
2 - d ln nxa
d ln u
d ln cxa
d ln u
2 (1
u 22)-d 1-nxa
d ln u
d ln cxa
d ln u
Здесь использованы обозначения: u =
V
X — логарифмический градиент плотности
атмосферы, Ke = K cos у, K — аэродинамическое качество, у — скоростной угол крена. Формула (2) обобщена на случай cxa Фconst.
Проверим точность этих аналитических оценок, сопоставляя их с результатами численного интегрирования точных уравнений движения на нижнем участке в предположении об экспоненциальной модели плотности атмосферы в диапазоне уменьшения значения u от 4 до 0.5. Начальные условия при u = 4 зададим в соответствии с формулами (1) и (2).
Результаты сопоставления при Ke = const, cxa = const изображены на рис. 1 — 4 для
Ke = 0.25, 0.5, 1, 2 в виде зависимостей от безразмерной скорости u значения продольной
перегрузки nxa и значения (—sin 0). Максимальная погрешность формул при определении
значения nxa изменяется от 6% при Ke = 0.25 до 2.5% при Ke = 2, а при определении значения
(—sin 0) — от 0.035 до 0.005, что можно признать удовлетворительным. По мере увеличения Ke
кривые, полученные в результате численного расчета, приобретают отчетливо выраженный колебательный характер. В результате, гладкие кривые, построенные по формулам (1) и (2), приближенно определяют «невозмущенное» решение уравнений движения, поэтому их удобно использовать в качестве номинальных зависимостей.
По мере уменьшения скорости аэродинамическое качество, как правило, возрастает. Примем, что значение Ke линейно изменяется в зависимости от скорости. В этом случае ошибки несколько возрастают. Обозначим Ke (u = 4) через K1, а Ke (u = 0.5) — через K2. Тогда, при постоянном значении cxa , убедимся, что указанная максимальная ошибка в определении значения nxa составляет 4% при K1 = 0.5, K2= 1; 13% — при K1 = 0.5, K2 = 2; 3% при — K1 = 1, K2 = 2 и
K2 = 3. Некоторые результаты расчетов проиллюстрированы на рис. 5 — 7.
Если возрастание аэродинамического качества сопровождается уменьшением угла
Рис. 1. Сравнение зависимостей пха (и) и -8т0 (и)
подсчитанных по формулам (1) и (2) (жирные линии) и полученных численным интегрированием (тонкие линии),
Ke = 0.25
Рис. 2. То же, что на рис. 1, но Ke = 0.5
атаки, то оно частично обуславливается уменьшением коэффициента сопротивления. Примем, что этот коэффициент линейно изменяется в зависимости от скорости. Рассмотрим случай, когда cxa (и = 0.5)/cxa (и = 4) = 0.5. Тогда при K1 = 0.5, K2 = 2 максимальная ошибка в определении nxa составляет 12%; при K1 = 1, K2 = 2 — 5%; при K1 = 1, K2 = 3 — 4%. Соответствующие результаты расчетов проиллюстрированы на рис. 8, 9.
0.475
0.47
0.465
0.46
0.455
0.4
—sin0
0.3
0.2
0.1
1 /
Г
/
0 1 U 4
Д. ■1 I
1 . .1
0 ; 3 4
Рис. 3. То же, что на рис. 1, но Ke = 1
Рис. 4. То же, что на рис. 1, но Ke = 2
Рис. 5. Сравнение зависимостей пха (и) и -8т0 (и),
подсчитанных по формулам (1) и (2) (жирные линии) и полученных численным интегрированием (тонкие линии), К = 0.5, к2 = 1, Сха (0.5)/Сха (4)= 1
п
X
0.8
0.7
0.6
0.5
0.4
-эшО
0.3
0.2
0.1
0 2 3 и 4
о
2 3
Рис. 6. То же, что на рис. 5, но К = 1, К2 = 2, Сха (0.5 VСха (4)= 1
Рис. 7. То же, что на рис. 5, но К = 1, К2 = 3, Сха (0.5 VСха (4)= 1
Рис. 8. То же, что на рис. 5, но К! = 1, К2 = 2, Сха (0.5)/Сха (4)= 0.5
Точность приближенных формул (1) и (2) можно оценить и другим путем: вычислением зависимости эффективного аэродинамического качества от скорости на основе решения обратной задачи динамики полета. Пренебрегая зависимостью гравитационного ускорения от высоты и считая модель плотности атмосферы экспоненциальной, можно показать, что указанная зависимость описывается следующей формулой:
(l - V 2 ) cos 0(v ) nx (V)
(nx (v) + sin0(V))v dsin0(V) nx (V)cos0(V) dV
(3)
где sin 0 определяется формулой (2).
Результаты расчетов по формуле (3) приведены в табл. 1.
Как видно, при не слишком малых значениях эффективного аэродинамического качества (Ke > 0.75) точность можно считать удовлетворительной.
Т аблица 1
Ke = 0.25 Ke = 0.5 Ke = 0.75 Ke = 1 Ke = 1.5
u Keref Keref Keref Keref Keref
4 0.218 0.491 0.760 1.023 1.537
3.5 0.207 0.478 0.748 1.012 1.527
Рис. 9. То же, что на рис. 5, но Ki = 1, K2 = 3, ^(0.5)/^(4) = 0.5
3 0.195 0.463 0.733 0.998 1.516
2.5 0.185 0.448 0.718 0.983 1.503
2 0.184 0.443 0.710 0.975 1.495
1.5 0.214 0.466 0.727 0.987 1.501
1 0.294 0.529 0.776 1.026 1.526
0.5 0.311 0.544 0.786 1.032 1.526
Имея в виду в дальнейшем исследовать некоторые свойства пространственного движения ВА, определим зависимость от скорости радиуса кривизны горизонтальной развертки траектории
тт-2 2 ґ\ 2 2 /л
V cos 0 и cos 0
gnya sin Y XKnxa sin Y
(4)
Форма горизонтальной развертки траектории, как будет показано ниже, напоминает логарифмическую спираль (рис. 10). Уравнение этой спирали в полярных координатах записывается в виде:
r = r0 exp (k ф),
где k = const, причем радиус кривизны определяется формулой
R = г>/1 + k2,
а пройденный путь от первой до второй точки — формулой
5 = (r - r2 )л/ 1 + k2 /k.
(5)
(6)
(7)
Отметим еще, что угол между вектором скорости и направлением на начало координат сохраняет постоянное значение, равное
ф = arctan 1/ к,
а расстояние между двумя точками определяется формулой
Лг = R1 sin ф1 + R2 sin ф2,
(8)
(9)
где ф1; ф2 — углы между прямой, соединяющей две точки, и направлениями вектора скорости в этих точках: знаки выбраны так, чтобы при малых расстояниях эти углы были положительны и стремились к нулю при Аг - > 0 (рис. 10).
В реальности горизонтальная развертка траектории отличается от логарифмической спирали: «местное» значение параметра спирали не является постоянным и описывается формулой
к = -
sin0(a,V-2&cos0) sin0(u2 -2Kenxa + 2cos0)
nxaK sin Y
(10)
Оценим характер изменения этого безразмерного параметра на нижнем участке траектории, подставляя в формулу (10) значения sin 0 и nxa в зависимости от безразмерной скорости, вычисленные по результатам интегрирования уравнений движения. Указанный параметр зависит не только от значения Ke, но и от значения sin y. Рассмотрим для примера случай у = 45°. В табл. 2 приведена определенная таким образом зависимость к (и) для различных значений Ke.
Таблица 2
и 4 3.5 3 2.5 2 1.5 1 0.5
Ke = 0.5 3.31 3.20 3.07 2.94 2.79 2.58 1.91 0.520
Ke = 1 1.73 1.72 1.63 1.51 1.37 1.19 0.805 0.243
.5 N e Ke 1.17 1.18 1.10 0.991 0.896 0.767 0.490 0.156
Ke = 2 0.884 0.902 0.818 0.733 0.678 0.559 0.356 0.115
Ke = 4 0.444 0.442 0.392 0.389 0.331 0.268 0.169 0.056
Если при том же значении эффективного аэродинамического качества рассмотреть другие значения угла крена, то в соответствии с формулой (10) оказывается, что функцию к (и) следует
умножить на ctgy. Как видно, эта функция примерно пропорциональна 1/Ke, убывает по мере
уменьшения скорости и при и > 2 изменяется не очень значительно. Отсюда следует, что горизонтальная развертка траектории не очень сильно отличается от логарифмической спирали. Незначительные нарушения монотонности при Ke = 1.5 и 2 объясняются колебаниями функций
nxa (и), sin0 (и), см. рис. 4. При подстановке в формулу (10) оценок этих функций (1) и (2)
значения к (и) почти не изменяются и
приобретают строгую монотонность.
В табл. 3 приведены также значения радиуса кривизны в километрах при у = 45°, подсчитанные по формуле (4), где принято 1/ X = 7 км. В реальности на малых высотах значение 1/X возрастает до 8 — 10 км. При других значениях угла крена, в соответствии с формулой (4), оказывается, что функцию R (и) следует умножить на ctgy. Значения этой функции быстро убывают при малых Ke и медленно убывают при больших Ke.
Рис. 10. Логарифмическая спираль
Т аблица 3
и 4 3.5 3 2.5 2 1.5 1 0.5
Ke = 0.5 153 120 90.5 64.7 42.0 22.5 7.32 0.93
Ke = 1 130 99.6 73.3 51.2 33.0 18.1 6.83 1.34
.5 N e Ke 122 93.0 67.8 47.3 30.7 16.9 6.84 1.52
2 II e Ke 119 90.0 65.5 46.0 29.8 16.3 6.91 1.61
Ke = 4 114 86.2 64.0 44.4 28.5 15.9 6.97 1.71
При использовании традиционной схемы управления боковой дальностью путем регулирования только угла крена ВА предусматривается выполнение перекладок по крену с изменением знака угла крена в зависимости от относительного угла азимута — угла между горизонтальной проекцией вектора скорости и направлением на желаемую точку приведения. Если продлить эту схему на диапазон малых скоростей, то количество перекладок существенно возрастает, поскольку с уменьшением скорости ВА становится слишком «податливым» в смысле возрастания угловых скоростей разворота по курсу. С учетом того, что при выполнении поворота по крену на немалый угол с конечной угловой скоростью (порядка 10°/с) возникает значительный импульс эффективной подъемной силы, управление боковым движением нарушает процесс управления продольным движением. Поэтому возникает потребность в ограничении количества перекладок, которая выливается в требование о том, что при достижении некоторой скорости на оставшемся интервале полета допускалось выполнение только одной перекладки по крену.
Рассмотрим движение гипотетического возвращаемого аппарата при постоянном значении угла крена у = 45° на нижнем участке траектории после выполнения последней перекладки
по крену. Зададим следующие значения для параметров ВА: К = с уа1сха = 1, сха8/ш = 0.001 м 2 /кг.
Начальные условия зададим на высоте 30 км, начальную скорость примем равной 1048 м/с, что соответствует значению и = 4, начальный угол наклона траектории зададим равным -6.72°. Начальные значения угла наклона траектории и продольной перегрузки при этом соответствуют вычисленным по формулам (1), (2). Конечное значение скорости, при котором происходит раскрытие парашюта, зададим равным 200 м/с.
Результаты расчета такой траектории приведены в табл. 4. В ней содержатся вычисленные с шагом 10 с значения следующих параметров: Ь, I, ф, Н, V, г, ге, я, se. Эти параметры обозначают соответственно продольную дальность, боковую дальность, угол курса, высоту, скорость, истинное горизонтальное расстояние между текущей и конечной точкой, оценку горизонтального расстояния между текущей и конечной точкой, истинный оставшийся путь до конечной точки, оценку оставшегося пути до конечной точки.
Таблица 4
t L 1 Ф H V r re s se
0 0 0 0 30.00 1048 59.76 72.15 70.39 69.30
10 9.90 0.37 4.46 28.73 951.5 51.25 61.81 60.48 59.93
20 18.79 1.43 9.31 27.39 858.8 43.63 52.36 51.53 51.47
30 26.63 3.08 14.63 25.98 770.6 36.86 43.85 43.51 43.89
40 33.44 5.22 20.49 24.51 687.4 30.89 36.30 36.37 37.11
50 39.20 7.74 26.96 22.99 609.9 25.67 29.65 30.08 31.09
60 43.96 10.53 34.16 21.41 538.9 21.13 23.87 24.56 25.74
70 47.74 13.49 42.22 19.80 474.9 17.20 18.88 19.76 21.00
80 50.58 16.49 51.32 18.15 418.2 13.80 14.64 15.62 16.79
90 52.56 19.45 61.70 16.49 368.8 10.86 11.07 12.06 13.06
100 53.72 22.25 73.67 14.82 326.3 8.30 8.13 9.02 9.76
110 54.16 24.80 87.58 13.17 290.0 6.06 5.73 6.42 6.87
120 53.96 27.00 103.8 11.56 259.0 4.06 3.76 4.20 4.41
130 53.23 28.74 122.8 10.00 232.3 2.27 2.10 2.31 2.36
140 52.11 29.93 144.8 8.52 209.4 0.668 0.631 0.669 0.671
145 51.53 30.26 155.9 7.88 199.0 0 0 0 0
Как отмечалось выше, в рассматриваемом случае радиус кривизны горизонтальной
развертки траектории убывает по мере спуска. В качестве второй точки принимается конечная
точка, определяемая достижением скорости 200 м/с.
Здесь грубая оценка расстояния между текущей и конечной точкой определяется формулой:
A/ «(R + R2)sinф^ (11)
Приближенная оценка оставшегося пути до конечной точки определяется формулой
As1 «Ar/cos ф1. (12)
Проанализируем вкратце полученные результаты. Оценка расстояния до конечной точки по упрощенной формуле (11) имеет удовлетворительную точность при скорости, не превышающей 500 м/c. Оценка оставшегося пути по формуле (12) имеет более высокую точность.
Задавая несколько значений угла крена у2 после выполнения перекладки, близких к исходному (в данном случае — 45°), вычисляем в процессе полета (при отрицательном значении
угла крена) значения оценок оставшегося расстояния и оставшегося пути по формулам (11) и (12) и приравниваем их соответственно к текущему расстоянию до требуемой конечной точки и к программному значению оставшегося пути, вычисляемому при данном значении у который определяется формулой:
Í
cos
0(V )vdv
nx (v) + sin0(V)’
(13)
где sin 0 определяется формулой (2).
В итоге получим два соотношения:
(R (y2 ) + Rk (y2 )) sin ф1 = ЛГ
Aí(y2) = Лг/cosф1.
(14)
Здесь значение Я (у 2) в соответствии с формулой (4) пропорционально 1/|8т у 2|, малое значение Як представляет собой радиус кривизны в конечной точке, а функция Аз (V, у2) определяется по формуле (13). Момент перекладки tp определяется условием, когда
выполняются оба равенства (14) при некотором значении у 2.
Для проверки применимости предложенного алгоритма определения последнего момента перекладки по крену и отыскания необходимого значения У2 были проведены расчеты управляемого траекторного движения космического летательного аппарата с умеренным аэродинамическим качеством типа «несущий корпус». Алгоритм управления на верхней части траектории был выбран подобным использованному на ВКС «Буран» [2]. Результаты статистического моделирования управляемого траекторного движения (Ы = 1000) для некоторой характерной точки внутри области достижимости с учетом расчетных возмущений (и без учета навигационных ошибок) проиллюстрированы на рис. 11 в виде отклонений продольной и боковой угловой дальности точки раскрытия парашюта, определяемой достижением высоты 6 км. Средняя величина промаха составила около 1.5 км.
км
1Л t . »
см * ••
-5 ,0 -2 5 Я£5.'*2 Я» ,5 5 • .0
км
Рис. 11. Статистическое моделирование. N = 1000. Влияние ветра, отклонений плотности атмосферы, разброса начальных условий и аэродинамического качества Отметим, что при моделировании траекторного движения учитывалось ограничение на конечную угловую скорость выполнения перекладок по крену с изменением его знака & = 10°/ с, в то время как при выводе соотношений (14) предполагалось, что перекладки выполняются мгновенно. Тем не менее, результаты можно признать удовлетворительными, поскольку при
неуправляемом спуске на парашюте с высоты 6 км под действием ветра отклонение точки посадки может заметно превысить указанную величину промаха. При необходимости можно усложнить логику определения момента выполнения последней перекладки с учетом конечной угловой скорости поворотов по крену.
Работа выполнена при поддержке Российского фонда фундаментальных исследований, проект № 06-01-08053-офи.
ЛИТЕРАТУРА
1. Ярошевский В. А. Вход в атмосферу космических летательных аппаратов. —
М.: Наука, 1988.
2. Ярошевский В. А. Алгоритмы управления траекторным движением космических летательных аппаратов на этапе спуска в атмосфере // Аэрокосмическая техника и технология. 1999. № 1.
Рукопись поступила 4/V 2006 г.
