Stability of the equilibrium position of a system of delay differential equations occurs in biological applications

Full Text

1. Введение и постановка задачи Дифференциальные уравнения с запаздыванием широко используются при построении математических моделей живых систем. В настоящее время имеется большое количество журнальных статей, обзоров и монографий, посвященных разработке и исследованию математических моделей динамики популяций в иммунологии, эпидемиологии, экологии и демографии на основе дифференциальных уравнений с запаздыванием. © Н.В. Перцев, М.В. Мулюков, А.Н. Пичугина, 2026 This work is licensed under a Creative Commons Attribution-NonCommercial 4.0 International License https://creativecommons.org/licenses/by-nc/4.0/legalcode 553 Важным аспектом изучения динамики некоторых популяций A1 и A2 является учет контактного взаимодействия индивидуумов. В результате осуществления контактного взаимодействия из пары индивидуумов популяций A1 и A2 образуется комплекс K, существующий некоторое конечное время. При распаде комплекса K индивидуум популяции A1 погибает, индивидуум популяции A2 либо погибает, либо возвращается обратно в эту популяцию. Типичным примером контактного взаимодействия индивидуумов являются межклеточные взаимодействия в иммунных процессах. Так, в частности, при моделировании врожденного иммунитета требуется построение уравнений, описывающих контактное взаимодействие зараженных клеток и NK-клеток (естественных клеток-киллеров). NK-клетки разрушают зараженные клетки в результате контактного цитолиза, длительность которого составляет 1-2 часа [10, с. 160]. При исследовании динамики ВИЧ-1 инфекции в лимфоузле необходимо учитывать контактные взаимодействия размножающихся CD4+ Т-лимфоцитов и зараженных CD4+ Т-лимфоцитов, производящих вирусные частицы. Контактное взаимодействие указанных клеток приводит к формированию неких «каналов» или «трубочек», по которым в еще не зараженную клетку переносятся отдельные компоненты, необходимые для сборки вирусных частиц. Длительность переноса компонент, необходимых для сборки вирусных частиц, может составлять от нескольких десятков минут до нескольких часов in vitro [11, с. 559]. Описанная выше специфика контактных взаимодействий индивидуумов требует привлечения соответствующего математического аппарата. В детерминированной постановке здесь могут быть использованы дифференциальные уравнения с распределенным запаздыванием. Рассмотрим один из вариантов указанной системы дифференциальных уравнений. Обозначим через x1(t), x2(t), y(t) численность популяций A1, A2 и численность комплексов K в момент времени t. Изучаемая система уравнений модели имеет вид (1.1) , (1.2) , (1.3) и дополняется начальными данными , (1.4) (1.5) Параметры модели имеют следующую интерпретацию. Популяции A1, A2 пополняются из внешних источников с постоянными скоростями ρ1 > 0, ρ2 > 0. Индивидуумы популяций A1 и A2 подвержены миграционному оттоку или гибели по естественными причинам с постоянными интенсивностями μ1 > 0, μ2 > 0. Контактное взаимодействие индивидуумов популяций A1, A2 осуществляется с постоянной интенсивностью γ > 0 в расчете на одну пару этих индивидуумов. Каждый образованный из пары индивидуумов комплекс K существует некоторое время η ∈ (0,ω], 0 < ω < ∞. Длительность η существования комплекса K задается функцией распределения R(s) на промежутке s ∈ [0,ω], такой, что R(0) = 0, R(ω) = 1. При распаде комплекса индивидуум популяции A1 погибает, тогда как индивидуум популяции A2 либо погибает (доля погибающих 1 - δ), либо возвращается обратно в эту популяцию (доля возвращающихся 0 < δ 1). При t = 0 под производными переменных x1(t), x2(t), y(t) понимаются их правосторонние производные, интеграл в (1.2), (1.3) понимается как интеграл Римана-Стилтьеса. В (1.4) начальные функции неотрицательны и непрерывны. В (1.5) начальное значение y(0) согласовано с формальным переходом от дифференциального уравнения (1.3) к интегральному уравнению для y(t) в форме свертки . (1.6) Уравнения (1.3), (1.6) и начальное условие (1.5) носят вспомогательный характер, а динамика популяций A1, A2 определяется уравнениями (1.1), (1.2) и начальным условием (1.4). Дополнительно отметим, что в уравнении (1.6) выражение интерпретируется как среднее время существования комплекса, состоящего из пары индивидуумов популяций A1, A2. В задачи настоящей работы входит: 1) исследование разрешимости задачи Коши для построенной системы уравнений; 2) нахождение неотрицательных положений равновесия системы уравнений модели и изучение их устойчивости; 3) анализ частных случаев системы уравнений модели (варианты модели, содержащие сосредоточенное запаздывание, и без запаздывания). 2. Глобальная разрешимость модели Модель (1.1), (1.2), (1.4) будем рассматривать как задачу Коши для систем функциональнодифференциальных уравнений с запаздыванием, исследованных в [7]. Решением задачи Коши (1.1), (1.2), (1.4) на конечном промежутке [0,τ) назовем непрерывные на промежутке [-ω,0] ∪ [0,τ) функции x1(t), x2(t), непрерывно дифференцируемые на промежутке [0,τ), удовлетворяющие начальному условию (1.4) и уравнениям (1.1), (1.2) для всех t ∈ [0,τ); при t = 0 под производными x1(t), x2(t) понимаем правосторонние производные этих переменных [3, с. 25]. Будем говорить, что задача Коши (1.1), (1.2), (1.4) является однозначно разрешимой на полуоси [0,∞) (глобально разрешимой), если эта задача имеет единственное решение на каждом конечном промежутке [0,τ). Предполагая неотрицательность и непрерывность x1(t), x2(t) на полуоси [0,∞), из (1.1), находим верхнюю оценку для x1(t) в форме константы: . (2.1) Используя (2.1), получаем верхнюю линейную оценку на слагаемые в (1.2), описывающие приток индивидуумов популяции A2: . (2.2) Опираясь на структуру уравнений (1.1), (1.2) и оценки (2.1), (2.2), устанавливаем, что для изучаемой задачи Коши выполнено сочетание предположений H2 и H3 из [7, с. 149]. Тогда, используя теорему 1 из [7, с. 153], получаем, что задача Коши (1.1), (1.2), (1.4) является глобально разрешимой и компоненты решения являются неотрицательными на полуоси [0,∞). 3. Положения равновесия и линеаризованная система уравнений Рассмотрим положения равновесия системы (1.1), (1.2) с неотрицательными компонентами. Учитывая, что , и полагая в системе (1.1), (1.2) x1(t) = x1 = const 0, x2(t) = x2 = const 0, приходим к системе ρ1 - μ1x1 - γx1x2 = 0, (3.1) ρ2 - μ2x2 - (1 - δ)γx1x2 = 0. Решение системы (3.1), (3.2) имеет вид (3.2) , (3.3) где - единственный положительный корень уравнения . (3.4) Заметим, что может быть выписан явно как положительный корень квадратного уравнения, вытекающего из (3.4) при δ < 1, или при δ = 1. От системы (1.1), (1.2) перейдем к системе линейного приближения в окрестности пары . Выполняя стандартные выкладки, введем переменные , для которых получаем систему уравнений (3.5) , (3.6) дополненную начальными данными . (3.7) Система (3.5), (3.6) имеет нулевое (тривиальное) положение равновесия , исследование устойчивости которого равносильно исследованию устойчивости положения равновесия системы (1.1), (1.2). Отбрасывая в (3.5), (3.6) нелинейные слагаемые и используя прежние обозначения переменных, получим линейную систему дифференциальных уравнений с распределенным запаздыванием (3.8) . (3.9) Обозначим . Опираясь на (3.8), (3.9), перейдем к линеаризованной системе (1.1), (1.2), записанной в форме . (3.10) 4. Устойчивость положения равновесия Для изучения устойчивости положения равновесия рассмотрим предварительно линейную систему (3.10), а потом перейдем к нелинейной системе (1.1), (1.2). 4.1. Анализ линейной системы. Обратимся к исследованию асимптотической устойчивости тривиального (нулевого) положения равновесия системы (3.10). Следуя [3, с. 50], частные решения системы (3.10) будем искать в форме x(t) = eztc, (4.1) где c = (c1,c2)T - постоянный вектор, - некоторое комплексное число. Подставляя в (3.10) выражение (4.1), приходим к характеристическому уравнению для z: , (4.2) где E - единичная матрица размера 2 × 2. Введем набор коэффициентов , . Отметим, что для коэффициентов (4.3) выполнены неравенства (4.3) q1 > q4 > 0, q2 > q3 > 0. Преобразуя (4.2), получаем характеристическое уравнение (4.4) . (4.5) Отметим, что система (3.10) является линейной автономной системой, в которой R(s) такова, что R(s) = 1 для всех s ω. Тогда функция, стоящая в левой части уравнения (4.5), является целой функцией переменной z, и уравнение (4.5) имеет либо конечное, либо бесконечное число корней z1, z2, ..., zk, ..., и Rezk → -∞ при k → ∞, см. [3, с. 51]. Опираясь на (4.3), (4.4), введем две константы: . (4.6) Лемма 4.1. Уравнение (4.5) не имеет корней z таких, что Rez 0. Доказательство. Пусть z = b+iϕ- предполагаемый корень уравнения (4.5) при. Подставляя z в (4.5), получаем: . (4.7) Из (4.7) следует, что От (4.8) приходим к неравенству . (4.9) Обозначим βb = (b2 + q1b + q2)2 - (q3 + q4b)2 = (b2 + (q1 - q4)b + q2 - q3)(b2 + (q1 + q4)b + q2 + q3) > 0, αb = (2b + q1)2 - 2(b2 + q1b + q2) - q42 = 2b(b + q1) + α0 > 0, где α0 > 0 приведено в (4.6). Используя указанные обозначения, запишем (4.9) в следующей форме: . (4.10) Полагая θ = ϕ2, получаем неравенство . (4.11) Поскольку в (4.11) коэффициенты αb и βb положительны, то θ либо не вещественно, либо отрицательно. Отсюда следует, что для каждого фиксированного неравенство (4.10) неразрешимо относительно ϕ. Приходим к противоречию, что завершает доказательство леммы. Лемма 4.2. Существует b∗ < 0 такое, что уравнение (4.5) не имеет корней z, удовлетворяющих неравенствам b∗ < Rez < 0. Доказательство. Пусть z = b + iϕ -некоторый корень уравнения (4.5) при b < 0. Подставляя z в (4.5), получаем: . (4.12) Из (4.12) следует, что От (4.13) приходим к неравенству . (4.14) По аналогии с (4.9), (4.10), перепишем (4.14) в форме , (4.15) где βb(+) = (b2 + q1b + q2)2 - e2|b|ω(q3 + q4b)2 = = (b2 + (q1 - e|b|ωq4)b + q2 - e|b|ωq3)(b2 + (q1 + e|b|ωq4)b + q2 + e|b|ωq3), α(+)b = (2b + q1)2 - 2(b2 + q1b + q2) - e|b|ωq42 = 2b(b + q1) + α0 + q42(1 - e2|b|ω). Пусть b → 0, b < 0. Тогда βb(+) → β0 > 0, α(+)b → α0 > 0, (4.16) где β0, α0 приведены в (4.6). По аналогии с леммой 4.1 получаем, что если для некоторого b < 0 верно βb(+) > 0, α(+)b > 0, то уравнение (4.15) неразрешимо относительно ϕ. В силу непрерывности функций и выполнения неравенств (4.16) получаем, что существует промежуток (b∗,0), b∗ < 0 такой, что уравнение (4.5) не имеет корней z таких, что Rez ∈ (b∗,0). Последнее утверждение завершает доказательство. Исходя из лемм 4.1, 4.2, устанавливаем, что все корни характеристического уравнения (4.5) таковы, что их вещественные части отрицательны и отделены от нуля. В соответствии с [3, с. 50- 54], приходим к следующему результату. Теорема 4.1. Тривиальное (нулевое) решение линейной системы (3.10) асимптотически устойчиво. 4.2. Анализ нелинейной системы. Перейдем к изучению устойчивости положения равновесия, используя систему (3.5), (3.6) с начальным условием (3.7). Применим теорему 33.2 из [4, с. 194] об устойчивости по первому приближению (см. также ссылки на указанную теорему в [3, с. 106-107] и в [2]). Опираясь на приведенные выше леммы 4.1, 4.2 и теорему 4.1, имеем, что линейная часть системы (3.5), (3.6) в форме системы (3.10) удовлетворяет первому из условий теоремы 33.2 из [4, с. 194]. Учитывая, что изучается локальная асимптотическая устойчивость положения равновесия , построим оценки нелинейных слагаемых из системы (3.5), (3.6) в достаточно малой окрестности нуля. Для дальнейшего изложения введем обозначение, используемое в [4, с. 157]. Пусть - значение решения системы (3.5), (3.6) с начальным условием (3.7) в некоторый момент времени t. Положим, что-отрезок решения, заданный соотношением . Отметим, что для xt можно было бы использовать и другую общепринятую запись, а именно: (см., например, [3, с. 87] и [8, с. 50]). Но тогда интегральные слагаемые в уравнениях (1.2), (1.3), (1.6) изучаемой модели следовало бы представлять немного в другой форме. Нелинейные слагаемые из системы (3.5), (3.6) обозначим через . Используем далее известные соотношения: если u, v - две вещественные переменные, то , | | | Тогда для достаточно малых |u| + |v| имеет место цепочка неравенств , (4.17) где α1, α2 - некоторые (достаточно малые) положительные константы. Следуя [4, с. 8, 151], положим . (4.18) Опираясь на (4.17), (4.18), оценим , предполагая малость , включая соответствующие выражения и для начальных функций (3.7). Имеем, что , (4.19) где константа β > 0 может быть выбрана достаточно малой. Оценка (4.19) обеспечивает выполнение второго из условий теоремы 33.2 из [4, с. 194] в окрестности нулевого положения равновесия системы (3.5), (3.6). В итоге приходим к окончательному результату. Теорема 4.2. Положение равновесия системы (1.1), (1.2) является асимптотически устойчивым. 5. Случай сосредоточенного запаздывания В приложениях обычно рассматриваются модели (1.1)-(1.5) в относительно простой форме. Одна из таких моделей для (1.1)-(1.5) представляет собой систему дифференциальных уравнений с запаздыванием, сосредоточенным в точке s = ω > 0: (5.1) , (5.2) (5.3) (5.4) (5.5) Уравнение (5.3) с учетом (5.5) сводится к интегральному . Система (5.1), (5.2) имеет положение равновесия , указанное в разделе 3. Система уравнений линейного приближения записывается в форме , (5.6) где матрица A приведена в разделе 3, матрица B такова: . Получаем, что характеристическое уравнение для системы (5.6) z2 + q1z + q2 - (q3 + q4z)e-ωz = 0 (5.7) представляет собой частный случай уравнения (4.5). Система вида (5.6), содержащая матрицы A и B с произвольными вещественными элементами, впервые изучалась в [1]. Современное исследование корней характеристического уравнения (5.7) для произвольных вещественных матриц A и B из (5.6) выполнено, например, в [5, 12]. Здесь коэффициенты q1,...,q4 уравнения (5.7) зависят от восьми элементов матриц A и B, что приводит к достаточно сложным соотношениям, обеспечивающим условия асимптотической устойчивости или неустойчивости тривиального (нулевого) положения равновесия системы (5.6). В настоящем исследовании существенно важно, что для системы (5.6) с матрицами A, B конкретного вида используются связи между элементами этих матриц, возникающие из изучаемой модели. Как показано выше, указанные связи между элементами матриц A, B (соотношения (4.4)) обеспечивают отрицательность вещественных частей корней уравнения (5.7) и их отделимость от нуля. Нетрудно заметить, что отброшенные при линеаризации системы (5.1), (5.2) слагаемые представляют собой квадратичные функции и имеют порядок выше первого относительно совокупности переменных , . Используя леммы 4.1, 4.2 и опираясь на [9, с. 148-149], приходим к следующему результату. Теорема 5.1. Положение равновесия системы (5.1), (5.2) является асимптотически устойчивым. 6. Дополнительный вариант модели В качестве простейшего дополнительного варианта модели приведем систему уравнений без запаздывания: (6.1) , (6.2) (6.3) (6.4) В уравнениях (6.2), (6.3) параметр r > 0 задает интенсивность распада комплексов K, состоящих из пар индивидуумов. Начальные данные (6.4) содержат неотрицательные константы. По аналогии с разделом 2 нетрудно показать, что задача Коши (6.1)-(6.4) является глобально разрешимой и компоненты ее решения являются неотрицательными на полуоси [0,∞). Система (6.1)-(6.3) имеет положение равновесия , где , а компоненты указаны в разделе 3. Переходя от (6.1)-(6.3) к системе уравнений линейного приближения и применяя к ней критерий Рауса-Гурвица, устанавливаем асимптотическую устойчивость положения равновесия. Параметр r, используемый в (6.2), (6.3), можно интерпретировать как параметр экспоненциального распределения длительности η существования комплекса, состоящего из пары индивидуумов, полагая, что R(s) заменяется функцией F(s) = 1 - exp(-rs), s ∈ [0,∞). (6.5) Примем, что в (6.4) y(0) = 0 и выполним ряд формальных преобразований. Интегрируя (6.3), находим, что . (6.6) Подставляя (6.6) в (6.2), (6.3) и учитывая (6.5), приходим к системе (6.7) , (6.8) (6.9) дополненной начальными данными . (6.10) Получаем, что для начальных данных (6.10) система (6.1)-(6.3) может быть записана в форме (6.6)-(6.9), близкой к исходной системе (1.1)-(1.3). В случае y(0) > 0 требуется учитывать дополнительные слагаемые в уравнениях (6.8), (6.9). Один из возможных способов построения изучаемых систем уравнений приведен в [6]. Если для системы (6.7)-(6.9) вместо (6.5) использовать произвольную функцию распределения F0(s), отличную от R(s), s ∈ [0,∞), то возникает еще один вариант модели, учитывающий неограниченную длительность η существования комплекса, состоящего из пары индивидуумов популяций A1, A2. Система (6.7), (6.8) с функцией F(s) ≡ F0(s) требует отдельного исследования. С биологической точки зрения интерес представляет случай, означающий, что среднее время существования комплекса является конечным.

About the authors

N. V. Pertsev

Sobolev Institute of Mathematics SB RAS

Email: pertsevnv@mail.ru
ORCID iD: 0000-0002-5674-7096
SPIN-code: 9626-6689
Scopus Author ID: 8964903300
ResearcherId: F-8086-2013
Novosibirsk, Russia

M. V. Mulyukov

Perm State National Research University; Perm National Research Polytechnic University

Email: mulykoff@gmail.com
ORCID iD: 0000-0001-5610-8035
SPIN-code: 3449-7824
Scopus Author ID: 56040106000
ResearcherId: AAQ-1475-2021
Perm, Russia

A. N. Pichugina

Saint Petersburg Electrotechnical University «LETI»

Author for correspondence.
Email: anna.pichugina@gmail.com
SPIN-code: 4196-6607
Scopus Author ID: 55653900600
ResearcherId: W-5290-2018
Saint-Petersburg, Russia

References

Supplementary files