Анализ устойчивости положений равновесия одной модели популяционной динамики Текст научной статьи по специальности «Математика»

МАТЕМАТИКА

Вестн. Ом. ун-та. 2009. № 2. С. 50-53.

УДК 519

Н.В. Перцев, Г.Е. Царегородцева

Омский филиал Института математики им. С.Л. Соболева СО РАН

А.Н. Пичугина

Омский государственный университет им. Ф. М. Достоевского

АНАЛИЗ УСТОЙЧИВОСТИ ПОЛОЖЕНИЙ РАВНОВЕСИЯ ОДНОЙ МОДЕЛИ ПОПУЛЯЦИОННОЙ ДИНАМИКИ*

В статье дается анализ устойчивости положений равновесия модели популяционной динамики.

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

1. Уравнения модели

Будем рассматривать популяцию индивидуумов, динамика которой определяется следующими факторами: 1) индивидуумы дают потомство, при этом репродуктивный потенциал может снижаться за счет воздействия вредных веществ; 2) индивидуумы погибают за счет самолимитирования и конкуренции, а также вследствие потребления вредных веществ, поступающих в составе ресурсов питания. Обозначим через Xj (t) численность популяции в момент времени t > 0 . Предположим, что изменение количества индивидуумов за бесконечно малый промежуток времени (t, t + dt), t > 0, dt > 0 определяется балансовыми соотношениями между процессами рождения и гибели индивидуумов, интенсивность которых зависит от Xj(t) и X2(t), где Х2 = x2 (t) - количество вредных веществ в момент времени t > 0 . Считаем, что индивидуумы нулевого возраста производятся за счет существующих индивидуумов независимо от их возраста и не поступают из внешнего источника. Примем, что количество индивидуумов нулевого возраста, появившихся за интервал времени (t, t + dt), равно b( x 2(t)) xj(t )dt, где b( x2) = eg (0( x 2)), в = const > 0 описывает максимальную скорость репродукции одного индивидуума без учета влияния вредных веществ. Функция g(u) учитывает снижение репродуктивного потенциала индивидуумов за счет неблагоприятного воздействия вредных веществ. Функция 6(x2) отражает скорость потребления вредных веществ в составе пищевых ресурсов. Примем, что количество индивидуумов, погибающих за счет самолимитирования и конкуренции за время (t, t + dt), dt ^ +0, независимо от их воз-

*Работа поддержана РФФИ, проект № 09-01-00098-а.

© Н.В. Перцев, Г.Е. Царегородцева, А.Н. Пичугина, 2009

раста равно yxf (t)dt, у = const > 0 - интенсивность взаимодействия индивидуумов. Потребление вредных веществ в составе ресурсов питания приводит к гибели индивидуумов за время (t, t + dt),

dt ^+0 в количестве G0(x2(t))x1(t)dt,

где G = const > 0 - коэффициент, отражающий степень токсичности вредного вещества.

Используя сделанные предположения, запишем, что

x1 (t) = b(x2 (t))x1 (t) - Yl (t) - G0(x2 (t))x1 (t) .

Примем, как и в работе [1], что количество вредных веществ в окружающей среде изменяется по правилу

x2 (t) = p - 0(x2 (t))x1 (t) - 8xi (t) , где p = const > 0 описывает скорость поступления вредных веществ, 8 = const > 0

- интенсивность распада вредных веществ.

Объединяя все полученные соотношения и добавляя начальные данные, получаем, что уравнения модели имеют следующий вид:

xi = Pg(0( x2 )) x1 - Yl -G0( x2) x1 , x2 = p -0(x2)x1 -8x2 , I1)

xi (0) = x° > 0, i = 1,2. (2)

Полагаем, что функция g = g(u) при u > 0 - непрерывная, невозрастающая, g(0) = 1, g(u) ^ 0, u ^ +ro, имеет непрерывную производную g' (u). Обозначим эти условия (g*). Функция 0(X2) -непрерывная, неубывающая, 0(0) = 0,

имеющая конечный предел 0 при Х2 ^ +го и непрерывную производную 0'(Х2). Условия на функцию 0 будем обозначать (0*).

2. Исследование свойств решений модели

