Моделирование систем дифференциальных уравнений с динамическими инвариантами Текст научной статьи по специальности «Математика»

УДК 517.925:519.2:519.6

DOI: 10.18698/2309-3684-2019-1-98117

Моделирование систем дифференциальных уравнений с динамическими инвариантами

© Е.В. Карачанская1'2

Дальневосточный государственный университет путей сообщения, Хабаровск, 680021, Россия

2Тихоокеанский государственный университет, Хабаровск,680035, Россия

В статье рассматривается метод построения систем дифференциальных уравнений, имеющих заданный набор гладких функций в качестве первых интегралов - глобальных инвариантов. Данный алгоритм позволяет строить как системы детерминистических дифференциальных уравнений, так и стохастических дифференциальных уравнений Ито. Стохастические уравнения могут быть диффузионными (с винеровскими возмущениями) или диффузионными со скачками, вызванными скачками пуассоновского процесса. Предложенный алгоритм опирается на предыдущие работы автора и не имеет аналогов. Приведено детальное описание алгоритма в фазовых пространствах размерности 2 и 3. Представляется MathCad-программа, автоматизирующая процесс построения систем детерминистических дифференциальных уравнений и систем стохастических дифференциальных уравнений Ито с винеровским процессом. Приведены примеры автоматизированного построения систем дифференциальных уравнений разных типов. Правильность построения систем уравнений проверяется с помощью численного решения полученного уравнения и определения значения функции, объявленной первым интегралом, в зависимости от решения этой системы (для детерминистических- уравнений) или реализации решения (для стохастических уравнений). Решение систем уравнений производится классическими методами - Эйлера, Эйлера-Маруямы и с использованием метода статистического моделирования Монте-Карло. Решения уравнений и значения функции - инварианта представлены графически. Применение представленной теории и разработанной программы MathCad продемонстрировано на примере построения программных управлений с вероятностью 1 (РСР1) для стохастической S1R-модели. Результаты работы могут быть использованы для построения модели динамической системы с инвариантами и дальнейшего исследования таких систем.

Ключевые слова: первый интеграл, система дифференциальных уравнений, стохастические уравнения, алгоритмизация процесса, MathCad

Введение. Модели реальных явлений всегда основываются на определённых ограничениях. Для моделей детерминированных процессов такие ограничения присутствуют, например, в форме закона сохранения. В условиях, когда уравнения, описывающие физическую систему, ещё не известны, законы сохранения являются эффективным способом её исследования, позволяющим иногда установить и сами уравнения.

Включение в рассмотрение неопределённости существенно изменяет подходы к моделированию. Однако и в этом случае динамическая система может иметь инварианты [1]. Задача

построения уравнений движения по заданным свойствам движения и исследование устойчивости этого движения является одной из обратных задач динамики. Решение этой задачи состоит в построении системы дифференциальных уравнений (ДУ), имеющих данную интегральную кривую, и установлении условий устойчивости (в том или ином смысле). Традиционно для первой части используется метод, предложенный Н. П. Еругиным [2], далее развитый Р. Г. Мухарлямовым, И. А. Мухаметзяновым и их учениками. Предложенные методы основаны во многом на линеаризации функций правой части, что не всегда даёт адекватный результат. В 2001 г. М. И. Тлеубергенов при рассмотрении обратных задач подошёл к построению систем СДУ Ито, содержащих случайную составляющую только в виде винеровского процесса [3]. В 1989 г. В. А. Дубко [4] предложил метод построения системы дифференциальных уравнений (СДУ) Ито на основе инварианта - функции, которая должна с вероятностью 1 сохранять постоянное значение на любом решении строящейся системы СДУ - первого интеграла этой системы уравнений.

Особенность метода построения систем ДУ, предложенного Н. П. Еругиным и развитого Р. Г. Мухарлямовым, И. А. Мухаметзяновым и их учениками состоит в том, что авторы восстанавливали систему на основе локальных первых интегралов (интегральных кривых), связанных с конкретными начальными условиями.

Предлагаемый в работе метод отличен от метода Н. П. Еругина и методов, на нём построенных, и не является ни обобщением, ни частным случаем описанных выше. В данном методе для построения системы ДУ используются глобальные первые интегралы при любых начальных условиях и вопрос об устойчивости инварианта снят - его значение сохраняется с вероятностью 1, подобно [4].

В работах [5-7] был предложен алгоритм построения систем СДУ Ито с пуассоновскими скачками (ОСДУ Ито) со свойством сохранения с вероятностью 1 инварианта, однако его применение достаточно сложно технически в связи с вычислением множества определителей векторно-функциональных матриц, связанных с матрицами Якоби. Однако данный алгоритм может служить основой для программного продукта, позволяющего упростить этот процесс получения системы дифференциальных уравнений, подбора функций, обеспечивающих удобство моделирования и визуализации полученных результатов.

Целью статьи является представление алгоритма и программной реализации аналитического построения моделей динамических систем, как детерминистических, так и стохастических, для которых семейство гладких функций определяет множество глобальных

