CC BY
64
2007
НАУЧНЫЙ ВЕСТНИК МГТУ ГА Серия Прикладная математика. Информатика
№ 120
УДК 629.735
МОДЕЛИРОВАНИЕ РИСКА СТОЛКНОВЕНИЙ ДЛЯ ВОЗДУШНЫХ СУДОВ, ДВИЖУЩИХСЯ ПО ОДНОМУ МАРШРУТУ.
КИНЕТИЧЕСКИЙ ПОДХОД
Предложена модель расчета риска столкновений для стационарного потока воздушных судов, двигающихся по одному маршруту на одной высоте. В модели такого класса впервые учтена кинетика функции распределения ошибок состояния воздушных судов, учтена ее нестационарность.
Одним из важнейших параметров, определяющих уровень безопасности полетов, является гипотетическая вероятность столкновения в расчете на один час полета воздушного судна (ВС). Рекомендованные ИКАО ограничения на эту величину должны учитываться как при разработке новых ВС и наземных средств обеспечения полетов, так и при планировании воздушного пространства. Актуальность этой тематики подтверждается большим числом публикаций последних лет, посвященных моделированию ситуаций, чреватых возможностью столкновений ВС. Здесь следует выделить цикл публикаций К. Махадхеби, завершившийся документом [1], подготовленным для утверждения в ИКАО. Представленные в этой работе материалы аккумулируют основные достижения в этой области исследований за более чем тридцатилетний период. Вместе с тем, документ отражает (в непереработанном виде) и те устоявшиеся подходы, которые, на наш взгляд, следует модифицировать. Речь здесь, в первую очередь, идет о представлении функции распределения ошибок состояния ВС и о технике ее обработки при вычислении фундаментальной (для задач этого класса) величины - плотности вероятности перехода в состояние катастрофы, порожденной столкновением. Мы не будем здесь останавливаться на критике опубликованных результатов (читатель при желании сам может сопоставить подходы, приведенные ниже и в работе [1]), а лишь представим наше виденье решения проблемы на примере наиболее «туманно» изложенной модели движения ВС по одному маршруту на одной высоте (раздел 5 работы[1]).
Анализ современных моделей расчета риска столкновений показывает, что основные проблемы возникают при попытке реконструкции роли функции распределения ошибок состояния воздушных судов в системе отсчета, связанной с другим ВС.
Действительно, согласно хорошо известной формуле Райса [2,3]
для Л (/)- плотности вероятности перехода системы из состояния «столкновения нет» в состояние «столкновение произошло» (плотности вероятности перехода в состояние катастрофы в момент времени I) неопределенна лишь функция (х, и). При ее задании дальнейшее реше-
ние задачи может вызывать лишь некоторые (связанные с вычислением значений интегралов) технические затруднения, не носящие принципиального характера. В работах [4,5] было показано, что /((х, и) есть функция распределения ошибок в состоянии одного из ВС конфлик-
В.Л. КУЗНЕЦОВ
Введение
1. Идеология построения модели
(1)
тующей пары в системе отсчета, связанной с реальным (вообще говоря, случайным) положением второго ВС, определяемая по формуле
Я (хй) = I^ ^ У1 (^^г)• у; (г1+ху, + й, г), (2)
где / и /2 - функции распределения ошибок, записанные в земной системе координат для первого и второго ВС, соответственно.
Таким образом, соотношение (2) сводит проблему вычисления /(х, й) к более простой задаче - определению функций / (г, V, г) г = 1,2, описывающих отклонение каждого из ВС конфликтующей пары от плановых параметров полета.
В наиболее часто встречающейся постановке задачи моменту времени г , для которого рассчитывается вероятность столкновения, предшествует момент времени г0г < г , в который г-й
экипаж передает информацию о состоянии ВС. Это информация не является абсолютно достоверной, т.к. всякой измерительной процедуре свойственны ошибки. Однако эти ошибки относительно просто оцениваются так, что функцию распределения ошибок р (г -г0г, V -, г0г)
можно считать известной. Здесь г0г и - сообщаемые диспетчеру параметры полета.
Таким образом, задача сводится к следующей постановке: в момент времени г0г известна функция р (г -г0г, V - у0г, г0г); необходимо вычислить функцию / (г, V, г) для момента г > г0г.
Если интервал времени (г - г0г) не слишком велик, и ВС не совершает за это время запланированных маневров, то его скорость V можно считать постоянной, т.е. его движение равномерно и прямолинейно. Такое приближение видится естественным и часто встречающимся в задачах моделирования риска столкновений [6,7].
В этом случае эволюция искомой функции распределения определяется простым кинетическим уравнением [8, 9]:
+ ? М. = 0. (3)
Э г Э ?
Как нетрудно убедиться непосредственной подстановкой, решением уравнения (3), удовлетворяющим начальному условию
/, ( 7 1 г0г ) = Рг (7 - 4 ,17 - У0г , г0г ) (4)
является функция
Я ( 7 1 г0г )=Р (7 - 17 ( г - г0г )-4 ,17 - У0г, 4 ) . (5)
Подставляя (5) в (2) и полученный результат в (1), находим решение ядра задачи о моделировании риска столкновений на малых временных интервалах.
Применим теперь очерченный подход к задаче расчета риска столкновений для ВС, движущихся по одному маршруту на одной высоте.
2. Постановка задачи
Рассмотрим стационарный поток ВС, двигающихся друг за другом по одному маршруту на одной высоте с одинаковой плановой скоростью. Длина маршрута в зоне ответственности диспетчерской службы равна Ь. Для каждой пары соседних ВС вероятность столкновения рассчитывается в соответствии со следующим сценарием. Первоначально в момент времени г0 = 0 соседние ВС с номерами г и (г +1) удалены друг от друга на расстояние г+1 и сообщают диспетчеру свои положение и скорость. Периодичность таких докладов Т. При г = Т ВС вновь сооб-
щают свои координаты и скорости. Если наблюдается недопустимое сближение ВС, диспетчер выдает рекомендации, получаемые на борту ВС в момент времени I = Т+т .
Полагая плотность распределения вероятности ^(я) для расстояния между соседними ВС в потоке ({я г+1} с {я} ) заданной, определить вероятность столкновений ВС на маршруте в расчете на час полетного времени.
3. Предварительные замечания относительно связи БЕ- и гауссового распределений
Изложенные выше соображения относительно алгоритма расчета плотности вероятности Л (^) показывают, что основные трудности связаны с техникой вычисления поверхностных и многократных интегралов в (1),(2).
Для описания отклонений ВС от плановых параметров полета в ИКАО принято так называемое ББ - распределение, имеющее в одномерном случае следующий вид
N
/“ (X) = -^ е1 . (6)
21
Использование напрямую этого выражения оказывается не слишком удобным для конкретных вычислений, поэтому мы воспользуемся здесь следующим искусственным приемом.
Представим ББ - распределение (6) в виде следующего интегрального преобразования:
где £(1 Л)=1е
О 1 ч_^Х_ „ 212 1
Его легко можно получить из формулы (3.325), приведенной в [10].
Соотношение (7) можно переписать также в виде
/ш (х)=Ь /о (х), (7а)
го, л 1 - Х21Ъ{Х0 )2 ~
где / (х) = ,—-е ' ' - гауссово распределение, а Ь - интегральный оператор с ядром
42к 1 х
£ (Л° ,1). С другой стороны, из (1), (2) нетрудно видеть, что Л (^) может быть также представлен в компактном виде
Л(0 = ¿{./¡“/Г), (8)
где (О - некоторый интегральный оператор.
Учитывая коммутируемость операторов й и Ь из (7 а) и (8), получаем:
Л(0=О ■ {— /?■— )=Ц-й ■ {/?-/?) . (9)
Соотношению (9) можно придать следующую интерпретацию. Плотность вероятности перехода при гипотезе о справедливости ББ - распределения (обозначаемая далее как Лш ^)) выражается с помощью интегральных операторов — Ь2 через величину Ло (^), соответствующую гипотезе о справедливости гауссовости распределения ошибок. Таким образом, имеем:
Лш (Г) = ЦЬ2 Ло (I). (9а)
х2
Достоинство соотношения (9а) заключается в том, что для случая гауссовых распределений в аналитике можно продвинуться значительно дальше.
Приведем здесь также связь функций распределения относительных ошибок для ББ- и га-уссового распределений
/“(г,у;!4 1<2)=//.../< <1Л1 ш(<, <) ш(< <<2))*/° (Г,у;Х1 <) (10)
0 г=1
Здесь введены следующие обозначения 1(1) =(<(1), 1(21), <:(1)), < = (<п, <12, <) - векторы параметров ББ - и гауссового распределения, <(1) ° ЛЯ, <(1) ° ° ° < и т.д.
4. Оценка плотности вероятности перехода в состояние катастрофы при гипотезе гауссовости начальных распределений
Рассмотрим два ВС, двигающихся по одному маршруту на одной высоте с дистанцией я . В предположении гауссовости распределений ошибок в положении каждого из них можно показать, что неопределенность состояния одного из них в системе отсчета, связанной со вторым, имеет вид:
1
26я5П(0,-Л?)
П єхр
1
4Х2
и --
(1 )2
Л2'
л к и
)
(X + и)2
4
(її)
Здесь {я} = {я, 0,0},і = 1,2,3, ст - три компоненты гауссового распределения ошибок ВС по
Е 1°-м па 42
скорости,
4(1 )2+М2
, 1° =^-41 ° +1 ° , а 1°і и 1° - параметры гауссовых распре-
2
делений ошибок в положении первого и второго ВС. Здесь мы полагаем, что параметры распределения для первого и второго ВС, вообще говоря, различны. Это требование вытекает из уравнения (9), где интегральные операторы Ь1 Ь2 действуют на функции распределения разных ВС.
Теперь с учетом (1) можно получить выражение для Ле (^)- плотности вероятности перехода в состояние катастрофы при гипотезе о гауссовости распределения ошибок. Полагая, что критическая область может быть смоделирована прямоугольным параллелепипедом, и вводя обозначение
Г 12£2,21
Ф(,; 1'11 = ЄХРГ1+4Я' 21»
1-X-, ( „( 1-X-,Л Л
ег/
2Я2
-1
(12)
находим:
Л° (,1° ,а)
І2 /
X-,
Л Ь3 (
йх21 ехр
я
V
+Е2Ф0; ^ Е2,12)АеХР“
Х3Ф(,, Ь3, Х3, Л3 )Ь1 ехр <
ь2
41 + М2)
%
4(1 + М2)
ь3 (
• -1ехр
0
^2 (
• -1ехр -
о V
4(1 +м2,2) J о
•2 Л
X
2 Л
4(132 +м32,2)
йхъ +
X,
4(12 +м32,2)
ёх3 +
х,
2
4(1 + м2,2)
г =1
0
г =1
2
я
Выражение (13) позволяет определить плотность вероятности перехода в состояние катастрофы, вызванной столкновением, как функцию времени I - продолжительности совместного полета ВС в предположении того, что в момент времени t0 = 0 столкновения заведомо не было.
5. Оценка плотности вероятности перехода в состояние катастрофы при гипотезе
о справедливости БЕ- распределения
Перейдем теперь к вычислению Лш (ї)- плотности вероятности перехода в состояние катастрофы, для чего воспользуемся соотношением (9а) и результатами, полученными в п.4. Произведение Ь1 •Ь2 в (9а) представляет собой шестикратный интегральный оператор, который можно представить в виде произведения трех двукратных операторов вида
(14)
0 0
Учитывая, что согласно (13) 1 и 1 фигурируют в выражении для Ла (ї) только в виде комбинации [(Да )2 + (1 )2], в (14) можно сделать замену переменных, соответствующую переходу в полярную систему координат:
1 =1а С0*Ф, , 1 =1а С08 р , І = 1,3
Поскольку Ла (ї) не зависит от величины (р, , то операцию интегрирования по этой переменной можно легко выполнить, понизив при этом в 2 раза общее число интегралов. Окончательно, для величины Лш (ї) получаем следующее выражение
Л “ (ї) = 2 Л,“ (ї),
І=1
где
Л.“ (ї) = 4 1
X
II 2 Ь2
2 е 122 • {^е
2
0
X
У
2 -— ¿з
2
¿1
Р &1&2&3
2 т2
2 —2 ^
52 Л
2 —1___________________________________
{йт 1-• е112 •е4—2+сті2ї2)
0 1
X
У
(15)
(16)
{ й— ^ Є 122 Ф(ї, ¿2, Ї2,—2 ) • Є 4(122 +"22ї2)
X
У
2 -— ¿3
2
12
Л3°Б (ї) = ~2
А
Р &1&2&3
а2 \
\
2 —1___________________
{йт • їг —-• е112 •е 4—2+Сті¥)
•' Л
X
У
:2 Л Л
0
X
{<— Ї32- 1-Ф(!,¿3,Ї3, —3)' е
ь32 Л
4(—32 +СТ32^)
3 13
(17)
4(—1 +01 ї )
е
(
4(—2 +СТ2 ї )
4(—3 +СТ3 ї )
е
0
0
0
Л
¿22 Л (
4— +^3 ї )
е
0
0
с
2
4(—2 +^2 ї )
0
0
Здесь Ф(^, I, Х,1) определено соотношением (12).
Заметим, что полученное в (15)-(18) выражение для Лш(^) зависит от ряда параметров:
^ - известного диспетчерской службе из докладов расстояния между ВС и Л,д — векторных параметров DE-распределения, определяющих стандартные отклонения значений координат и скоростей ВС от докладываемых. Поэтому для корректности вместо Лш (^) будем использовать
обозначение Лш ^ $,Л,&). На рис.1 представлены некоторые графики для Лш ^ $,Л,&), иллюстрирующие эту параметрическую зависимость.
Рисунок. Графики зависимостей плотности вероятности столкновения пары ВС - Лш ^ я,Л,&) при следующих значениях варьируемых параметров задачи: а) Л1(^)- при ^ = 30ёг , 1. = 25г ; Ь) Л2(7)- при ^ = 30ёг , 1 = 20г ; с)Л3(^)- при ^ = 29,5ёг , 1 = 25г
6. Оценка риска столкновений для всего маршрута в зоне ответственности диспетчерской службы
Пусть протяженность рассматриваемого маршрута в зоне ответственности диспетчерской службы равна Ь, поток ВС, движущихся по маршруту, стационарен. Очевидно, что в этом предположении вероятность столкновения ВС за час полетного времени в зоне ответственности Ь равна вероятности столкновения, рассчитанной для некоторого «цуга» ВС (уединенной группы ВС численностью К, двигающихся друг за другом по одному маршруту на одной высоте и с одной скоростью) в течение одного часа. Рассчитаем эту вероятность.
Нетрудно видеть, что вероятность отсутствия столкновений в цуге за время ^ - (1 — Р"^ ^))
равна произведению вероятностей отсутствия столкновений между любой из пар цуга, поскольку именно эти события являются независимыми
(1 - Рш(I)) = П(1 -('■ ^))
(19)
Из формул (15)-(18) видно, что плотность вероятности столкновений пары ВС катастрофически быстро спадает с увеличением 5, так что Лш(I 5,1,0)>>Лш(I 25,1,0’) . Поэтому в (19) можно оставить только пары ближайших ВС, т.е.
N—1
(1—Р1 (I)) = П(1 — Л»5«.1))»1 — 2рТш(I. 5,,,+1) (20)
Полагая число ВС в цуге достаточно большим N ■
>> 1, где (s) = I s • w(s)ds
среднее расстояние между ВС, выражение для вероятности столкновения в цуге - Р"^ (I) можно
записать в следующем виде
: А+i']xiр?ш(t,s)•w(s)ds = L+1']x\dt'\ds•w(s)s,l,s). (21)
p..® (t)=
ooa V '
Поскольку ограничения накладываются на величину риска столкновений, численно равного вероятности столкновения за час полетного времени в расчете на одно ВС, а информация у диспетчеров обновляется по истечении интервала времени Т, то выражение для ожидаемого числа столкновений в рассмотренной модели может быть записано в виде следующего соотношения
T T ^
Nax = T | dt\dS W(S) •L(t1S,1,S)
(22)
Здесь Т0 - 1 час полетного времени или 3600 секунд. Это необходимо учитывать, поскольку на рисунке единичный интервал времени равен одной секунде.
Заключение
2=1
О
О
Предложена новая модель расчета столкновений ВС, двигающихся по одному маршруту на одной высоте. Основное отличие развиваемого подхода от известных заключается в привлечении к решению задачи аппарата кинетической теории, что впервые позволило отказаться от приближения независимости позиционных и скоростных ошибок при описании состояния ВС, а также корректно учесть динамику вероятностных распределений этих величин.
Показана связь между плотностями вероятности перехода в состояние катастрофы, вызванной столкновением ВС при гипотезах о справедливости гауссового и DE-распределений. Полученное выражение для риска столкновений отличается от известных ранее [1],[11].
В заключение автор выражает благодарность доктору технических наук В.Б. Спрыскову за полезные обсуждения.
ЛИТЕРАТУРА
1. Mehadhebi K. A unified framework for collision risk modeling in the airspace planning manual document with further applications. ICAO SASP-WG/WHL/9-WP/12, attachmenta.
2. Reich, P.G. Analysis of long-range air traffic systems - separation standards. Part I, II, III // The Journal of Institute of Navigation, 1966, vol.19, p.88,169,331.
3. Hsu, D.A. The evaluation of aircraft collision probabilities at intersecting air routes // The Journal of Navigation, 1981, vol.34, N1, p.78.
4. Kuznetsov V.L. Markov collision risk model. ICAO SASP-W6/WHL7/-IP/2, Monreal, May, 2005.
5. Кузнецов В.Л. Марковская модель оценки риска катастроф на воздушном транспорте // Научный Вестник МГТУ ГА, серия Эксплуатация воздушного транспорта и ремонт авиационной техники. Безопасность полетов, №91, 2005.
6. Anderson D. A longitudinal distance-based collision risk model based on reliability theory installment 1 - crossing angles less then 90 degrees. SASP-WG/WHL/3-WP/3.
7. Fujita M., Nagaoka S., Amai O. Safety assessment prior to implementation of 50NM longitudinal separation minimum in R220 and R580. SASP-WG/WHL/9-WP/14.
8. Кузнецов В.Л. Кинетический подход к описанию эволюции неопределенности состояния воздушного судна в задаче расчета риска катастроф // Научный Вестник МГТУГА, серия Прикладная математика. Информатика, № 105, 2006.
9. Кузнецов В.Л., Соломенцев В.В. К задаче моделирования риска столкновений воздушных судов / Статья в настоящем Научном Вестнике.
10. Грандштейн И.С., Рыжик И.М. Таблицы интегралов, сумм, рядов и произведений. М.: Наука, 1971.
11. Руководство по методике планирования воздушного пространства для определения минимумов эшелонирования. ICAO Doc 9689.- Монреаль, ИКАО, 1998.
COLLISION RISC MODELING FOR AC MOVING ON SAME RATES - KINETIKAL APPROACH
Kuznetsov V.L.
A collision risk model for a stationary aircraft flow, moving on same rates, is suggested. A kinetic of distribution of aircraft position in models of such class is take account for the first time.
Сведения об авторе
Кузнецов Валерий Леонидович, 1949 г.р., окончил МГУ им. Ломоносова (1972), доктор технических наук, заведующий кафедрой прикладной математики МГТУ ГА, автор более 90 научных работ, область научных интересов - методы математического моделирования в задачах распространения излучения в пространственно неоднородных случайных и периодических средах, безопасность полетов.