Рассмотрим систему (1), (2). Предположим, что выполнены условия (g*), (0*). Система (1), (2) является автономной. Опираясь на теорему существования и единственности для нормальной системы уравнений [2], получаем, что (1), (2) имеет единственное решение на некотором промежутке 0 < t < t0, 0 < 10 <<Х). Ис-

пользуя схему доказательства из [1], можно показать, что решение (1), (2) существует, единственно и неотрицательно на промежутке 0 < t < да . Для анализа поведения решений системы (1), (2) исследуем устойчивость положений равновесия * / * * \

x = (х^, x2) с неотрицательными компонентами, опираясь на матрицу системы линейного приближения в точке x *

(

A =

pg(0) - 27х? - G0 х0 '(pg'(0) - g)

Л

-0

-0 X * -8

v1

и ее собственные числа, det(A — AE) = 0 . СЛУЧАЙ 1.

Положение равновесия:

х* = 0 х* = p/8

(3)

Здесь явно получаем, что

л,=-г< о, ъ =в еО)) -сеО).

О о

Положение равновесия (3) асимптоти-

Р Р

чески устойчиво, если в(е(—)) < се(-).

ОО

При

выполнении

неравенства

Р Р

> се(—) оно неустойчиво. В

ОО случае Л2 = 0 требуются дополнительные исследования.

СЛУЧАЙ 2.

Точка положения равновесия:

* * х1 > 0, х2 > 0 . (4)

Для нахождения положений равновесия (4) решаем систему уравнений:

(pg(0(Х2 )) -Yci - G0(Х2 ))Х1 = 0 Р -0(Х2)Х1 8Х2 = 0.

(5)

Выражаем хх > 0 из системы (5), приходим к уравнению Е(х2) = @(х2) , где

Е(Х2) = Р -0х2 ,

0(х2) = У~е(х2)(в(е(х2)) -се(х2)) , х2 > 0.

Покажем далее, что характер пересечения графиков функцийЕ(х2) , Q(х2) определяет асимптотическую устойчивость (неустойчивость) соответствующих точек пересечений, которые задают положения равновесия вида (4).

52

Н.В. Перцев, Г.Е. Царегородцева, А.Н. Пичугина

F

Рис. 1. Типичный вид функций Е(Х2) и Q(Х2)

На рисунке 1 показан пример графиков функций Е(Х2), Q(Х2), которые пересекаются в нескольких точках. Это оз*

начает, что Х вида (4) может существовать, причем не одно. Пусть положение равновесия вида (4) существует. Исследуем его на устойчивость. Составляем матрицу А в точке (4). Для нахождения ее собственных чисел решаем уравнение

Л + аЛ + Ь = 0, где а = ус* + в'(х* )х* + 5 ,

Ь = х*(у(в X х*) Х* +5) +

+ в( х*)в'( х*)(в '(в( х*))-а)).

Для определения знаков ReЛг. применим критерий Рауса-Гурвица для уравнений второго порядка. Получаем, что Re Л2 < 0 О а > 0, Ь > 0 . В рассматриваемом случае а > 0 . Для исследования знака Ь используем расположение графиков функций Е (Х2) и ^ Х2) и геометрический смысл производной, поскольку видно, что Ь = ухх ^'(Х2) — Е'(Х2)) . Как следствие, замечаем, что среди точек пересечения графиков функций Е (Х2) и Q( Х2), приведенных на рисунке 1, положения равновесия, соответствующие точкам А, С, Е, будут асимптотически устойчивы, положения равновесия В, О, Н - неустойчивы.

Покажем, что в положительной части фазовой плоскости для системы (1), (2) выполняются условия критерия Дюлака отсутствия предельных циклов [3]. Обозначим

к(х1 , х2) = в(в(Х2))х1 — ^с2 — ав(Х2)Х1,

Ь(х1, х2) = р — в(х2)х1 — 5с2.

В качестве множителя B(x1, x2), фигурирующего в условии критерия, возьмем B(x1, x2) = X-1. Тогда

5

— (B(X^ X2)K(X1, X2)) +

dx1

d

+ — (B(X1 , X2 )L(X1 , X2 )) = N(X1 , X2 ) ,

dX2

где

N(X1, x2 ) = -y- gQ(x2 ) - Q'(x2 ) - ¿X-1 < 0, при X1 > 0 , x2 > 0 . Следовательно, в положительной части фазовой плоскости не существует замкнутых контуров, составленных из траекторий системы (1), (2) (нет предельных циклов).