первых интегралов строящейся системы дифференциальных уравнений.

Исходные обозначения и предварительные положения.

Обозначим: {ё}.}п=0={ё0,ёх,...,еп} — ортонормированный базис в

™ -гт-р, п \—7 / ч ды (t, x) ^ -А ды и, x) ^ _ _

х Ж , Vы(t, x) = —1—- е0 + > —1—-— обобщенный градиент дt 1=1 дхг

функции ы(г, х), х е Еп с — случайное событие; w(t) — винеровский вещественно- значный процесс (ВП), у(Аг, Ау) — пуассоновская мера (ПМ), т.е. случайная величина с распределение Пуассона, имеющая свойства меры; Я(у) — пространство значений вектора у е К4. Рассматриваемые ниже w(t) и у(Аг, Ау) являются взаимно-независимыми.

Рассмотрим задачу Коши для ОСДУ Ито

¿х1 (г) = а (г, x(t+ Ьа (г, x(t))dw, (г) + (г, x(t), у)у(Л, dу), (1) x(0) = x0, I = 1,..., п,

с условиями: а (г, x) е С/Х, Ьц (г, x) е С^Х, gi (г, x, у) е С1^, и функция

x(t, у) е Еп — ее решение, зависящее от начального условия. Указанные условия являются дополнительными к стандартным условиям существования и единственности решения системы (1) [8].

Определение [9]. Неслучайную функцию ы (г, x) будем называть первым (детерминированным) интегралом (ПИ) системы ОСДУ (1), если она с вероятностью 1 сохраняет постоянной значение для любой реализации случайного процессаx(t,x0,с) являющегося решением этой системы: ы(г, x(t, x0, с)) = ы (0, x0).

Поскольку первый интеграл динамической системы формально зависит от времени, его можно будем считать динамическим инвариантом.

Следует отметить, что ПИ может быть определен только в пространстве размерности п > 2, поскольку наибольшее число независимых ПИ на единицу меньше размерности пространства [1,4] (также этот вывод можно получить, опираясь на второе утверждение Теоремы 2, сформулированной ниже). Понятие ПИ системы СДУ имеет важное значение, например, для построения модели динамической системы, обладающей инвариантами [7, 10] и построения программных управлений с вероятностью 1 (РСР1) [5, 6, 11, 12] и др.

Применение СДУ к моделированию реальных систем открывает перспективу для использования первого интеграла таких уравнений

для решения подобных задач, связанных со стохастическими системами.

Для случая только винеровских возмущений (классическое уравнение Ито) в [1] исследование свойств ПИ как детерминированной функции было связано с установлением независимости такой функции от реализаций ^(1). При добавлении пуассоновских возмущений (обобщённое уравнение Ито) такое требование приводит к следующим условиям.

Теорема 1. Пусть случайный процесс x(í) является решением

системы ОСДУ (1). Функция и , x)е С/х2 есть первый интеграл системы (1) тогда и только тогда, когда она удовлетворяет условиям:

• ^^ + (,, ^ - 2 Ь„ (^ ^ ] = 0;

д1 дх, 2 дх,

К Ц, x)

оди(1, x, с) _

дх,

= 0 для всех к = 1, т;

• для любых у е М4 во всей области определения процесса: и(?, x, со) - и(?, x + g(?, x, у), со) = 0.

Получение данного результата стало возможным после пприменения обобщенной формулы Ито-Вентцеля - специального «правила дифферен цирования» для случайной сложной функции [9,

13].

Построение системы стохастических дифференциальных уравнений Ито по заданному набору функций. Построение системы п дифференциальных уравнений по известному множеству его первых интегралов опирается на векторное произведение в пространстве Мп. В п -мерном пространстве с базисом{е, }П=1 любой

вектор а, ортогональный множеству векторов т определяться следующим образом:

-(1)

}п

-]>]=1 *(п-1)

будет

а = Л -[т

(1)

-(п-1)

Т ] е

(

Л- det

(1)

ЧТ1

( п-1)

2

(1)

(п-1)

п

(1)

(п-1)

Л е М

Это произведение отлично от нуля в случае, если векторы линейно

ч,)

независимы. Введем вместо векторов т" обобщённые градиенты V/, (/, x), тогда линейная независимость указанных векторов

пространства [0, Т ] х Мп соответствует независимости функций /} (г, х), у = 1п.

Построение системы ОСДУ Ито по заданному набору функций. Определим вид системы ОСДУ Ито с начальными данными, имеющей известный СПИ.

Теорема 2 [3]. Пусть х : М+ ^ Мп, п > 2, и(г, х) е С^2, х(г, х0) — решение ОСДУ Ито:

ёх(г) = А(г, х(г))ёг+В (г, х(г))<^(г) + Г g(г, х(г), у)у(ёг, ёу)

^ (у)

х(0) = х0, г > 0.

(2)

Пусть ^ (г, х), у = 1, п -1, и рр (г, х(-,у)), р = 1, п - 2, — произвольные функции из С/Х и С^Х'У, независимые с функцией и(г, х), и пусть Н (г, х) — матрица вида

