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

Том ХЫН

УЧЕНЫЕ ЗАПИСКИ ЦАГИ 2012

№ 3

УДК 629.735.33.015.075

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

А. В. БОБЫЛЕВ, В. А. ЯРОШЕВСКИЙ

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

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

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

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

V =

р соэ ( + фдВ )-сшдЗ

тё

_ /

1 _ Р эт ( + фдв ) + судЧ8

тё

(1)

БОБЫЛЕВ Анатолий Владимирович

кандидат технических наук, начальник сектора ЦАГИ

ЯРОШЕВСКИИ Василий Александрович

доктор технических наук, член-корреспондент РАН, советник дирекции ЦАГИ

где V — путевая скорость самолета; Р — тяга двигателей; m = 100 000 кг — масса самолета; q — скоростной напор; S = 168 м2 — характерная площадь; g = 9.81 м/с2 — ускорение свободного падения; cxa = 0.105, cya = 0.5 — коэффициенты сопротивления и подъемной силы; f = 0.05

(среднее значение между сухим бетоном и заснеженной полосой) — коэффициент трения качения по ВПП; Фдв — угол установки двигателя; угол наклона ВПП принят равным нулю. Будем

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

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

P = P0 (1 - kvVw),

где для рассматриваемого самолета коэффициент kv = 0.0021 с/м, Р0) = 270 000 н. Наряду с этим вариантом предварительно рассматривается более упрощенный вариант Р = const = 250 000 н.

Все параметры, входящие в уравнение (1), за исключением скорости V и тяги двигателей P в процессе разбега, будем считать постоянными. Это допущение справедливо, по крайней мере, до момента достижения скорости принятия решения V . Именно к этому моменту необходимо уточнить параметры движения самолета. Следует отметить, что V — воздушная скорость самолета. Таким образом, уравнение (1) можно преобразовать к следующему виду:

V =

(Р/mg - f) - B (cxa - fcya )], (2)

здесь В = qS/mg — безразмерный параметр.

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

V =

отр

2mg = 78.9 м/с. (3)

pSc

У отр

В процессе разбега измеряются параметры: L (расстояние от точки начала разбега), скоростной напор q, связанные перегрузки nx и Пу. Считается, что путевая скорость Vизвестна точно.

Измеренное значение скоростного напора qизм можно записать в виде:

^зм = 05pVW, (4)

здесь Vw = V - W — воздушная скорость самолета; W — скорость ветра; р = 1.225 кг/м2 — плот-

ность атмосферы. Измеренные значения перегрузок:

Пхизм = pmg - f - B(cxa - fcya), (5)

Пу изм = cyaB. (6)

Вначале рассмотрим вариант Р = const. Длина разбега L в зависимости от скорости V может быть получена по формуле:

L (V ) = - S ( m f ) ‘n

pS (cxa fcya )

1 PS ( cxa - fcya )V 2

2mg (Р - f)

(7)

Для номинальных условий длина разбега при достижении скорости VI составляет I (ух ) = 1193 м, г (ух ) = 35.5 с, Ь (Котр ) = 1860 м, t (Котр ) = 44.76 с.

Поскольку определяющим моментом является достижение воздушной скоростью Vw скорости принятия решения VI, то скорость V целесообразно использовать в качестве независимой переменной. Четыре измеряемые параметра движения самолета могут быть использованы для уточнения четырех наиболее важных параметров движения самолета. В качестве таких параметров выберем Ж, P, m и / Отклонения этих параметров от номинальных значений составят вектор состояния 5Г, причем параметры ДР, Дm и Д/ целесообразно сделать безразмерными, разделив их

на номинальные значения. Таким образом, вектор состояния 5Г включает Ж, ДР, Дт и Д/ .

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

<7ном (V) пхном (У), пу ном (V) и Аюм У). Затем продифференцируем измеряемые параметры

по элементам вектора состояния, чтобы получить линейную систему уравнений:

Производные измеряемых параметров по элементам вектора состояния образуют матрицу

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

, остальные

дпх 2^ном^(Сха /суа) дпх Р

(8)

дЖ mgV ’ дДР mg ’

(9)

дпу _ 2^

дЖ ~

(10)

ю

ю

(11)

тёпх ном

^пхном ’ дД/ I пхном V

mg )

(12)

При P ф const уравнение (2) принимает вид:

V = :

Р ( - )/тё - / - В (Сха - /суа )]. (13)

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