3. Вычислительный эксперимент Целью вычислительного эксперимента являлось выявление наиболее характерных режимов поведения решений системы (1), (2). Для проведения расчетов использована программа Maple, библиотека DEtools и функция DEplot. В процессе вычисления варьировались параметры, отвечающие за наличие и устойчивость положений равновесия.

ПРИМЕР 1.

Подбираем функции и параметры таким образом, чтобы система (1), (2) имела ровно одно ненулевое положение равновесия:

Q(x2) = (2exp(x2 / 20) - 2)(1 + exp(x2 / 20))-1,

g (u) = 0.5(exp-1 (u) + (u +1)-1),

P = 2, у = 0.002, g = 0.1, p = 16, 8 = 2.

При выбранном наборе параметров система (1), (2) имеет ненулевое положение равновесия (x1 = 984, x2 = 0.3125) ,

которое является асимптотически устойчивым, кроме того, существует положение

равновесия (x1 = 0, x2 = 8) , которое является неустойчивым. Результаты численного решения модели при различных начальных данных

n1: (x1 (0) = 400, x2 (0) = 6),

n2: (x1 (0) = 90, x2(0) = 0.3),

n3: (x1 (0) = 25, x2 (0) = 9) представлены на рисунке 2.

О 2 4 В 8

Рис. 2. Динамика численности популяции в зависимости от начальных данных

ПРИМЕР 2.

Рассмотрим случай, в котором существует устойчивое положение равновесия (3). Подберем функции и параметры таким образом, чтобы система (1), (2) имела три положения равновесия:

в(х2 ) = (3 ехр(х2 / 20) — 3)(1 + ехр(х2 / 20))—, Е (и) = 0.5(ехр—1 (и) + (и +1)—1), в = 2, у = 0.002, а = 8, р = 20, 5 = 4. Положения равновесия ( х1 = 894,

х2 = 0.281, (х1 = 0, х2 = 5)) - асимптотически устойчивы.

( х* = 49.076, х2* = 2.606) - неустойчиво. Исследуем поведение решений при различных начальных данных:

п1: ( х1 (0) = 49.076, х2 (0) = 2.606), п2: (х1 (0) = 49.076, х2 (0) = 2.607) .

Рис. 3. Динамика численности в случае двух устойчивых положений равновесия

Из рисунка 3 видно, что небольшое отличие в начальных данных приводит к совершенно разным результатам. В случае п1 наблюдается выход популяции на ненулевой стационарный уровень, тогда как в п2 - вырождение популяции.

ПРИМЕР 3.

Приведем набор функций и параметров модели, для которых система (1), (2) имеет пять положений равновесия. Пусть в( х2 ) = тт(2(вт( х2 ) + х2 ),20)

Е (и) = 0.5(ехр—1 (и / 20) + (и +1)—1),

в = 2, у = 0.02, а = 0.1, р = 110, 5 = 11.

На рисунке 4 изображено расположение графиков функций Е (Х2) , Q(Х2), которые пересекаются в четырех точках. Как отмечено выше, точки пересечения этих графиков и начальные данные определяют пове-деиие решений модели при t —> +00 .

Рис. 4. Графики функций F(x2) и Q(x2)

4. Заключение

Результаты исследования указывают на достаточно сложную динамику популяции, развивающуюся в условиях воздействия вредного вещества. Эффект воздействия этого вещества проявляется в существовании нескольких устойчивых стационарных режимов динамики популяции, одним из которых может являться режим, связанный с вырождением популяции. Предложенная модель позволяет получить условия, при которых возможно снижение численности популяции до наименьшего положительного уровня или ее полное вырождение.

ЛИТЕРАТУРА

[1] Пичугина А.Н. Интегродифференциальная мо-

дель популяции, подверженной воздействию вредных веществ // Сибирский журнал индустриальной математики. Т. 7. № 4 (20). 2004. С. 130-140.

[2] Понтрягин Л.С. Обыкновенные дифференциальные уравнения. М.: Наука, 1974.

[3] Баутин Н.Н., Леонтович Е.А. Методы и приемы

качественного исследования динамики систем на плоскости. М.: Наука, 1976.