Н (г, х) =

ди(г, х) ди(г, х)

ди (г, х)

дг дх1 " дхп

т х) д/г(г, х) д/(г, х)

дг дх1 дхп

д/п-1 (г, х) д/п-1 (г, х) д/п-1 (г, х)

дг

дх,

дхп

Если и (г, х) является СПИ системы (2), тогда:

1. матрица В (г, х) состоит из столбцов Вк (г, х) = Ьк (г, х)ег. к = 1, т, вида Вк (г, х) (г, х) • det М (г, х)|, где М (г, х) — минор, соответствующий элементу/п0 матрицыН(г,х);

2. А(г, х) имеет вид: А(г, х) е| R(г, х) +1 [дВк(г,х)]• Вк(г,х) 1, где

[ 2 дх ]

Я (г, х) — матрица-столбец, компоненты г. (г, х), . = 1, п, которой определяются следующим образом: С-1 (г, х) • ёй Н(г, х) = е0 + г (г, х)ег.; С (г, х)— алгебраическое дополнение элемента е0 матрицы Н (г, х),

ёе! С (г, х) ф 0 и [ Вк (г, х);

дВк (г, х)

дх

] — матрица Якоби для векторной функции

0

п

3. коэффициент g(г, x,у) = gi (г, x, у)е. при пуассоновской мере определяется представлением, g (г, x, у) = у (г, x, у) - x где у (г, x, у) — решение системы дифференциальных уравнений, построенных на основе векторного произведения

ду(-,у)

ду

ди(г, у(-,у)) , у(-,у))

д^(г, у (-,у)) °ух

^2

д^(г, у(- у))

^2

дФп-2(г, у(-,у)) ^п^ у(-,у))

°У2

^ y(•, у))

^п

дPl(г, у(-,у))

^п

^-2^ У (-,у))

дУп

удовлетворяющее начальному условию у (г, x, 0) = x.

Доказательство теоремы приведено в [7].

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

Описание алгоритма в М2. Набор исходных данных: функция и (г, x) — первый интеграл; дополнительная функция ^г, x) линейно

независимая с и (г, x) и произвольная функция ^0(г, x) . В М2 положим x = ( х У)*.

Построение системы обыкновенных ДУ (без случайных возмущений). Вычисление векторного произведения на основе частных производных функций и (г, x) ^г, x) в соответствии с теоремой 2

det

ди (г, х, у) ди (г, х, у) ди(г, х, у)

дг

дИ(г, х, у) дг

х

дИ(г, х, у ) х

у

дИ(г, х, у)

ду

(3)

= а0 (г, х, у)в0 + ах (г, х, у)ех + а2 (г, х, у)ё2.

Определение алгебраических дополнений элементов первой строки определителя (3) и компонент вектора - коэффициента при (1г:

2

п

0

2

4(г, х, у) =

_ ах(г,х,у)

С1(г,х, у)

, х, у) =

_ а20,х, у)

С1(г,Х, у)

Векторная форма коэффициента при детерминистического уравнения:

Л(г, х(г), у (г )) =

"4(г, х(г), у (г)) Л2(г, х(г), у (г))_

Вывод результата ОДУ: Система обыкновенных дифференциальных уравнений с первым интегралом и (г, х, у), дополнительной функцией ^г, х, у) и произвольной функцией д0(г, х, у) (без случайных воздействий:

йх(г)

йг йу(г) йг

= Л(г, х(г), _у(г)).

Построение системы СДУ Ито с винеровскими возмущениями. Система строится с использованием предыдущего шага.

Построение векторного произведения на основе частных производных функции и (г, х):

Б(г, х, =

с2

ди(г, х, у ) ди (г, х, у)

дх

ду

4(г, х, у) й2(г, х, у)

Определяем компоненты вектора В1(г, х, у) — столбца матрицы при винеровской составляющей:

_ ^ x, у)

В1 (г, х, у) = д0 (г, х, у) П(г, х, у) =

Ъ2(г, x, у)

Определяем стохастическую добавку в коэффициент сноса: дЪ1 (г, х, у) дЪ1 (г, х, у)

НР(г, х, у) =

дх ду

дЪ2 (г, х, у) дЪ2 (г, х, у)

Ъl(г, х у) Р1(г, x, у)

Ъ2(г, х у)_ _Р2(г, x, у)_

дх ду

Определение компонент коэффициента сноса - коэффициента при

йг:

4и(г, х, у) = 4(г, х, у) +

Ап(г, х, у ) = А2(г, х, у) + Коэффициент сноса: А0 (г, х(г), у(г)) =

Pl(г, x,у) 2

P2(г, X, у) 2 .

А01 (г, х(г), у (г)) _ А02(г, х(г), у (г))_

Коэффициент при винеровской составляющей для скалярного ВП: В(г, х(г), у (г )) = В^г, х(г), у(г))

Вывод результата СДУ1: Вывод результата СДУ1: система стохастических уравнений Ито с одномерным винеровским процессом w(г), имеющая первый интеграл и (г, х(г), у (г)), с дополнительной функцией ^г, х, у) и произвольной функцией д0(г, х, у):

