CC BY
50
УЧЕНЫЕ ЗАПИСКИ ЦАГИ
Том XXIV 1993 № 3
УДК 629.735.33.015.073 629.735.33.015.077
ОЦЕНКА ЭКСТРЕМАЛЬНЫХ ОТКЛОНЕНИЙ ПАРАМЕТРОВ ПОЛЕТА ЛЕТАТЕЛЬНЫХ АППАРАТОВ*
В. П. Кузьмин, В. А. Ярошевский
Рассматриваются задачи, в которых успешность полета определяется выполнением ограничений на параметры траектории в фиксированный момент либо в течение всего времени полета. Приведены результаты исследований по разработке методики оценки вероятности успешного полета летательного аппарата при действии случайных возмущений. Обсуждаются проблемы численного решения различных задач. Данные, полученные по предлагаемой методике, сравниваются с другими аналитическими и численными оценками для модельных примеров. Приводятся результаты решения прикладных задач.
Задачи об оценке экстремальных отклонений параметров полета от номинальных значений приобрели большую актуальность за последнее десятилетие в связи с внедрением вероятностных норм в практику авиационной науки и техники.
Рассматриваются два типа задач. Примером задачи первого типа («локальная» задача) является анализ процесса автоматической посадки. По существующим нормам [1] основное требование, предъявляемое к системе управления посадкой неманевренного самолета, сводится к тому, чтобы вероятность серьезного летного происшествия при автоматической посадке не превышала 10-6 -10-8. Возможность возникновения летного происшествия зависит в первую очередь от значений параметров полета (вертикальная скорость Уу, дальность Ь, боковое смещение 2 от оси полосы и т. д.) в момент касания самолетом взлетно-посадочной полосы (Н= 0). Поэтому очень важно оценить вероятность выхода этих параметров за допустимые пределы, например, оценить значения
Р[Уу(Н = 0)<¥ут[а]. Р[ЦН = 0) < Ьт1п]
•Доклад, прочитанный на VII Всесоюзном съезде по теоретической и прикладной механике. — Москва, 20 августа 1991 г.
Ко второму типу (задачи о выбросах) относятся задачи об оценке вероятности выхода кинематических параметров полета за допустимые пределы на определенном конечном интервале времени. Например, при полете самолета в зоне турбулентной атмосферы на интервале времени | г0, ] требуется оценить вероятность выхода угла атаки самолета а за
допустимые пределы, т. е. оценить значение
Р\а(т)>
атах Ь Т е[(0>(к 1-
Особенностью задач обоих типов являются малые значения искомых вероятностей, что препятствует применению при математическом моделировании широко распространенного метода Монте-Карло: число потребных реализаций N оказывается непомерно большим, достигает миллионов и даже миллиардов.
Линейные оценки, основанные на гипотезе о гауссовском распределении выходных параметров, также, как правило, оказываются несостоятельными, поскольку в экстремальных условиях, например при воздействии больших ветровых возмущений, органы управления летательными аппаратами выходят на упоры, либо скорости их перемещения выходят на ограничения, к тому же часть возмущений имеет негауссовское распределение, особенно в области больших уклонений от математического ожидания.
В то же время малость искомых вероятностей облегчает применение асимптотических оценок для определения экстремальных отклонений параметров полета.
Наряду с вероятностными подходами к решению поставленных задач на практике используются и «гарантированные» подходы, когда для имеющихся случайных возмущений (например, заданных в виде совокупности случайных параметров) задаются предельные значения и в качестве расчетных случаев рассматриваются наиболее неблагоприятные сочетания предельных возмущений с соответствующими знаками [2]. Не отрицая в принципе гарантированного подхода, авторы статьи считают все же более обоснованным и отвечающим природе вещей вероятностный подход.
При разработке приближенной методики решения поставленных задач были использованы результаты, нашедшие свое отражение в ряде монографий [3 — 10].
При постановке задачи первого типа в конечномерном случае считается, что случайные возмущения сведены к конечному числу независимых гауссовских нормированных случайных параметров с'; (/ = !,/?):
Рис. 1. К оценке вероятности в конечномерном случае
Выходная переменная у является непрерывной дифференцируемой функцией от с,-, вероятность превышения переменной у значения у„ определяется интегралом
Р[У>У*] = \р(с)с1с1.. .с!с„, (1)
£2
взятым по области П, расположенной за поверхностью уровня у = у*. Если эта вероятность достаточно мала и поверхность У = У* расположена достаточно далеко от начала координат, то основной вклад в величину инте1рала вносит окрестность точки поверхности с», ближайшей к началу координат, поскольку в области О функция р(с) быстро убывает по мере удаления от точки с,. Тогда естественно предположить,
что при замене поверхности уровня У = У* плоскостью, касательной к поверхности в точке с, (рис. 1), можно грубо оценить в первом приближении искомую вероятность
где й= £<£
и=1 )
Экспоненциальный множитель согласуется с оценками, строго обоснованными в [4].
Следовательно, для оценки вероятности в первом приближении (Pj) достаточно решить задачу на условный экстремум:
П
R2 = min V с? при условии у(с) = у+ (3)
с‘ £г
П
или эквивалентную задачу у» = max у(с) при условии V с} = R2, где R
С/ « является функцией заданной вероятности Р: Ф(Л) = Р.
Решение такой задачи может быть основано на сочетании стохастических и градиентных методов поиска, причем их трудоемкость существенно возрастает с увеличением п — размерности пространства случайных параметров [11].
Для уточнения этой оценки в случае нелинейной функции у(с) следует учесть кривизну поверхности уровня, причем наиболее существенна форма поверхности уровня в окрестности точки с,. Поворачивая
систему осей С1С„ так, чтобы ось q была направлена в точку с*, запишем уравнение поверхности уровня в виде
2Л1=2*=2 где аи > -1.
Тогда уточненная оценка искомой вероятности определяется формулой
-И? > Д’.] * Рг =------------, (4)
2 [<1е1(Е + А)]1/2 ’ ^
где Е — единичная матрица, А — матрица, составленная из элементов Щк [12].
Оценка (4) может привести к значительным по1решностям, если некоторые из значений = -1 (вырожденный случай).
Другой метод уточнения оценки (2), приемлемый и в случаях, близких к вырожденным, связан с применением специальной модификации метода Монте-Карло — метода существенной выборки при использовании сферических координат, с помощью которого удается приближенно вычислить многомерный интеграл (1) [12, 13].
Применение описанной методики иллюстрируется результатами расчетов вертикальной скорости Уу в момент касания самолетом взлетно-посадочной полосы (рис. 2) с несколько упрощенной моделью возмущений, характеризуемой восемью случайными параметрами. Здесь
Рис. 2. Сравнение оценок вероятности первого приближения (штриховая линия) и второго приближения (сплошная линия) с результатами численного интегрирования прочности вероятности (+)
штриховая линия соответствует оценкам вероятности Р[Уу <Уу*] по первому приближению (2), пестиками обозначены результаты расчета по второму приближению (4), сплошная линия — результаты расчетов с применением методики вычисления многомерных интегралов в сферических координатах. Как видно, крестики располагаются близко к сплошной кривой.
В полной модели случайных возмущений, действующих на самолет в процессе автоматической посадки, число независимых случайных параметров следовало бы увеличить до сотни и более, при этом трудоемкость вычислений катастрофически возрастает.
Для преодоления этой трудности предложено разделить всю совокупность возмущений на «главные» (не более 10 — 15) и «второстепенные» возмущения. Такое разделение оказывается возможным во многих практических случаях. Последовательность проведения расчетов сводится к следующему [13]. Сначала учитываются только главные возмущения и решается задача (2) на условный экстремум, что дает возможность
наити точку с», значение К, а также определить градиент в
указанной точке. Далее фиксируется найденная точка с* («наихудшая» совокупность главных возмущений) и проводится статистическое моделирование процесса с варьированием второстепенных случайных возмущений в пределах 100 — 200 реализаций, что позволяет найти дисперсию Д*, выходной координаты у, обусловленную влиянием второстепенных возмущений. После этого найденное значение К корректируется по формуле
(5)
1 +
2
при этом искомая вероятность с учетом всей совокупности возмущений определяется соотношением
В итоге такой двухэтапной процедуры расчетов удается ограничиться суммарным числом реализаций, не превосходящим нескольких сотен или тысячи, для сколь угодно малых значений вероятности.
Применение описанной методики для расчетов процесса автоматической посадки неманевренного самолета с использованием полной модели возмущений проводилось в сотрудничестве с Московским институтом автоматики и электромеханики (В. В. Смирнов, О. Б. Кербер, В. Н. Мазур). В качестве главных возмущений рассматривались две составляющие систематического горизонтально направленного ветра их
и и1, а также атмосферная турбулентность (продольные и боковые порывы ветра IVх и Ж?). Для того чтобы описать случайные функции, определяющие атмосферную турбулентность, конечным числом случайных параметров, эти функции заменялись каноническими разложениями [7], например
где йк — нормированные независимые гауссовские случайные параметры; <рк (г) — детерминированные координатные функции. Каноническое разложение гауссовской функции может быть выполнено бесконечным
Р[У>У.]« Ф(Ю-
(6)
П
о
М 2Аг
t
Рис. 3. Координатные функции канонического разложения порывов ветра
числом способов. В частности, для случая, когда корреляционная функция для переменной И/’х(() является экспонентой, координатные функции изображены на рис. 3.
Некоторые результаты расчетов изображены на рис. 4 и 5. На рис. 4 приведено изменение по дальности высоты нижней кромки шасси #ш, угла тангажа 9, вертикальной скорости Уу и возмущения Н\. на экстремальной траектории, соответствующей £тах - 900 м. Штриховым линиям на рис. 4 соответствует изменение параметров на траектории без возмущений (1¥х = 0). На рис. 5 приведены относительные значения предельных отклонений различных параметров траектории самолета _ х* - х
х* ~ ~ в сравнении с гипотетическими законами распределения — гау>
ах
совским (сплошная линия) и законом распределения порывов ветра (штриховая линия). На рис. 5 использованы обозначения: у — угол крена, у/ — угол курса, I — боковое смещение от осевой линии ВПП.
Из приведенных результатов следует, что нелинейность задачи приводит к невозможности использования для оценок предельных отклонений какого-либо одного распределения для всех координат. Отметим, что, как следует из результатов работы [181, для линейной задачи оценки предельных отклонений должны лежать между оценками, соответствующими двум гипотетическим распределениям, т. е. лежать между сплошной и штриховой линиями на рис. 5.
В бесконечномерном случае характерная задача первого типа сводится к анализу решений стохастического уравнения, рассматриваемого в смысле Стратоновича [3]:
— = /(*,/) +С(*,г)£(0, ...
<и (7)
< + Д^>= Е6(А1),
где £(/) — вектор независимых белых шумов; /(х,/) и <7(х,/) — векторная и матричная функции; х(/) — вектор фазовых координат.
Если начальное значение х(/0) = х0 фиксированно и требуется оценить вероятность Р\хх^к) > х,] * Рх = Ф(/?), где х, — первый компонент вектора х, то в первом приближении, по аналогии с конечномерным случаем, задача сводится к эквивалентной детерминированной вариационной задаче [4, 14] по определению вектора £(/), доставляющего
минимум значению
т !к » = 1 /0
Щ'.м/с
5,0 - -
2,5 Л I \ 1 1
и 0
-2,5 - - \
-5 П 1 1 1 1 >/ 1
3, и -500 -250 С 250 500 750
1*,У„
♦
2 О -2
* -500 -250 О 250 500 7501,м
Рис. 4. Изменение по дальности параметров экстремальной траектории посадки, соответствующей
Рис. 5ч Сравнение вероятностей превышения характерных значений различных параметров траектории в момент приземления самолета с нормальным распределением (сплошная линия) и распределением порывов ветра (штриховая линия)
при условии х(/0) = х0, х1(^) = х*, где т размерность вектора белого шума.
Тогда
Р[х,(/*)>х.]«. Л = *(*)• (8>
Если начальные условия для вектора х0 не фиксированны, а имеют гауссовское распределение, то в выражение для Я2 добавляется квадратичная форма от значений Д х,0, характеризующих отклонения компонентов вектора х0 от их математических ожиданий.
В ряде случаев удается решить эквивалентную вариационную задачу аналитически; например, для скалярного уравнения
— + /(*) = В €(*) (9)
<а
при /0 = - оо, /(0) = 0, х, > 0 оказывается, что
g 4 = 2/(х) при /(х) > О, g4=0 при Ддс)<0;
для уравнения
^ + 2*^ + /(*) = *{М (К»
«г»2 ви при Ч = -°°, к> О, /(0) = 0, — (х) > 0
<1 X
оказывается, что
е 4 1 Лх
41
По аналогии с конечномерным случаем можно уточнить оценку искомой вероятности, если рассмотреть вторую вариацию функционала и воспользоваться методом соседних экстремалей [15], что эквивалентно формуле (4) для конечномерного случая:
Р2 =
[с!е1(Е + А)]1/2
р
о
5
Я
-----РГР1!Р0----------?2=Р2/РЙ
Рис. 6. Сравнение оценок вероятности первого (сплошные линии) и второго (штриховые линии) при-
ближений
где матрица А определяется в результате интегрирования некоторой системы обыкновенных дифференциальных уравнений после решения вариационной задачи [14].
На ряде модельных примеров удается сопоставить точные результаты с оценкой искомой вероятности по первому и второму приближению. Например, на рис. 6 приведены результаты такого сопоставления для случая, когда в уравнении (9) /(х) является линейной функцией с
Задачи второго типа, связанные с оценкой вероятности выбросов случайных функций на определенном интервале времени, оказываются значительно сложнее задач первого типа. Здесь удалось получить необходимые оценки только для случаев, когда динамическая система, возбуждаемая белым шумом, имеет первый или второй порядок.
Для динамической системы первого порядка вида (9) искомая вероятность Р\х{г) > хД т е [0,*] определяется на основе решения соответствующего уравнения Фоккера — Планка — Колмогорова [3] методом Фурье, для чего необходимо определить спектр собственных значений цк одномерной краевой задачи:
При достаточно больших значениях х, и / эта вероятность определяется простейшей формулой:
насыщением: /(х) =
і*[х( т) > х* ] = 1 - Ьке~^к{.
Р[х(т) > х*] * 1 - е~М
(П)
практически независимо от начального распределения переменной х, а значение //0 определяется соотношением [3]:
-1
(12)
где р„(х) — стационарное распределение переменной х, которое нетрудно найти как стационарное решение уравнения Фоккера — Планка — Колмогорова. Данный результат относится к простейшему примеру недифференцируемой случайной функции х(0у дисперсия производной которой бесконечна (егх =<х>).
В случае дифференцируемой случайной функции х(() (ах <ао), которая характерна для систем более высокого порядка, возбуждаемых белым шумом, например вида (10), при достаточно больших х* и / искомая вероятность по-прежнему описывается приближенной формулой (11), где значение Ио является наименьшим собственным значением многомерной краевой задачи, для которого не удается получить простых аналитических выражений. Для гауссовских процессов доказано [16], что при х* -> оо число выбросов на заданном интервале времени приближается к И— числу пересечений функций х(0 уровня х* изнутри в единицу времени
Однако в случае больших, но конечных значений х, оценка (13) может оказаться существенно завышенной, поскольку она эквивалентна предположению о том, что последующие выбросы не зависят от предыдущих.
Сказанное иллюстрируется результатами численных экспериментов, проведенных для динамической системы:
при различных к > 0, и > 0, в ходе которых определялась вероятность
выбросов Р[|х(т)| > х, ], г е [0,?] на значительных интервалах времени Л
При п = 1 (линейная система), = 1, ах = 1 наибольшее количество выбросов наблюдается при значениях к, несколько меньших еди-
(13)
о
—■ + 2к^+ \х\п &щпх = 2 у[к £(/)
(14)
оценок параметра для различных значений к
ницы, а при больших и малых значениях к частота выбросов уменьшается (рис. 7), хотя оценка (13) в данном случае дает один и тот же результат, не зависящий от к:
МО »Н=-е-^/\ (15)
Л
где Я = х, (в данном случае следует удвоить оценку (12), поскольку фиксируются пересечения двух уровней ± Я). Причина заключается в том, что при больших к переменная х является выходом цепи из двух апериодических звеньев с сильно различающимися постоянными времени (2к и 1/2к), на вход которой подается белый шум. Тогда роль звена с малой постоянной времени оказывается второстепенной, переменная х(0 становится «близкой к недифференцируемой», пересечения функцией х($ любого из уровней происходят «с дрожью», т. е. сериями (рис. 8, а). В другом предельном случае, при очень малых к, пересечения уровней функцией х(0 также происходят сериями, поскольку процесс х(0 является колебательным, слабодемпфированным, и переменная хф, достигнув однажды уровня х„, во многих случаях повторяет эти пересечения (рис. 8, б, в).
X ♦
2 О -2 -Н *
2 О -2
*
2 О -2
Рис. 8. Характерные реализации случайного процесса д ля различных
значений параметра к
В результате решения уравнения (14) при и =1 в случае больших и малых к можно аппроксимировать решениями уравнения первого по-
йгх
рядка: в первом случае для этого достаточно отбросить член в
Л2
уравнении (14), а во втором — перейти к эквивалентному уравнению для изменения амплитуды колебаний согласно методике [3]. Тогда, применяя формулу, аналогичную (12), можно получить оценки
Ао * *е-*2/2 (16)
к4 2 я
при больших к,
Но * кЮе~*гП <17)
при малых к.
Результаты расчетов по этим формулам также приведены на рис. 7 (штриховыми и штрихпунктирными линиями соответственно). Как видно, совокупность формул (14), (15) и (16) с достаточной полнотой
отражает «точные» результаты, полученные в численном эксперименте.
Рис. 9. Сравнение численных и анали- Рис. 10. Сравнение численных и ана-
тических оценок параметра р0 для литических оценок параметра
различных к и п = 0 для различных к и п = 3
По мере роста К, как это следует из теоремы, доказанной в [16], пределы применимости формулы (15) расширяются.
Аналогичные численные эксперименты были проведены для случаев, когда в уравнении (14) принималось п — 0 и п = 3; соответственно для этих случаев были получены формулы, аналогичные формулам (15), (16) и (17). Результаты сопоставления приведены на рис. 9 и 10.
Полученные результаты позволяют в ряде случаев с достаточной точностью оценить вероятность выбросов угла атаки летательного аппарата с нелинейными аэродинамическими характеристиками при полете в турбулентной атмосфере [17].
ЛИТЕРАТУРА
1.Air worthiness technical manual, chapter 7, section 8, part 3.
2. P о w e 11 R. W., Stone N. W. Analysis of the space shuttle oibiter entry dynamics from Mach 10 to Mach 2,5 with the November 1976 Flight control system//NASA TP — 1667. — 1980.
3. Стратонович P. Л. Избранные вопросы теории флуктуации в радиотехнике//М.: Советское радио. — 1961.
4. Вентцель А. Д., Фрейдлин М. И. Флуктуации в динамических системах под действием малых возмущений//М.: Наука. — 1979.
5.Диментберг М. Ф. Нелинейные стохастические задачи механических колебаний//М.: Наука. — 1980.
б.Черноусько Ф. Л., Колмановский В. Б. Оптимальное управление при случайных возмущениях//М.: Наука. — 1978.
7.Пугачев B.C. Теория случайных функций и ее применение к задачам автоматического управдения//М.: Физматгиз. — 1962.
8. Т и х о н о в В. И. Выбросы случайных процессов//М.: Наука. — 1970.
9. Казаков И. Е., Мальчиков С. В. Анализ стохастических систем в пространстве состояний//М.: Наука. — 1983.
10. Малышев В. В., Кибзун А. И. Анализ и синтез высокоточного управления летательными аппаратами//М.: Машиностроение. — 1987.
11. Кузь мин В. П., Ярошевский В. А. Оценка предельных отклонений параметров траектории самолета при автоматическом заходе на по-садку//Ученые записки ЦАГИ. — 1985. Т. 15, № 2.
12. Кузьмин В. П.,Ярошевский В. А. Приближенный метод оценки вероятности больших отклонений параметров траектории самолета при действии случайных возмущений/Дехническая кибернетика. — 1989, № 2.
13. Кузьмин В. П. Оценка вероятности безопасного полета летательного аппарата при действии возмущений//Ученые записки ЦАГИ. — 1990. Т. 21, № 1.
14. Кузь мин В. П., Ярошевский В. А. Асимптотические оценки вероятностей больших случайных отклонений фазовых координат динамических систем от средних значений//Ученые записки ЦАГИ. — 1990. Т. 21, № 3.
15. Брайсон А., Хо Ю-Ши. Прикладная теория оптимального управ-
ления//М.: Мир. — 1972. .
16. Крамер Г., Лидбеттер М. Стационарные случайные процессы//М.: Мир. — 1969.
17. Н е м ы ки н А. С., Я р о ш е в с к и й В. А. Оценка амплитуды колебаний угла атаки летательного аппарата с нелинейными характеристиками демпфирования при наличии атмосферной турбулентности//Ученые записки ЦАГИ. — 1984. Т. 15, № 3.
18.Кузьмин В. П., Парышева Г. В. Определение статистических характеристик движения самолета при автоматическом заходе на посадку// Ученые записки ЦАГИ. — 1985. Т. 16, № 6.
Рукопись поступила 12 /XII 1991г.