дпх ^Р) 2“ном^ (Сха - /суа ) дпх р (1 — kvV)

dW mg mgV ’ SAP mg

(14)

дт = (ха - /Суа)-/>° (1 к"У) , д/ = _/(1 -п ). (15)

дДт т^ у ' mg дД/ 4 '

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

дЬ Г pkVV - 2“номЛ' (Сха - /Суа)Л дЬ \ -УР0 (1 - к^) (1б)

дЖ Г ^пхном , дДР ,0 mgnxном ,

dL tV

L_ Г I/0 (1 -kvV)-(c*a - fcya)]dt dL_ = - Vf Г - cyaqS (17)

dAm : mgnXHOM ’ dAf tn nx ном Г

mg )

Итак, в момент времени ^ производится измерение перечисленных выше параметров, результаты которого можно представить в следующем виде [3]:

5'г = С5Г + бёк, (18)

где 5'г — вектор отклонений от номинальных значений; 5"гк — вектор ошибок измерений, которые имеют нулевые математические ожидания и образуют корреляционную матрицу Я. Будем предполагать, что компоненты вектора ошибок измерений некоррелированы между собой, тогда матрица Я будет диагональной. В работе принято, что диагональные элементы матрицы Я имеют значения

Лп = 100 н2/м4, Я22 = Я33 = 0.0001, Л44 = 1 м2.

Что касается точности измерения дальности при разбеге по ВПП, то в настоящее время

возможно использование Дифференциальной глобальной системы позиционирования (БОР8), которая обеспечивает указанную точность [4].

Оценку вектора состояния перед измерением в момент времени и - 0 будем обозначать

51-, а соответствующую корреляционную матрицу К_. После измерения оценка вектора состояния и корреляционная матрица приобретают новые значения 5Т+ и К + из условия минимума среднеквадратической ошибки [3]:

51+ = 51+ М (5ё - С51), (19)

К + = (стЯ-1С + К-1 )-1, (20)

М = К-Ст (Я + СК-Ст )-1. (21)

Параметры дискретного алгоритма представляют собой рекуррентные соотношения. Для вычисления элементов вектора 8r+ необходимо задать начальные значения вектора 8F , которые можно принять равными нулю, и начальные значения элементов корреляционной матрицы K _. Выпишем начальные значения диагональных элементов матрицы K _.

К11_ = 1 м2/с2, К22_ = К33_ = К44_ = 0.001.

Изменение элементов матрицы G в процессе разбега самолета по ВПП при P = const показано на рис. 1 — 4. Обращает на себя внимание тот факт, что все элементы в процессе разбега увеличиваются по модулю, за исключением элементов G22, G23 и G24. Элемент

dnx = _P_ dP mg

Предположим, что в процессе разбега самолета по ВПП отклонения параметров W, AP Am и А/ от номинальных значений составляют:

G22 =

• = const.

Ж = -1 м/с, АР = 0.05, АШ = -0.05, А/ = 0.1.

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

ч

О

-2.5

-5.0

-7.5

О

ч

0.2

0.1

о

-0.1

-0.2

-0.3

Рис. 1. Элементы матрицы G

т'\1 ) 2 0 5 *» 0 б ч 0 ' Ю Г. м

V Ч ч Ч k G4I \

ч\ °и ' \ \ ч > \

\ \

п

22

] 0* 7 О-—-.. Р. 4 ) 5 0 6 _ rz 0 7 К»

о "V,. "-"ч # ч

> Ч

м/с

Рис. 2. То же, что на рис. 1

о

-0.01

-0.02

-0.03

-0.04

-0.05

^21

“‘1 (Т 1 г---? 0.__4 ■'*1 . 1 1 О» 0 6 0 7 (] 8 1

Сз, 1

Г м/

. -■

• ' —'

Рис. 3. То же, что на рис. 1

О.

>ч

2000

1000

-1000

-2000

#

043 ✓ 1 ** * ✓ * О* ✓

„ •* -*• С/ 44 •—

1 ) 2 ?1. ^4 1) 5 • 0 6 0 7 \

о." ' -ч * ч ^ Ч

ч ч

Рис. 4. То же, что на рис. 1

№, м 100

50

0

-50

-100

-150

Рис. 5. Измерения (Ш = -1 м/с, АР = 0.05, Ат = -0.05, А/ = 0.1)