дх(г) йу(г)

= А0 (г, х(г), _у(г ))ёг+В(г, х(г), _у(г)) • w(г).

Коэффициент при винеровской составляющей для двумерного ВП w(г ) = К(г )^(г ))*:

ВВ(г, х(г), у (г)) =

"ь,(г, х(г), у (г)) Ь1(г, х(г), у(г)) Ь(г, х(г), у (г)) Ь2(г, х(г), у(г))

Вывод результата СДУ2: система стохастических уравнений Ито с двумерным винеровским процессом, имеющая первый интеграл и (г, х, у), с дополнительной функцией ^г, х, у) и произвольной функцией д0 (г, х, у) :

ёх(г) Ф (г)

= А (г, х(г), у(г))ёг + ВВ(г, х(г), у (г)) •

^(г)" w2(г).

Построение системы ОСДУ Ито. К сожалению, полностью автоматизировать процесс построения системы СДУ Ито с винеровскими возмущениями и пуассоновскими скачками невозможно, поскольку требуется решить систему уравнений с частными производными вида, неизвестного заранее. Кроме того, не всегда построенная на основании теоремы 2 система будет иметь решение в явном виде, а именно явный вид решения указанной выше системы позволяет определить коэффициент для слагаемого при ПМ. В некоторых случаях описанную проблему можно решить за счет подбора дополнительных и произвольной функций.

В случае, когда матрица g(г, х, у, у) =

определена,

gl(г, хy, у)

_ g2(г, x, У, у) _

то она вводится дополнительно. Тогда в алгоритм добавляется следующий шаг.

Вывод результата СДУ1Р: система стохастических уравнений Ито со скалярным винеровским процессом w(г) и пуассоновской составляющей, имеющая первый интеграл и(г, х, у) , с дополнительной функцией Н(г, х, у) и произвольными функциями д0(г, х, ^)и ^(г, х, у) (последняя функция вводится из соображений удобства моделирования):

¿х(г)

йу(г)

= А (г, х(г), _у(г ))ёг + В(г, х(г), _у(г)) • w(г) +

+1В( )g(г, х, у, у) , ¿у).

Для двумерного ВП результат аналогичен предыдущему.

Описание алгоритма в М3. Отличием от алгоритма в случае М2 является набор исходных функций и возможность учета трехмерного винеровского процесса. Введение исходных данных: функции и. (г, x) — первого интеграла, дополнительной функции Н} (г, x) ,

линейно независимых с и. (г, x) и произвольной функции ^0(г, x), . + у = 3. Процесс построения систе-мы ДУ аналогичен ранее изложенному.

Программная реализация построения системы ДУ в М2 и М3.

Программа построения системы дифференциальных уравнений на основе инвариантов в М2 и М3 реализована в MathCad. Для построения системы дифференциальных уравнений входными данными являются: функции и. (г, x) — инвариант (ПИ), дополнительные функции Н^ (г, x)

и произвольная функция (г, x). Общее число инвариантов и дополнительных функций равно размерности фазового пространства (п = 2,3). В результате получаем систему обыкновенных ДУ (без стохастической части), систему СДУ Ито с винеровским процессом -скалярным и двумерным. Для решения системы детерминистических и стохастических ДУ используются известные метод Эйлера и стохастический метод Эйлера соответственно [14]. В силу того, что вспомогательная система ДУ, построенная согласно п.3 Теоремы 2, может быть любого вида, нет возможности ее численного решения «универсальным» методом, поэтому построение системы СДУ с пуас-соновскими скачками программно не реализуется, т.е. может решаться

только непосредственно. Кроме того, система не всегда имеет решение в явном виде. Если вспомогательная система ДУ имеет решение в явном виде, то возможно построение системы СДУ со скачками и далее находится ее численное решение с помощью метода Монте-Карло. Поскольку целью программной реализации является аналитическое построение системы ДУ, а численное решение построенных систем уравнений используется только для демонстрации правильности метода построения, то вопросы сходимости и т. п. не рассматриваются.

Классический метод Эйлера. Пусть дана задача Коши для системы уравнений первого порядка:

^ = / (г, х), х(г )|г=0= хо,

йг

где функция / определена на некоторой области V е [0, Т] х Мп. Решение ищется на интервале изменения параметра г: г е [0, Т]. На этом интервале введем узлы: 0< г1 < ... < гы = Т, N — число точек разбиения, зависит от шага разбиения. Численное решение в узлах гг, которое обозначим через хг, определяется по формуле:

х = х,._1 + И/ (г _1, х,._1), г = 1,2,3,., N,

