Breeze написал(а):
- Rob, я с нетерпением жду Вашу модель боевого столкновения F-22 c любым самолётом на Ваш выбор
Начнем...
Ниже предложены приемы составления дифференциальных уравнений, описывающих процессы военных действий, позволяющие выполнить смешанное аналитико-имитационное моделирование исследуемого процесса. В качестве практических примеров разработана модель процесса ВБ двух истребителей (специально для
Breeze) и процесса изменения численности ЛА на ТВД (если не надоест).
Немного теории (для тех, кому действительно интересно).
Дифференциальные модели являются одним из способов аналитического описания процессов военных действий. Они позволяют, зная сходные состояния объектов и интенсивности переходов, определить состояния объектов в заданный момент времени и, таким образом, исследовать динамику процесса.
Основным достоинством дифференциальных моделей (как и прочих аналитиче-
ских моделей), в случае их адекватности, удовлетворительной для исследования, является простота записи (по сравнению с имитационными моделями) и малое время счета. Благодаря этому дифференциальные модели нашли широкое использование в исследовании различных процессов, начиная от описания функционирования отдельных устройств, до боевых действий группировок войск, включающих в свой состав большое количество применяемых сил и средств. В последнем случае, как правило, используются модели Ланчестера.
Решение дифференциальных уравнений в простейших моделях осуществляется анали-
тически, а в случае сложных моделей типа упомянутых "экспресс-моделей" – численно. При этом, в соответствии с теорией дифференциального исчисления, при выполнении дифференцирования не допускаются разрывы функций и скачкообразные изменения параметров. Однако реальные исследуемые процессы могут быть разрывными, иметь скачкообразное изменение параметров и оцениваемых показателей, в том числе стохастическое. Учет в аналитической модели этих особенностей процесса позволил бы повысить ее адекватность при сохранении малого времени счета.
Введение названных уточнений в структуру дифференциальной модели возможно
при численном методе ее решения благодаря тому, что решение дифференциальных
уравнений численными методами по своем характеру соответствует имитационном
моделированию – прежде всего, выполняется пошагово во времени. В связи с этим
уточнение модели заключается в том, что на дискретную основу численного решения
дифференциальных уравнений следует наложить логические и стохастические функции, детализирующие описание исследуемого процесса. В результате может быть получена желаемая гибридная модель, сочетающая достоинства аналитических и имитационных моделей.
Далее предлагаются специальные приемы, позволяющие с позиций имитационного
моделирования провести указанную корректировку дифференциальных уравнений.
1.Приемы численного дифференцирования с использованием элементов имитации.
Будем исходить из того, что численное решение дифференциальных уравнений выполняется одношаговыми методами. Сущность этих методов состоит в пошаговом решении дифференциального уравнения dу(t)/dt = f(t, у(t)) при начальном условии у(t0) = у0. Вместо исходного дифференциального уравнения в точках tj = t0 + j•∆t (j = 0, 1, …) решается разностное уравнение
уj+1 = уj + ∆t•ϕ(tj, уj). (1)
Функция ϕ задается формулами типа Рунге-Кутта, связывающими ∆t, t , у с позиций минимизации погрешности метода вычислений. В простейшем случае при ϕ = f имеет место метод Эйлера. Для повышения адекватности модели исследуемому процессу при пошаговом решении разностного уравнения (1) необходимо следовать определенным рекомендациям и приемам.
1.1.Выбор шага вычислений.
Величина ∆t, с одной стороны, определяет точность численного решения дифференциального уравнения. Так, при вычислении по методу Эйлера, погрешность вычислений в j-м узле составляет ∆t2•f′′(ξj)/2, tj-1 ≤ ξj ≤ tj. При использовании метода Рунге-Кутта точность может быть существенно выше. В случае сложных моделей, включающих совокупности сцепленных друг с другом по входным и
выходным переменным систем дифференциальных уравнений, при малом шаге вы-
числений время счета может составить десятки минут и часы. Следовательно, выбор шага вычислений является результатом компромисса между точностью и временем вычислений.
Кроме того, моделируемый процесс, как в приведенном ниже примере, может иметь определенную дискретную структуру, с которой также следует согласовывать значение ∆t. Значение ∆t, как показано ниже, влияет и на структуру моделируемых стохастических процессов, что также необходимо учитывать.
1.2. Приемы моделирования детерминированных скачкообразных изменений
процесса.
Будем различать случай скачкообразного изменения процесса в заданный момент (моменты) времени и не временной случай скачкообразного изменения процесса
– в зависимости от достигнутого значения параметра, учитываемого в дифференцируемой функции. При необходимости описать скачкообразное изменение процесса в заданный момент времени условие скачка дифференцируемой функций f(t, у(t)) задается с помощью одной из логических функций, реализованных в среде моделирования, например, функции if [А, х1, х2], которая возвращает х1, если логическое условие А верно, и х2 в остальных случаях. Для этих целей так же удобно использовать функцию Хевисайда Ф(х), возвращающую 1, если х ≥ 0, и 0 – в остальных случаях.
Так, если до момента скачка tск дифференциальное уравнение, описывающее процесс, задавалось функцией f1(t, у(t)), а начиная с этого момента – функцией f2(t, у(t)), то модель процесса имеет вид dу(t)/dt = if [t < tск, f1(t, у(t)), f2(t, у(t)].
Решение разностного уравнения, соответствующего данному дифференциальному
уравнению, будет выполняться пошагово с переходом с функции f1 на f2 начиная с j-го шага, на котором значение tj достигнет tск. Аналогично, с помощью функции Хевисайда, можно записать:
dу(t)/dt = f1(t, у(t))[1-Ф(t- tск)] + f2(t, у(t)[Ф(t- tск)].
За счет вложений в логическую функцию количество скачков может быть любым. На
рисунке 1 показан результат моделирования процесса, развитие которого до момента t1 описывается функцией f1(t, у(t)) = tу1/3,
в интервале [t1, t2] – функцией f2(t, у(t)) = tу-1/2,
а с момента t2 – функцией f3(t, у(t)) = tу1/6.
Модель процесса имеет вид:
dу(t)/dt = if [t ≤ t2, if (t < t1, tу1/3, tу-1/2), tу1/6], у(0) = 1.
Рисунок 1 – Моделирование процесса с временными скачками.
Если описанное выше дифференцирование со скачкообразным изменением процесса
в заданные моменты времени, возможно, вообще говоря, как численными, так и ана
литическими методами, то учет скачкообразного изменения процесса в зависимости от достигнутого значения некоторого его параметра или от достигнутого значения самой дифференцируемой функции, возможен только при численном ифференцировании.
В основе учета не временного скачкообразного изменения исследуемого процесса
лежит оценка на каждом шаге вычислений достигнутого значения переключающего параметра (или самой функции) и принятие решения на "переключение" процесса по схеме, описанной выше. Для удобства моделирования целесообразно изменение переключающего параметра описывать дифференциальным уравнением в одной системе уравнений с дифференцируемой функцией.
Пример 1. Исследуемый процесс у(t) описывается дифференциальным уравнением
dу(t)/dt = tу1/3, у(0) = 1.
При достижении уровня у(t) = 100 процесс скачкообразно изменяется и развивается в соответствии с уравнением dу(t)/dt=tу1/6.
Модель процесса задается системой уравнений:
dу(t)/dt = tу1/3
Ф(100 - у2(t)] + tу1/6
Ф(100 - у2(t)], dу2(t)/dt = tу21/3.
Результат численного решения системы
приведен на рисунке 2 - а).
Пример 2. Исследуемый процесс у(t) описывается дифференциальным уравнением
dу(t)/dt = tу1/3, у(0) = 1.
Изменение переключающего параметра описывается уравнением s(t) = t3/3. При достижении параметром значения s(t) = 50 процесс останавливается и возобновляется в виде dу(t)/dt = tу1/6, у(0) = 1 при s(t)≥100.
Модель процесса задается системой уравнений:
dу(t)/dt = tу1/3
Ф(50 - s(t)] + tу1/6
Ф(100 - s(t)], ds(t)/dt = t2.
Результат численного решения системы приведен на рисунке 2 - б).
Рисунок 2 – Моделирование скачка процесса: а) в зависимости от достигнутого им значения; б) в зависимости от достигнутого значения переключающего параметра.
С использованием описанных приемов возможно моделирование различных ком-
бинаций скачков временного и не временного характера.
1.3. Приемы моделирования стохастических изменений процесса.
Пошаговый характер численного дифференцирования дает возможность использования на отдельных шагах вычисления разностного уравнения (1) датчиков случайных чисел. Эта возможность позволяет имитировать различные стохастические особенности исследуемого процесса – от его регулярного стохастического изменения (имитация различных шумов), до моделирования случайных событий в отдельных точках процесса (типа отказа устройства, поражения цели и т.п.).
На основе датчика равномерно распределенной случайной величины стандартными
приемами метода Монте-Карло могут моделироваться различные непрерывные случай-
ные величины с известными законами распределения (методом обратных функций) и дискретные случайные события с известными вероятностями (методом розыгрыша попадания равномерно распределенной величины на интервал, пропорциональный вероятности события).
1.3.1. Моделирование регулярных стохастических изменений.
Порядок моделирования рассмотрим на основе отдельных примеров.
Пример 3. Результаты моделирования случайных непрерывных процессов различ-
ной природы с шагом ∆t = 1 представлены на рисунке 3 (верхняя часть рисунка). Модели описываются следующими уравнениями:
а) случайный процесс, возрастающий на каждом шаге на величину rnd(1), где rnd(х) –
значение равномерно распределенной на интервале [0, х] случайной величины:
dу1(t)/dt = rnd(1), у1(0) = rnd(1);
б) гауссов случайный процесс, со среднеквадратичным отклонением на каждом шаге
σу = 10 и средним mу = 0. Значения случайной величины вычисляются по методу об-
ратных функций с использованием встроенной функции gnorm(α, mу, σу), возвращаю-
щей α% - ную квантиль нормального распределения с параметрами mу, σу:
dу2(t)/dt = gnorm(rnd(1), mу, σу), у2(0) = gnorm(rnd(1), mу, σу);
в) случайный процесс, возрастающий на каждом шаге на случайную величину,
имеющую показательный закон распределения с параметром λ = 0.1. Значения случай-
ной величины вычисляются по методу обратных функций:
dу3(t)/dt = - ln[1 - (rnd(1))] /λ, у3(0) = - ln[1 - (rnd(1))] /λ.
Рисунок 3 – Реализации, полученные при моделировании случайных непрерывных процессов различной природы.
Приведенные выше модели случайных процессов объединяет порядок пошагового
вычисления значений процесса. Здесь случайные величины представляют собой константы дифференциальных уравнений. Вследствие этого, не зависимо от метода численного решения, его сущность в соответствии с разностным уравнением (1) заключается в пошаговом весовом суммировании вычисляемых на каждом шаге значений случайной величины. Весовой коэффициент каждого слагаемого равен ∆t.
В нижней части рисунка 3 по результатам расчетов показаны значения случайных ве
личин уj,1, уj,2, уj,3, (j = 0, 1, …), вычисленных на каждом шаге моделирования соответствующего процесса и суммируемых при построении процесса с весовым коэффициентом ∆t.
Отмеченное взвешивание случайных величин определяет зависимость степени
"случайности" процесса от шага моделирования. Использование чрезмерно малого шага в заданном временном диапазоне моделирования приводит к "сглаживанию" процесса, снижению величины случайного отклонения на каждом шаге моделирования. В связи с этим при моделировании конкретных случайных процессов необходимо согласовывать выбор параметров процесса и шага моделирования. Аппарат численного решения дифференциальных уравнений в совокупности с встроенными функциями среды моделирования позволяет сформировать в решаемой системе уравнений самостоятельный датчик случайных чисел с определенными характеристиками, не имеющий описанного выше кумулятивного эффекта. Например, случайный процесс, близкий к гауссовому, синтезирован в виде результатов численного решения следующего дифференциального уравнения:
dу(t)/dt = [runif(1, m1, m2) – ceil(y)] /∆t, у(0) = runif(1, m1, m2),
где runif(n, m1, m2) – вектор n независимых случайных чисел, каждое из которых имеет равномерное распределение на интервале [m1, m2]; ceil(х) – наименьшее целое, не меньшее х. Одна из реализаций процесса, полученных при m1 =-1, m2 = 2, приведена на рисунке 4.
Рисунок 4 – Реализация синтезированного случайного процесса.
Формирование случайного процесса может осуществляться за счет аддитивного или
мультипликативного наложения на детерминированную функцию значений случай-
ной величины, рассчитанных в отдельных точках. Одним из вариантов формирования случайного процесса может быть также использование детерминированной функции в качестве параметра закона распределения случайной величины, определяемой пошагово методами Монте-Карло (розыгрышем событий или методом обратных функций).
Соответствующие процессы приведены в следующем примере.
Пример 4. В качестве детерминированной составляющей используется функция, яв-
ляющаяся решением рассмотренного ранее уравнения dу(t)/dt = tу1/3
при у(0) = 100. Шаг моделирования ∆t = 1.
а) аддитивная смесь с равномерно распределенной в интервале [0, 10] случайной величиной:
dу1(t)/dt = tу11/3 + rnd(10);
б) мультипликативная смесь с равномерно распределенной в интервале [0, 2] слу
чайной величиной:
dу2(t)/dt = tу21/3 •rnd(2);
в) гауссов процесс с динамически изменяющимся средним mу(t) = у(t) и постоян-
ным среднеквадратичным отклонением σу =50:dу3(t)/dt = gnorm(rnd(1), tу31/3, 50).
Результаты моделирования приведены на рисунке 5.
Рисунок 5 – Реализации, полученные при моделировании случайных процессов с детерминированной составляющей.
1.3.2. Моделирование отдельных случайных событий и нерегулярных стохастических изменений процесса.
Моделирование может осуществляться по схеме, описанной в пункте 1.2, в зависимости от момента времени, или от значения переключающего параметра (или самой функции). В отличие от описанного в пункте 1.2 детерминированного случая, решение об изменении процесса должно приниматься на основе розыгрыша случайного события с известной вероятностью. Событие считается свершившимся, если значение равномерно распределенной величины попадает на интервал, пропорциональный вероятности p рассматриваемого события. Математически это условие записывается логической функцией if[p < rnd(1), х1, х2].
Пример 5. Рассматриваются два гауссовых процесса с известными параметрами,
описывающих траектории движения цели и наводимой на нее ракеты. В момент встречи ракеты с целью происходит поражение цели с вероятностью p. Необходимо оценить распределение по времени случайных моментов поражения цели.
Модель процесса обстрела цели опишем следующей системой дифференциальных
уравнений :
dу1(t)/dt = gnorm(rnd(1), tу11/3, 50), у1(0) = 100,
dу2(t)/dt = gnorm(rnd(1), tу21/5, 50), у2(0) =1000,
ds(t)/dt = if [у2 ≤ у1, if (p < rnd(1), 0, 2000), 0], у(0) = 10.
Первое и второе уравнения описывают соответственно траектории ракеты и цели.
Третье уравнение представляет собой индикаторную функцию, предназначенную для отображения результатов розыгрыша случайного события поражения цели с вероятностью p в момент пересечения траекторий. Две реализации моделирования процесса обстрела цели для р = 0.8 приведены на рисунке 6. На первой реализации разыгрываемое событие поражения цели произошло после пересечения траекторий – цель не поражена. На второй реализации событие свершилось в точке пересечения траекторий – цель поражена.
Рисунок 6 – Реализации моделирования процесса обстрела цели: а) цель не поражена; б) цель поражена в точке встречи в момент t = 23.
Приступим к практике…..
Т.к. полная математическая модель ВБ достаточна обширна, рассмотрим частный случай - ситуацию, столь любимую
Breeze-ом, а именно – F-22 (слева) при сближении с самолетом противника F-15 (справа) обнаруживает его на дальности, когда противник его еще не видит, и производит пуск УРВВ.
2. Пример моделирования процесса наведения ракеты УРВВ на
истребитель противника.
2.1. Описание моделируемого процесса.
Моделирование осуществляется в координатах R-H "дальность-высота" (рисунок 7).
Положение цели на момент начала моделирования считается заданным. Траектория цели определяется в общем случае динамически изменяющимися характеристиками ее вектора скорости – абсолютной величиной Vц(t) и углом Θц(t), составляемым с направлением горизонтальной оси координат. Траектория ракеты определяется методом наведения ее на цель, который позволяет рассчитать в каждый момент времени после пуска ракеты потребное расстояние между ракетой и целью D(t) и угол визирования цели из местоположения ракеты ϕ(t).
Рисунок 7 – Схема расположения истребителей.
В начале системы координат размещены истребитель F-22. Его РЛС осуществляет циклическое, с периодом ∆t, измерение параметров траекторий цели. Наведение ракеты осуществляется в соответствии с методом пропорциональной навигации. Пуск ракеты по цели производится в момент t = 0.
Уравнения метода наведения:
dD(t)/dt = Vц(t) cos(Θц(t) - ϕ(t)) - Vр(t) cos(Θр(t) - ϕ(t)),
dϕ(t)/dt = [Vц(t) sin(Θц(t) - ϕ(t)) - Vр(t) sin(Θр(t) - ϕ(t))] /D(t),
D(0) = Rц0, ϕ(0) = arctg(Hц0/Rц0), Θр(t) = kϕ(t) + c,
где k, c – константы метода, реализованные в конкретной системе вооружения.
Скорость цели Vц(t) оценивается РЛС F-22 с погрешностью, распределенной по
нормальному закону с известным значением дисперсии σVц2. Скорость ракеты Vр(t) имеет переменный характер – возрастает по закону Vр.н(t) с момента старта на разгонном участке полета и снижается по закону Vр.п(t) на пассивном участке, после окончания работы двигателя. Момент окончания работы двигателя ракеты является равномерно распределенной в интервале [tдв.1, tдв.2] случайной величиной.
Требуется оценить ошибки наведения ракеты на цель.
2.2 Построение модели и результаты моделирования.
Дополним исходную систему уравнений зависимостями, описывающими изменение координат цели Rц(t), Нц(t) и ракеты Rр(t), Нр(t):
dRц(t)/dt = Vц(t) cosΘц(t), dНц(t)/dt = Vц(t) sinΘц(t),
dRр(t)/dt = Vр(t) cosΘр(t), dVр(t)/dt = Vр(t) sinΘр(t),
Rц(0) = Rц0, Hц(0) = Hц0, Rр(0) = 0, Vр(0) = 0.
Для учета в модели случайной погрешности измерения скорости цели и случайного
характера окончания работы двигателя ракеты используем приемы, описанные в при
мерах 4,5. В соответствии с этим, в исходной системе уравнений метода, а также в уравнениях для расчета значений Rр(t), Vр(t) произведем замену:
Vц(t) → gnorm(rnd(1), Vц(t), σVц),
Vр(t) → if [t ≤ runif(1, tдв.1, t дв.2), Vр.н(t), Vр.п(t)].
Реализации моделирования процесса наведения ракеты на цель, при движении ее по разным траекториям, приведены на рисунке 8. При моделировании использовались средние значения параметров РЛС и УРВВ .
Рисунок 8.
1) Реализации моделирования процесса наведения ракеты: а) – на прямолинейно движущуюся и б) – на маневрирующую по закону синуса цель;
2) Гистограмма промахов ракеты по маневрирующей цели.
Моделирование показывает, что в точке встречи от реализации к реализации имеет место разброс координат траектории ракеты, обусловленный действием случайных факторов.
Гистограмма разности координат ∆R= Rц - Rр в точке встречи ракеты с маневрирующей по закону синуса целью получена
по сорока реализациям моделирования процесса.
Применимость построенной модели не ограничивается приведенными результата-
ми.
Модель позволяет исследовать разнообразные условия сближения истребителей и наведения ракет на цель, в том числе случайные факторы, связанные с ошибками измерения различных координат цели и ракеты. Незначительные доработки модели по схеме, приведенной в пункте 1.3.2. позволяют оценивать не только траекторные характеристики , но и возможности поражения цели.
Пример моделирования процесса изменения численности ЛА воюющей группировки
приводить не буду, т.к. не актуально, да и надоело уже…..
П.С. - Выражаю огромную благодарность д.т.н. Горевичу, за предоставленный материал.