На рис. 7, 8 представлено изменение диагональных элементов матрицы М, являющейся матрицей коэффициентов усиления дискретного алгоритма. Элемент Мп вначале увеличивается примерно по линейному закону, достигает максимума (по модулю) при t«1.5 с, а затем уменьшается до малых значений по закону, близкому к экспоненциальному. Элемент М22 (0) =1.1 и при t > 0 также уменьшается до малых значений примерно по экспоненциальному закону. Этот элемент определяет характер изменения оценки тяги двигателя. Элемент М 33 достигает максимума (по модулю) лишь при t « 17.5 с, он определяет характер изменения оценки массы самолета

Л<7, II/? л

Л<7 _ Ь „Д. 4/'? ' Л '■- V' ** т? &Ч--

и*.г- ■ г&.“/ь ' & 5 1 Дя..& А: V V \у’

V ’ ь_.л 5 2 ~ ^ 0 2 5 3 о ■ 5 1,

* л «■»

АЬ ** ч V 4

Рис. 6. То же, что на рис. 5

М,

.1 1 1 1

\ \ ч N к. 1Г:

Л / А- (и I 5 2 0 2 5 1 ми в * ?'— /,

1.0

0.5

Рис. 7. Диагональные элементы матрицы М

М..

-0.001

-0.002

1 Г

1 5 2 ) 2 5 3 ) : 5 /, с

у/44

Рис. 8. То же, что на рис. 7

и свидетельствует о том, что постоянная времени в оценке массы самолета примерно на порядок превышает постоянную времени в оценке тяги. Элемент М44 в процессе разбега самолета меняет знак (при ?« 10 с) и достигает наибольших значений лишь к окончанию разбега (при ?« 28 с), что свидетельствует о трудностях уточнения параметра А/ .

Изменение диагональных элементов корреляционной матрицы К показано на рис. 9, 10. Элементы Кц, К22 и К33 к моменту окончания разбега уменьшаются до малых значений,

а элемент К44 уменьшился по сравнению с К44 (0) примерно в 3 раза.

к.

*1

1.0

0,8

0.6

0.4

0.2

и„

1"'

5 10 15 20 25 30 35

Рис. 9. Диагональные элементы матрицы К

Г, с

0.102

хЮ-2

0.051

1 V ч \ .Ки \

X \ К22 \ V' ч \ ^3

10 15 20 25

Рис. 10. То же, что на рис. 9

30

35

I. о

Рис. 11. Апостериорная оценка Ш

На рис. 11, 12 представлено изменение апостериорных оценок всех параметров в зависимости от скорости разбега. Значение скорости ветра определяется быстро и с большой точностью. Этого и следовало ожидать, поскольку скоростной напор не содержит информации о других определяемых параметрах. Оценки тяги и массы самолета определятся достаточно хорошо, а оценка А/ лишь к окончанию разбега переходит через ось абсцисс и принимает значение нужного знака. Совершенно очевидно, что для повышения точности оценок необходимо увеличить частоту измерений.

0.1

0.05

-0.05

л,,.,