где И = г{ - гг-1 — шаг разбиения.

Метод Эйлера - Маруямы. Рассмотрим задачу Коши для системы стохастических дифференциальных уравнений Ито

йх(г) = Л(г, х(г ))йг+в (г, х(г ))<^(г), х(г) |г=0= х0.

Численное решение определяется по формуле:

х = х..-1 + ИА(г-1,х-1)+4нв(г .-1,х-Х,, г = 1,2,3,.,N,

где И = гг - гг-1 — шаг разбиения, £ — значение случайной величины, имеющей нормальное распределение £ ~ N(0,1).

Стохастический метод Эйлера для уравнения со скачками.

Численное решение задачи Коши для системы СДУ Ито со скачками

йх(г) = Л(г, х(г))йг+в(г, х(г))dw(г) + Г g(г, х(г), у)у(йг, йу),

^ (у)

х(г) |г=0 = х0.

определяется по формуле:

х = х-1 + ИА(гг-1, хг-1 ) + ^В(гг-1, х.-1 К + ^-1, г =1,2,3, . , N,

где И = ti - ti_1 — шаг разбиения, £ — значение случайной величины, имеющей нормальное распределение £ ~ N(0,1), — слагаемое, определяющее величину случайного скачка в момент времени ti, разыгрываемое, например, по методу Монте-Карло:

Г8^,■-1,^-1), г-1 < Л_1;

п_1 = 10

[0, иначе,

Г_1 — значение случайной величины, равномерно распределённой на промежутке (0,1). Число Л}-1 также может быть случайным, любой природы, в том числе равномерно распределённым.

Примеры построения моделей динамических систем с инвариантами. Построение систем ДУ производится MathCad-программой. Будем фиксировать только результаты.

Пусть функция и(г, х, у) = ^ + е^ + х - у интерпретируется как первый интеграл (инвариант) при начальном условии х(0) = 0, у(0) = 1, дополнительная функция И^, х, у) = у - х2 +1 и д^, х, у) = х — произвольная.

Система без возмущений. Строим систему детерминистических уравнений, для которой функция и^, х, у) = е1 + е- + х - у — первый интеграл. Система уравнений имеет вид:

dx(t) = et - е^ +1 dy(t) = 2х^)(е^ - et) -1 (4)

dt 2x(t) -1 dt 2х(г) -1

Результаты численного моделирования решения системы уравнений (4) методом Эйлера приведены на рис. 1. Указанная функция является инвариантом для данной системы уравнений и сохраняет значение 1.

Возмущения в виде одномерного винеровского процесса. Построенная на основе функций и ^, х, у) и И^, х, у) система СДУ имеет вид:

, , ч Г е^^ - е^ +1 х(0 ^ , , , Л dx(t ) = -+ —— dt - х^ )dW(t),

. 2х^)-1 2 ,

(5)

,/ч Г 2 x(t )(e-t - et) -1 /ч

(]у^)=--^-— + dt - x(t).

. 2x(t) -1 2 )

Результаты численного моделирования решения системы уравнений (5) методом Эйлера приведены на рис. 2. Использованный метод применяется к линейным уравнениям, а стохастическое

уравнение (5) нелинейное, то возникающие погрешности отразились на значениях одной из данных функций - дополнительной функции. При любых реализациях решения СДУ функция и (г, х(г), _у(г)) = ег + е- + х(г) - у (г) сохраняет постоянное значение, равное 1.

Рис. 1. Численное моделирование решения системы уравнений (4) и значений функции и (г, х(г), у (г)) в MathCad

:= У. + Ы ■ + 7Й1 ■ I ' »-2. «2 шопп(51ер51,0,1)

(Ч

(Ч_

иу, • • •

0.4 0.5 0.6 0.7 0.8 0.9

Начальное условие

'( хО

У0 :=

ЦУ. = 1

К)

Рис. 2. Реализация численного моделирования решения системы уравнений (5) и значений функции и (г, х(г), у(г)) в MathCad

Диффузионный процесс со скачками. Используется уже построенная система (5). Обозначим: х = (х,у)* и z = (г2)*. Коэффициент g (г, х, у) при ПМ определяется равенством g (г, г,у) = г (г, х,у) - х, где г(г, х,у) — решение системы ДУ из теоремы 2, удовлетворяющее начальному условию г (г, х, 0) = х. Для построения ОСДУ необходимо получить явное решение

соответствующей системы ОДУ в аналитическом виде. Приведем полностью процесс получения этого решения. В соответствии с утверждением теоремы 2 определим частные производные функции и (г, z):

Си (г,г) _1 си (г, г) _

— 1,---1.

Теперь строим систему ДУ:

dz(t, x, у)

ду

dz2

f dZj(t, x, у) ^ ду dZ2(t> Х,У) . дУ .

M

v-Ъ'

Её решение составляют функции г1 (г, х, у) — -у + С1 (х), (г, х, у) — -у + С2 (х). Начальное условие , х, у) |у—0 — х определяет функции С1 (х) — х и С2 (х) — у. Тогда г1 (г, х, у) — -у + х, г2(г, х, у) — —у + у. С учетом g (г, г,у) — г (г, х, у) - х, получаем:

g(г, х, у, у) — (—у, —у)* Соответственно, диффузионная система со скачками имеет вид:

dx(t) =

f e - e-t+1

x(t)

2 x(t) -1 2

Л

- +

dt - x(t)dw(t) - F yv(dt, dy),

JR (y)

dy(t) =

2x(t)(e-t - el) -1 x(t) 2

(6)

2 x(t) -1

dt - x(t)dw(t) -F yv(dt, dy).

JR(y)

Результаты численного моделирования решения уравнения (6) представлены на рис. 3. На любых траекториях решения уравнения (6) функция u(t, x(t), y (t)) = e* + e- + x(t) - y(t) является первым интегралом данной системы.

Применение: управляемая SIR-модель. Построение системы дифференциальных уравнений с заданным первым интегралом дает возможность определения управлений в динамической системе.

Определение. Программным управлением с вероятностью 1 (PCP1 - Programmed Control with Probability 1) будем называть такое управление в стохастической системе, которое с вероятностью, равной 1, обеспечивает сохранение постоянного значения для заданных функций, зависящих от положения системы в любой момент времени.

Рассмотрим применение PCP1 и изложенной ранее теории.

Рис. 3. Численное моделирование решения системы уравнений (6) и значений функции и(/, х(0, у (()) в MathCad

Распространение инфекционных заболеваний представляет собой сложное явление с множеством взаимодействующих факторов. Ключевая роль математической эпидемиологии заключается в создании моделей распространения патогенов. Эти модели служат в качестве математической основы для понимания сложной динамики распространения заболевания. Основа таких моделей - модель SIR (Susceptible - Infected - Removed), предложенная в 1927 г. У. Кермаком и А. МакКендриком. В дальнейшем на её основе было построено множество моделей для различных заболеваний, например,

в [16].

Введем обозначения: S (t ) — доля здоровых, но восприимчивых к инфекции людей; I (t) — доля инфицированный (больных), распространяющих инфекцию людей; R (t) — доля обладающих иммунитетом к болезни (изначально невосприимчивых,а также переболевших); X — интенсивность заражения, л — интенсивность выздоровления. Соответственно, классическая SIR-модель имеет вид:

^ = -XS (t )I (t )

ddtt) = XS (t )I (t ) -лI (t) (7)

dt

^=¡S (t ) I (t ), dt

Будем считать, что при этом общая численность людей - величина постоянная:

S (t ) +I (t ) + R (t ) = 1.

(8)

В классическую модель внесём сильные возмущения в виде непрерывного, вызванного одномерным винеровским процессом, и скачкообразного, сопоставленного с пуассоновским:

dS (t ) = -AS (t ) I (t )dt + b1 (t, S (t), I (t), R(t ))dw(t ), dl (t ) = (AS (t ) I (t ) - /I (t ))dt + b (t, S (t ), I (t ), R(t ))dw(t ) + +jR(^)g2 (t, S (t ), I (t ), R(t ), y)v(dt, dy),

dR(t ) = /S (t )I (t )dt + b (t, S (t), I (t ), R(t ))dw(t)

(9)

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

dS (t ) = (-AS (t ) I (t ) + s (t, S (t ), I (t), R(t )))dt +

+ (b (t, S (t ), I (t ), R(t )) + ^ (t, S (t ), I (t ), R(t )))dw(t ) dI (t ) = (AS (t ) I (t ) - /I (t ) + s2 (t, S (t ), I (t ), R(t )))dt + + (¿2 (t, S (t ), I (t ), R(t )) + й (t, S (t ), I (t ), R(t )))dw(t ) + , (io)

+ f )( g 2 (t, S (t ), I (t ), R(t ), y) + % (t, S (t ), I (t ), R(t )))v(dt, dr)

jR(f)

dR(t ) = (/S (t ) I (t ) + s3 (t, S (t ), I (t), R(t )))dt + + (¿3 (t, S (t ), I (t ), R(t )) + (t, S (t ), I (t ), R(t )))dw(t )

Таким образом, нужно определить PCP1 для системы (10), чтобы функция

и (t, S (t ), I (t ), R(t )) = S (t ) +1 (t ) + R(t ) -1 (11)

будет всегда принимать значение 0, т.е. была первым интегралом.

Положим: первый интеграл — функция и (t, x, y, z ) = x + y + z -1,

дополнительные функции — v(t, x, y, z) = 2e-t + x и h(t, x, y, z) = y, произвольная функция — q (t, x, y, z) = x . Соответственно, начальные условия: x(0) = 1, y(0) = 0, z(0) = 0 . Построенная система дифференциальных уравнений имеет вид:

v(dt, dy) (12)

dx(t ) "2e-i " " 0 " " 0 "

dy(t ) = 0 dt + x(t ) dw(t ) + f JR(y) x(t )y

dz(t ) -2e-t - x(t )_ - x(t )y_

Результат моделирования решения построенной системы без скачков - на рис. 4, при наличии скачков - на рис. 5.

Рис. 4. Численное моделирование решения системы уравнений (12) без скачков и значений функции и(/, х(/), у (/)) в MathCad

Рис. 5. Численное моделирование решения системы уравнений (12) со скачками и значений функции и(/, х(/), у (/)) в MathCad