/*| д/ *• I

¥ V У4 у» АР *v'./vV' .'■'-■S'. ,\l >/ * •••••• /V5* Г1

* О---2 0 Ат г » "'«■Ч 0*4 * -*4 .f\ 0* 3 0 6 0 7 ) г ,, мУ

V *v'”

Рис. 12. Апостериорные оценки АР, Ада, Af

Для получения статистических характеристик процесса определения неизвестных параметров было проведено статистическое моделирование в объеме 1000 реализаций при А( = 0.2 с. Результаты расчетов для двух точек траектории разбега (в момент достижения скорости принятия решения ¥1 и при достижении скорости отрыва Уотр) представлены в табл. 1.

Т аблица 1

Результаты статистического моделирования, P = const, At = 0.2 с

V X Xmin Xmax X CT

W, м/с -1.0664 -0.9508 -1.0113 0.01638

Vi АР 0.0287 0.0574 0.0428 0.00427

Am -0.0598 -0.0263 -0.0442 0.00526

Af -0.0025 0.0945 0.0457 0.01513

W, м/с -1.0467 -0.9710 -1.0099 0.01289

V отр АР 0.0371 0.0519 0.0442 0.00238

Am -0.0605 -0.0299 -0.0475 0.00398

Af 0.0196 0.1172 0.0705 0.01410

В первом столбце выписаны значения скорости разбега, при которых фиксировались апостериорные оценки неизвестных параметров, во втором столбце — определяемые параметры Ж, АР, Ат и А/, в третьем и четвертом — минимальные и максимальные значения этих параметров Хт;п, Хтах, в пятом — математическое ожидание X, в шестом — среднеквадратическое отклонение определяемых параметров с . Из таблицы следует, что скорость ветра определяется с высокой точностью, а все остальные параметры имеют значительные отклонения от истинных значений. Причем точность определения параметров к моменту окончания разбега заметно повышается по сравнению с моментом достижения скорости ¥1, особенно это заметно по оценке коэффициента трения качения А/.

Для того чтобы повысить точность определения неизвестных параметров, было проведено статистическое моделирование с увеличенной частотой измерений при А( = 0.1 с. Результаты расчетов представлены в табл. 2.

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

Результаты статистического моделирования, P = const, At = 0.1 с

V X Xmin Xmax X CT

W, м/с -1.0471 -0.9696 -1.0117 0.01136

V1 AP 0.0317 0.0547 0.0430 0.00327

Am -0.0593 -0.0340 -0.0474 0.00415

Af 0.0204 0.1151 0.0639 0.01502

W, м/с -1.0360 -0.9808 -1.0098 0.00889

V отр AP 0.0390 0.0495 0.0440 0.00165

Am -0.0600 -0.0410 -0.0503 0.00308

Af 0.0381 0.1251 0.0840 0.01168

В дальнейшем статистическое моделирование проводилось с учетом уменьшения тяги двигателей в процессе разбега самолета, объем статистических испытаний составил 1000 реализаций. Результаты расчетов для At = 0.1 с показали, что точность определения параметра Af так же, как

и для случая Р = const, At = 0.2 с (табл. 1), недостаточна. Удовлетворительные результаты удалось получить, повысив частоту измерений до 20 Гц. Результаты расчетов приведены в табл. 3.

Таблица 3

Результаты статистического моделирования, P Ф const, At = 0.05 с

V X X . min X ^-max X CT

W, м/с -1.0258 -0.9993 -1.0119 0.00412

V1 AP 0.0345 0.0516 0.0434 0.00266

Am -0.0550 -0.0386 -0.0472 0.00274

Af 0.0194 0.1140 0.0646 0.01530

W, м/с -1.0236 -0.9987 -1.0098 0.00317

V отр AP 0.0426 0.0496 0.0461 0.00128

Am -0.0558 -0.0436 -0.0502 0.00200

Af 0.0563 0.1408 0.0955 0.01324

Одна из основных задач, связанных с разбегом по ВПП, состоит в выявлении нештатных ситуаций, часто приводящих к аварии. В качестве примера рассмотрим невыпуск по какой-либо причине закрылков, из возмущений сохраним горизонтальный ветер Ж = -1 м/с. Результаты расчетов представлены в табл. 4.

Таблица 4

Результаты статистического моделирования, P Ф const, At = 0.05 с

V X Xmin Xmax X CT

W, м/с -1.0250 -0.9977 -1.0104 0.00394

V1 AP 0.2512 0.2697 0.2611 0.00264

Am 0.2444 0.2640 0.2540 0.00286

Af -0.0159 0.0826 0.0354 0.01520

W, м/с -1.0191 -0.9979 -1.0087 0.00292

V отр AP 0.2500 0.2565 0.2534 0.00106

Am 0.2434 0.2582 0.2518 0.00224

Af -0.0345 0.0466 0.0044 0.01369

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

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

ЛИТЕРАТУРА

1. Глубокая М. Г. Современное состояние вопроса решения проблемы безопасности на этапе взлета // Искусственный интеллект. 2005. №3, с. 270 — 380.

2. Глубокая М. Г. Метод построения процесса взлета по функции эффективной взлетной массы // Ученые записки ЦАГИ. 2009. Т. XL, № 1, с. 82 — 91.

3. АлексеевК. Б., БебенинГ. Г., ЯрошевскийВ. А. Маневрирование космических аппаратов. — М.: Машиностроение, 1970, с. 402 — 407.

4. T a n d a l e M., Bowers R., Valasek J. Robust trajectory controller for vision based probe and autonomous aerial refueling // AIAA-2005-5868, AIAA Guidance, Navigation, and Control Conference and Exhibit San Francisco, CA. — 2005, p. 1 — 16.

Рукопись поступила 14/VII2011 г.