Как видим, в системе (12), в отличие от (9), присутствует скачкообразная компонента. Поэтому необходимо ввести фиктивно скачкообразную компоненту с управлением , £^), I^), Я^)) в последнее уравнение системы (10). Теперь можно определить управления я ^, £ ^), I ^), Я({))),'\ = 1,2,3, из сопостравления коэффициентов систем (10) и (12):

s1(t, 5(г), I(г), Щ(г)) = 2е-г - ХБ (г)I(г); ^ (г, Б (г), I (г), щ(г)) = -Ъх (г, Б (г), I (г), щ(г)); s2 (г, Б (г), I (г), щ(г)) = ¿I (г) - ХБ (г) I (г); й(г, Б (г), I (г), щ(г)) = х(г) - Ъ2(г, Б (г), I (г), щ(г)); 61 (г, Б (г), I(г), щ(г)) = Гх(г) - g2 (г, Б (г), I (г), щ(г), г);

^(г, Б (г), I (г), Щ(г)) = -2е-г -¿Б (г) I (г); й (г, Б (г), I (г), щ(г)) = - х(г) - Ъъ (г, Б (г), I (г), щ(г));

6(г, Б (г), I (г), щ(г)) = -Гх(г).

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

Работа выполнена в рамках проекта №»24С/2019 Хабаровского краевого фонда поддержки научных исследований.

ЛИТЕРАТУРА

[1] Дубко В.А. Первый интеграл системы стохастических дифференциальных уравнений. Киев, Институт математики АН УССР, 1978, 21 с.

[2] Еругин Н.П. Построение всего множества дифференциальных уравнений, имеющих заданную интегральную кривую. Прикладная математика и механика, 1952, т.16, № 6, с. 658-670.

[3] Тлеубергенов М.И. Об обратной задаче восстановления стохастических диффе-ренциальных систем. Дифференциальные уравнения, 2001, т. 37, № 5, с. 714-716.

[4] Дубко В.А. Вопросы теории и применения стохастических дифференциальных уравнений. Владивосток, ДВНЦ АН СССР, 1989, 185с.

[5] Карачанская Е.В. Построение программных управлений динамической системы на основе множества ее первых интегралов. Современная математика. Фундаментальные направления, 2011, № 12, с. 125-133.

[6] Чалых Е.В. Построение множества программных управлений с вероятностью 1 для одного класса стохастических систем. Автоматика и телемеханика, 2009, т. 70, № 8, с. 110-122.

[7] Карачанская Е.В. Построение множества дифференциальных уравнений с заданным множеством первых интегралов. Вестник Тихоокеанского государственного университета, 2011, № 3 (22), с. 47-56.

[8] Гихман И.И., Скороход А.В. Стохастические дифференциальные уравнения. Киев, Наукова Думка, 1968, 354 с.

[9] Карачанская Е.В. Обобщенная формула Ито - Вентцеля для случая нецентрированной пуассоновской меры, стохастический первый интеграл и первый интеграл. Математические труды, 2014, т. 17, № 1, с. 99-122.

[10] Averina T., Karachanskaya E., Rybakov K. Statistical analysis of diffusion systems with invariants. Russian journal of Numerical Analysis and Mathematical Modelling, 2018, vol. 33, pp. 1-13.

[11] Карачанская Е.В., Петрова А.П. Применение программного управления с вероятностью 1 для некоторых задач финансовой математики. Математические заметки СВФУ, 2018, т. 25, № 1, с. 25-38.

[12] Карачанская Е.В. Построение программных управлений с вероятностью 1 для динамической системы с пуассоновскими возмущениями. Вестник Тихоокеанского государственного университета, 2011, № 2 (21), с. 51-60.

[13] Карачанская Е.В. Доказательство обобщенной формулы Ито - Вентцеля с помощью дельта - функции и плотности нормального распределения. Математические заметки СВФУ, 2014, т. 21, № 3, с. 46-59.

[14] Кузнецов Д.Ф. Стохастические дифференциальные уравнения: теория и практика численного решения. Санкт-Петербург, Изд-во политехнического университета, 2010, 816 с.

[15] Kermack W.O., McKendrick A. Contributions to the Mathematical Theory of Epidemics. Proc. Royal Society, 1927, A 115, pp.700-721

[16] Романюха А.А. Математические модели в иммунологии и эпидемиологии инфекционных заболеваний. Москва, БИНОМ. Лаборатория знаний, 2015, 256 с.

Статья поступила в редакцию 21.09.2019

Ссылку на эту статью просим оформлять следующим образом: Карачанская Е.В. Моделирование систем дифференциальных уравнений с динамическими инвариантами. Математическое моделирование и численные методы, 2019, № 1, с. 98-117.

Карачанская Елена Викторовна — канд. физ.-мат. наук, доцент, доцент кафедры ин-формационных технологий и систем, Дальневосточный университет путей сообщения; доцент кафедры высшей математики Тихоокеанского госуниверситета. Область интересов: стохастические динамические системы с инвариантами. ORCID ID 0000-0003-0815-3688, Web of Science Researcher ID V-8176-2019. e-mail: [email protected]

Modeling of a differential equations system with dynamical invariants

© E.V. Karachanskaya1'2

1Far-Eastern State Transport University, Khabarovsk, 680021, Russia 2Pacific National University, Khabarovsk, 680035, Russia

This article is devoted to the construction method of differential equations system with smooth functions as invariants. This method allows us to construct both deterministic differential equations system and Ito's stochastic differential equations system. These stochastic differential equations are diffusion equations or jump-diffusion ones. The algorithm under consideration is based on the author's previous research and has no parallel. We study the realization of the algorithm in R2 and in R3 in detail. We represent a MathCad-

program for construction of every type of differential equations system described above. We also show some examples of differential equations system with the given invariant constructed automatically be the computer program. Relevance of our research is guaranteed by numerical solution of the created system of differential equations. These systems are solved using well-known numerical techniques: Euler method, Euler-Murayama method and Monte-Carlo procedure. The numerical solution of a system of differential equations and the corresponding values of the invariant function are represented in diagram form. Programmed control with probability 1 for a stochastic SIR-model is an application of the presented theory and the MathCad-program. Results of this research can be used to construct models of dynamical systems with invariant and to further explore these systems.

Keywords: first integral, differential equations system, stochastic differential equation, al-gorithmization for computer, MathCad.

REFERENCES

[1] Dubko V.A. Pervyj integral sistemy stokhasticheskikh differentsial'nykh uravnenij [A First integral for stochastic differential eqyations]. Kiev, Inctitut matematiki AN USSR Publ., 1978, 21 p.

[2] Erugin N.P. Prikladnaya matematika i mekhanika — Journal of Applied Mathematics andMechanics,1952, vol. 16, no. 6, pp. 658-670.

[3] Tleubergenov M.I. Differencial'nye uravneniya — Differential equation, 2001, vol. 37, no. 5, pp. 751-753.

[4] Dubko V.A. Voprosy teorii i primeneniya stohasticheskih differencial'nyh uravnenij [Questions on a theory and applications of stochastic differential equations]. Vladivostok, DVNC AN SSSR Publ., 1989, 185 p.

[5] Karachanskaya E.V. Construction of programmed controls for a dynamic system based on the set of its first integrals. Journal of Mathematical Sciences, 2014, vol. 199, no. 5, pp. 547-555.

[6] Chalykh E.V. Constructing the set of program controls with probability 1 for one class of stochastic systems. Automation and Remote Control, vol. 70, no. 8, pp. 1364-1375.

[7] Karachanskaya E.V. Vestnik Tihookeanskogo gosudarstvennogo universiteta — Proc. PNU, 2011, №3 (22), pp. 47-56.

[8] Gihman I. I., Skorohod A.V. Stochastic Differential Equations. Berlin, Heidelberg, Springer-Verlag, 1972, 354 p.

[9] Karachanskaya E. The generalized Ito-Venttsel' formula in the case of a noncen-tered Poisson measure, a stochastic first integral, and a first integral. Siberian Advances in Mathematics, 2015, vol. 25, iss. 3, pp 191-205.

[10] Averina T., Karachanskaya E., Rybakov K. Statistical analysis of diffusion systems with invariants. Russian journal of Numerical Analysis and Mathematical Modelling, 2018, vol. 33, pp. 1-13.

[11] Karachanskaya E.V., Petrova A.P. Matematicheskie zametki SVFU — Math. Notes of NEFU, 2018, vol. 25, no. 1, pp. 25-38.

[12] Karachanskaya E.V. Vestnik Tihookeanskogo gosudarstvennogo universiteta — Proc. PNU, 2011, №2 (21), pp. 51-60.

[13] Karachanskaya E.V. A Proof of the Generalized Ifo-Wentzell Formula via the Delta-Function and the Density of Normal Distribution. Yakutian mathematical journal, 2014, vol. 21, no. 3, pp. 40-51.

[14] Kuznetsov D.F. Stohasticheskie differencial'nye uravneniya: teoriya i praktika chislennogo resheniya [Stochastic differential equations: theory and practice of numerical calculation]. St. Petersburg, Peter the Great St. Petersburg Polytechnic University Publ., 2010, 816 p.

[15] Kermack W.O., McKendrick A.G. Contributions to the Mathematical Theory of Epidemies. Proc. Royal Society, A, 1927, vol. 115, pp.700-721.

[16] Romanyuha A.A. Matematicheskie modeli v immunologii i epidemiologii in-fekcionnyh zabolevanij [Mathematical Modeling of Immunosenescence: Scenarios, Processes and Limitations]. Moscow, BINOM. Laboratoriya znanij Publ,, 2015, 256 p.

Karachanskaya E.V., PhD in math., an associated professor of a department of informational technologies and systems, Far-Eastern State Transport University, an associated professor of high mathematics department of Pacific National University, Khabarovsk. Research area: stochastic dynamical systems with invariants. ORCID ID 0000-0003-08153688, Web of Science Researcher ID V-8176-2019. e-mail: [email protected]