Математическое моделирование производственно-экономических процессов и систем

Применения теории графов в сетевом планировании и управлении

Выбираем вариант задания:

Номер варианта определяется по формуле:

N+1–25 N/25 . N – последние 3 цифры номера зачётной книжки. Квад-ратные скобки обозначают целую часть числа, заключённого в эти скобки.

В данном случае: 479+1–25 479/25 =5.

ОСНОВНЫЕ ЭТАПЫ ВЫПОЛНЕНИЯ РАБОТЫ

1 ЭТАП

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

Сетевое планирование и управление программами включает три основ-ные этапа: структурное планирование, календарное планирование и оперативное управление.

Этап структурного планирования начинается с разбиения программы на чётко определённые операции. Затем определяются оценки продолжительности операции, и строится сетевая модель, сетевой график, стрелочная диаграмма, ка-ждая дуга (стрелка) которой, отображает работу.

Вся сетевая модель в целом является графическим представлением взаимосвязей операций программы. Построение сетевой модели на этапе струк-турного планирования позволяет детально проанализировать все операции и вне-сти улучшения в структуру программы ещё до начала её реализации. Однако ещё более существенную роль играет использование сетевой модели для разработки календарного плана выполнения программы.

Конечной целью этапа календарного планирования является построение календарного графика, определяющего моменты начала и окончания каждой опе-рации, а также её взаимосвязи с другими операциями программы. Кроме того, ка-лендарный график должен давать возможность выявлять критические операции (с точки зрения времени), которым необходимо уделять особое внимание, чтобы за-кончить программу в директивный срок. Что касается некритических операций, то календарный план должен позволять определять их резервы времени, которые можно выгодно использовать при задержке выполнения таких операций или с по-зиции эффективного использования ресурсов.

Заключительным этапом является оперативное управление процессом реализации программы. Этот этап включает использование сетевой модели и ка-лендарного графика для составления периодических отчётов о ходе выполнения программы. Сетевая модель подвергается анализу и в случае необходимости кор-ректируется. В этом случае разрабатывается новый календарный план выполне-ния остальной части программы. В связи с этим потребовалась разработка эффек-тивных методов планирования, обеспечивающих оптимизацию всего процесса осуществления проекта. Один из них, метод критического пути /МКП/, рассмат-ривается ниже.

Каждый проект можно представить в виде некоторого графа. Для этого необходимо знать: 1.

Перечень всех операций проекта. 2.

Время, необходимое для выполнения каждой операции. 3.

Перечень операций, непосредственно предшествующей каждой операции проекта.

Для представления операции используется ориентированная дуга, на-правление которой соответствует процессу реализации проекта во времени. От-ношение упорядочения между операциями задается с помощью событий. Собы-тие определяется как момент времени, когда завершаются одни операции и начи-наются другие. Начальная и конечная точки, любой операции описываются парой событий, которые обычно называют начальным и конечным событием. Операции, выходящие из некоторого события, не могут начаться, пока не будут завершены все операции, входящие в это событие. Таким образом, каждая операция пред-ставляется ориентированной дугой, а каждое событие – узлом /вершиной/. Граф, изображающий отношения предшествования между операциями проекта, называ-ется сетевым графиком.

(а) (б)

Рис. 1

На рис. 1 /а/ приведён типичный пример графического изображения операции /i j/ c начальным событием i и конечным событием j. Из рис. 1 /б/ видно, что для возможности начала операции /3.4/ требуется завершение операций /1.3/ и /2.3/.

Заметим, что в случае, когда две или большее число операций допусти-мо выполнять одновременно ПОЯВЛЯЕТСЯ неоднозначность определения опе-раций через события. Чтобы исключить такую неоднозначность вводят фиктив-ные операции, которые не требуют затрат ни времени, ни ресурсов.

Пример такого случая приведён на рис. 2 /а/ где операции А и В имеют одинаковые начальные и конечные события. Чтобы исключить такую неодно-значность вводится фиктивная операция /D/ обозначается пунктиром /рис. 2б/. Фиктивные операции позволяют также правильно отображать логические связи, которые без их помощи нельзя задать на сети. Предположим что в некотором проекте операции А и В должны непосредственно предшествовать операции С а операция Е непосредственно предшествует только В. Представление указанных условий даёт рис. 3 в котором используется фиктивная операция.

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

Алгоритм нумерации событий

Шаг 1. Событию, в которое не входит ни одной дуги присваивается но-мер 0.

Шаг 2. Любому неперенумерованному событию, для которого все предшествующие события перенумерованы, присваивается следующий номер.

Повторять шаг 2 до тех пор, пока все события не будут перенумерова-ны.

Таблица 1 Номера ра-бот Каким работам предше-ствует Продолжительность работ Потребность в трудовых ресурсах 1 3,4,5 6 1 2 6 8 0 3 10 10 2 4 7,8,9 6 1 5 6 8 0 6 7,8,9 10 5 7 10 6 1 8 11,12, 10 3 9 14 12 7 10 11,12 8 8 11 13 6 3 12 14 14 5 13 15 8 8 14 16 4 4 15 6 7 16 8 1

Соответствующий этому примеру сетевой график приведён на рис. 4.

Рис. 4.

Обозначим сетевой график проекта через G (X, A), где:

Х= 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 –множество событий,

A= /0,1/, /0,2/, /1,2/, /1,3/, /1,4/, /2,3/, /3,4/, /3,5/, /3,6/, /4,5/, /5,6/, /5,7/, /6,8/, /7,9/, /8,10/, /9,10/ – множество операций.

2 ЭТАП

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

Расчёт критического пути включает два этапа. Первый этап называется прямым проходом, а второй - обратным проходом. Для каждого события на пер-вом этапе вычисляется число, представляющее наиболее ранний возможный срок его наступления, а на втором, число, представляющее наиболее поздний допус-тимый срок его наступления.

Обозначим через Е /j/ - наиболее ранний возможный срок наступления j – го события. Пусть d I j – продолжительность операции, соединяющей i – е и j- е события. Поскольку j – е событие не может произойти пока не будут завершены все операции, ведущие к j – му узлу, то наиболее ранний возможный срок его на-ступления будет вычисляться по следующему алгоритму.

Алгоритм расчёта наиболее ранних возможных сроков наступления со-бытий

Шаг 1. Полагается Е/0/=0.

Шаг 2. Для j = 1. 2. . . . n вычисляется:

E ( j ) = max E ( i ) + d i j

i: ( i j ) э А

где максимум берётся по всем операциям, завершающимся в j – м узле и выходящим из любого предшествующего i - го узла.

Обозначим теперь через L (i) наиболее поздний срок наступления i – го события, не влияющий на время завершения всего проекта. Начиная с завершаю-щего события движемся в обратном направлении через каждое предшествующее событие. Вычисления осуществляются в этом случае по следующему алгоритму.

Алгоритм расчёта наиболее поздних допустимых сроков наступления событий.

Шаг 1. Полагается L(n)=E (n).

Шаг 2. Для i = n – 1, n – 2, …, 0 вычисляется:

L ( i ) = max L ( j – d i j )

j: ( i j ) э А

где минимум берётся по всем операциям, начинающимся в i – ом узле и входящим в любой j – й узел.

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

Операция считается критической, если задержка её начала приводит к увеличению срока всего проекта. Критическая операция /i, j/ удовлетворяет сле-дующим трём условиям:

E (i) = L (i)

E (j) = L (j)

E (j) – E (i) = L (j) – L (i) = d ij.

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

Для нашего примера расчёт выглядит следующим образом:

E(0)=0,

E(1)=max 0+6 =6;

E(2)=max 0+8,6+8 =14;

E(3)=max 6+6,14+10 =24;

E(4)=max 6+10,24+6 =30;

E(5)=max 30+8,24+10 =38;

E(6)=max 24+12,38+14 =52;

E(7)=max 38+6 =44;

E(8)=max 52+4 =56;

E(9)=max 44+8 =52;

E(10)=max 52+6,56+8 =64.

L(i)=max L(j–d i j) ;

L(10)=E(10)=64;

L(9)=max 64-6 =58;

L(8)=max 64-8 =56;

L(7)=max 58-8 =50;

L(6)=max 56-4 =52;

L(5)=max 50-6,52-14 =44;

L(4)=max 44-8 =36;

L(3)=max 52-12,44-10,36-6 =40;

L(2)=max 40-10 =30;

L(1)=max 36-10,40-6,30-8 =34;

L(0)=max 34-6,30-8 =28.

Критический путь включает операции /0,1/, /1,2/, /2,3/, /3,4/, /4,5/, /5,6/, /6,8/, /8,10/.

Минимальное время реализации проекта равно 64.

3 ЭТАП

Теперь необходимо вычислить резервы времени для некритических операций. Очевидно, что резерв времени критической операции должен быть ра-вен нулю.

Рассмотрим некоторую некритическую операцию /i, j/. Какое макси-мальное количество времени можно выделить для её выполнения без задержки своевременного окончания выполнения этого проекта? Операция /i, j/ может на-чаться не ранее Е /i/ и должна закончиться не позднее L (j). Таким образом, без задержки окончания проекта на выполнение операции /i, j/ можно выделить не более L (j) – E (i) единиц времени. Следовательно, при выполнении этой операции можно допустить максимальную задержку L (j) – E (i) – d ij > 0. Величина L (j) – E (i) – d ij называется полным резервом времени операции (i, j).

Какое максимальное количество времени может быть выделено для вы-полнения операции (i, j) без введения дополнительных временных ограничений на последующие операции? Для соблюдения этого условия операция (i, j) должна быть закончена к моменту времени E (j).

Поскольку операция (i, j) может начаться не ранее E (i), на её выполне-ние без введения дополнительных ограничений на последующие операции модно выделять не более E (j) – E (i) единиц времени.

Величина E (j) – E (i) - d ij называется свободным резервом времени операции (i, j).

Свободный резерв времени равен максимальной задержке выполнения операции (i, j), не влияющий на выполнение последующих операций.

Какое максимальное количество времени может быть выделено для вы-полнения операции (i, j) без введения дополнительных временных ограничений на любую операцию проекта? Для выполнения этого условия операция (i, j) должна начаться в момент времени L (i) и закончиться к моменту времени E (j), следовательно, на выполнение операции (i, j) в этом случае можно выделить не более E (j) – L (i) единиц времени.

Величина E (j) – L (i) - d ij называется независимым резервом времени операции (i, j). Независимый резерв времени равен максимальной задержке, кото-рую. Можно допустить при выполнении операции (i, j) без введения дополни-тельных временных ограничений на любую другую операцию проекта. Отрица-тельное значение независимого резерва означает, что любая задержка с выполне-нием операции приведёт к дополнительным ограничениям на выполнение других операций.

Для нашего примера, рассмотренные резервы времени операции проек-та приведены в таблице 2.

Таблица 2

Операция i j d ij E (i) L (i) E (j) L (j) Полный резерв Свободный резерв Независимый резерв 1 0 1 6 0 28 6 34 28 0 -28 2 0 2 8 0 28 14 30 22 6 -22 3 1 4 10 6 34 30 36 20 14 -14 4 1 3 6 6 34 24 40 28 12 -16 5 1 2 8 6 34 14 30 16 0 -28 6 2 3 10 14 30 24 40 16 0 -16 7 3 4 6 24 40 30 36 6 0 -16 8 3 5 10 24 40 38 44 10 4 -12 9 3 6 12 24 40 52 52 16 16 0 10 4 5 8 30 36 38 44 6 0 -6 11 5 7 6 38 44 44 50 6 0 -6 12 5 6 14 38 44 52 52 0 0 -6 13 7 9 8 44 50 52 58 6 0 -6 14 6 8 4 52 52 56 56 0 0 0 15 9 10 6 52 58 64 64 6 6 0 16 8 10 8 56 56 64 64 0 0 0

4 ЭТАП

Конечным результатом выполняемых расчётов является календарный график.

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

Ограничимся в работе определением потребностей в трудовых ресурсах для каждого единичного интервала времени лишь в двух крайних случаях.

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

Процедура решения этих задач состоит в следующем. Для каждой опе-рации (i, j) в качестве момента её начала выбирается E (i) в первом случае и L (i) во втором. Затем распределяются все ресурсы, требуемые для её выполнения. По-сле окончания распределения ресурсов устанавливается их ежедневная суммарная потребность и вычёркивается график ежедневной потребности в ресурсах, назы-ваемых ресурсным профилем проекта. В заключение необходимо сравнить мак-симальные потребности в ресурсах для обоих случаев.

Для нашего примера на рис. 5 приведены: календарный график реализа-ции проекта /рис. 5а/ и два графика ресурсных профилей проекта для описанных крайних случаев /рис. 5б и 5в/.

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

РИС. 5 а.

РИС. 5 б.

РИС. 5 в.

На рис. 5 б. показана потребность в рабочей силе при условии выбора в качестве календарных сроков некритических операций начала их ранних сроков (так называемый ранний, или левый календарный план), а на рис. 5 в - потреб-ность в силе при выборе наиболее поздних сроков (так называемый поздний, или правый, календарный план). Пунктирной линией поставлена потребность крити-ческих операций, которая должна быть обязательно удовлетворена, если нужно выполнить программу в минимально возможный срок. (Отметим, что для опера-ции (0, 2) и (1, 2) ресурсы рабочей силы не требуются.)

Построение приведённых графиков проводилось на основе информации помещённой в таблицах один и два. Из рисунков 5 б и 5 в видно, что при раннем календарном графике максимальная потребность в ресурсах составляет 18 чело-век, а при позднем 20 человек.

Применение методов Рунге-Кутта для моделирования системы типа Хищник-Жертва Курсовой

ОГЛАВЛЕНИЕ

Введение'3

1. Взаимодействие хищник-жертва. Модель Базыкина.'4

2. Система типа хищник-жертва. Модель Maк-Артура.'4

3. Система типа хищник-жертва. Модель Холлинга-Тэннера.'4

4. Система типа хищник-жертва. Модель Вольтерра.'4

5. Практическая часть'4

5.1 Описание моделей'4

5.1.1 Задание 3. Система типа хищник-жертва. Модель Базыкина.'4

5.1.2 Задание 8. Система типа хищник-жертва. Модель Maк-Артура.'4

5.1.3 Задание 12. Система типа хищник-жертва. Модель Холлинга-Тернера.'4

5.1.4 Задание 21. Система типа хищник-жертва. Модель Вольтерра.'4

5.2 Программная реализация'4

5.2.1 Задание 3'4

5.2.2. Задание 8'4

5.2.3 Задание 12.'4

5.2.4 Задание 21.'4

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

Список литературы'4

ВВЕДЕНИЕ

Итальянский математик Вито Вольтерра в 1930-е годы в серии своих работ разработал основные положения математической популяционной биологии. При разработке моделей взаимодействующих популяций учитывались процессы, происходящие в биоценозах. Основной из них — взаимоотношение между популяциями, основанное на распределении трофического ресурса между ними. Математические модели строились на основе систем обыкновенных дифференциальных уравнений с использованием принципа парных взаимоотношений для описания скоростей изменения численности популяций. Значительная часть исследований Вольтерра относится не столько к моделированию динамики реальных биологических популяций, сколько непосредственно к решению математических задач. Тем не менее, роль Вольтерра как одного из основных разработчиков математических моделей взаимодействующих популяций, огромна. Последователями Вольтерра были предложены математические модели хищник-жертва, учитывающие их взаимодействие на территории, влияние антропогенного давления. Принципы построения математических моделей для взаимодействующих популяций стали применять и в задачах медицины, экономики и социальных систем. Новый интерес к задачам популяционной биологии возник в 1960-е годы, поскольку в этот период стали заметны экологические последствия деятельности человека. Возникла необходимость в прогнозе возможных изменений в биосфере, вызванных техногенным воздействием на нее. Накопленные данные полевых наблюдений зачастую не вписывались в разработанные к тому времени математические модели. Одним из значимых исследований стала работа А. Н. Колмогорова. В отличие от используемого Вольтерра принципа парных взаимоотношений в ней был предложен иной подход для описания межпопуляционных взаимодействий — предлагалось в математические модели вводить трофические функции общего вида, качественно отражающие характер взаимодействия, как между популяциями, так и внутри популяций.

После опубликования работы А. Н. Колмогорова происходило дальнейшее развитие моделей Вольтерра. Одним из последователей теории Вольтерра был А. Д. Базыкин.

1. Взаимодействие хищник-жертва. Модель Базыкина.

Базыкин А. Д. предложил свой вариант модели «хищник-жертва», близкий по форме к модели А. Н. Колмогорова, и исследовал решения для различных трофических функций '1'. Обобщенная модель А. Д. Базыкина представлена системой дифференциальных уравнений:

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

Функция обращается в ноль при , а функция — при . Предполагается, что при , или функции и не отрицательные и не убывающие функции своих аргументов, обращающиеся в ноль при или . С учетом этого в '6' предлагается считать, что , а (a – постоянная). При этом и . С учетом этих условий в окрестности точки , функции и раскладываются в ряд по степеням своих аргументов:

В модели Вольтерра:

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

Скорость гибели хищника в отсутствие жертвы должна быть пропорциональной его численности. Поэтому квадратичным слагаемым в разложении функции C(y) в ряд по малому параметру можно в первом приближении пренебречь и считать, что .

Таким образом, в основная модель хищник-жертва, обобщающая модель А. Н. Колмогорова, представлена системой уравнений:

Эта система уравнений заменой переменных приводится к виду:

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

Стационарные точки. Система уравнений (1.5), поскольку и , имеет стационарные точки:

, и ,

В первой стационарной точке собственными значениями матрицы Якоби правой части уравнений (1.5):

будут и . Поэтому, поскольку эта стационарная точка будет неустойчивой. Во второй стационарной точке собственными значениями матрицы Якоби будут и . Поскольку функция положительно определенная, а не убывающая функция, то будет отрицательной величиной при выполнении неравенств и вторая стационарная точка в этом случае будет устойчивой.

Остальные стационарные точки находятся как решение системы алгебраических уравнений:

Значение в стационарной точке определяется из второго уравнения в (1.8) по найденному значению из уравнения:

Поскольку и , то это уравнение имеет как минимум два корня:

В точке

а в точке

Таким образом, в окрестности точки функция является положительной и возрастающей. А в окрестности точки функция будет принимать отрицательные значения и будет возрастающей при выполнении неравенства:

или:

Это неравенство при условии, что , является условием существования стационарной точки, в которой .

При малых значениях параметра a или больших значения параметра γ неравенство (1.13) может не выполняться. Тогда система уравнений (1.8) будет иметь только стационарные точки, при этом будет устойчивой только стационарная точка , . Это означает, что в рамках модели (1.5) при малых скоростях потребления (или переработки) хищником жертвы (малые значения a) или при высокой смертности хищника по сравнению с рождаемостью жертвы (большие значения ), хищник может погибнуть. При этом, поскольку параметр обратно пропорционален емкости среды, то выживаемость хищника с ростом емкости среды жертвы обеспечивается и ростом скоростью переработка хищником жертвы.

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

и, соответственно:

(1.15)

b – константа.

Для случая функции:

«не тривиальные» корни уравнения (1.8) удовлетворяют кубическому уравнению:

которое на промежутке (0,1) может, как не иметь вещественных корней, так и иметь один или три вещественных корня. Для случая a=20, b=20, γ=0,2 на рис. 1.1 отражена зависимость корней уравнения (1.17) от параметра a.

Собственные значения матрицы Якоби правой части уравнений в устойчивых стационарных точках могут быть как отрицательными, так и комплексно сопряженными с отрицательной вещественной частью. В окрестности стационарных точек с малыми значениями и при больших значениях параметра a (при больших скоростях уничтожения жертвы хищником) могут возникать колебания. На рис. 1.2 для случая a=80 отражена зависимость функций и от времени ( ). При этом значении параметра a система уравнений (1.8) имеет три стационарные точки, в которых (рис. 1.1). Рис. 1.2 соответствуют случаю, когда система из окрестности неустойчивой стационарной точки переходит в окрестность устойчивой стационарной точки (рис. 1.1 — переход из точки A в точку B). Численное интегрирование осуществлялась в среде математического пакета Matlab с применением встроенных функций ode***. Результаты (рис. 1.2) получены с относительной и абсолютной точностью равными 0,000001. При точностях равными 0,00001 (в Matlab это точность по «умолчанию») численные методы, реализованные в ode***, строят периодические колебания (рис. 1.3).

Рис. 1.1

Рис. 1.2

Рис. 1.3

Для случая функции:

«не тривиальные» корни уравнения (1.8) удовлетворяют уравнению:

которое на промежутке (0,1) может, как не иметь вещественных корней, так и иметь до пяти вещественных корней. Для случая a=20, b=20, γ=0,2 на рис. 1.4 отражена зависимость корней уравнения, лежащих на промежутке (0,1), от параметра a. Собственные значения матрицы Якоби правой части уравнений в устойчивых стационарных точках могут быть как отрицательными, так и комплексно сопряженными с отрицательной вещественной частью.

К какой стационарной точке будет стремиться решение уравнений (1.5) зависит от выбора начальных данных. Так, например, при значениях a=210 (a=20, b=20, γ=0,2) и начальных данных и решение стремится в стационарную точку , – на рис. 1.4 отражена зависимость функций и от времени.

Трофические функции у взаимодействующих популяций могут изменяться под влиянием внешних факторов. Это может быть изменение внутреннего метаболизма особей, изменение пространственного распределения популяций, изменение свойств среды обитания под влиянием антропогенного давления. Антропогенная нагрузка может привести к изменению качества трофических ресурсов или к их уничтожению, к физическому поражению особей. В иных моделях не учитывалось, что трофические функции могут изменяться во времени. Учет этого фактора в моделях взаимодействующих популяций позволяет объяснить более широкий спектр явлений, наблюдаемых в природе. Временной фактор можно учесть, считая что a, γ, b1 и b2 в уравнениях является не параметрами, а функциями времени. Для случая значение a из области «два корня» со временем приближается к области «четыре корня», рис. 1.4), на рис. 1.6 отражена зависимость функций и от времени. Этот результат соответствуeт случаю, когда после длительного совместного существования хищника и жертвы на первый взгляд устойчивого (рис. 1.5 – зона «иллюзия» устойчивости) численность обеих популяций резко уменьшается. При этом хищник погибает, а жертва постепенно восстанавливает свою численность (рис. 1.5).

Pис. 1.4

Pис. 1.5

Pис. 1.6

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

Математические модели хищник-жертва, предложенные в работах А. Д. Базыкина, объясняют существование нескольких стационарных состояний в системе хищник жертва как устойчивых, так и не устойчивых. Одновременно с этим объясняют и возможное «резкое» исчезновение хищника. Основными факторами, определяющими, соотношение численностей, являются удельная скорость роста жертвы, удельная скорость гибели хищника, скорость уничтожения жертвы хищником и скорость переработки жертвы хищником.

2. Система типа хищник-жертва. Модель Maк-Артура.

Розенцвейг и Мак-Артур разработали метод анализа взаимоотношений хищника и жертвы при более реалистичных типах поведения хищника.

Рассмотрим уравнения:

где f(x) — скорость изменения х в отсутствие хищников;

ф(х, у) — интенсивность хищничества;

k — эффективность превращения жертвы в хищника;

е — смертность хищника.

Предположим, что численность хищника ограничивается только численностью жертвы, а в остальном уравнения имеют довольно общий вид. Если допустить, что скорость, с которой данная особь хищника поедает жертву, зависит только от плотности популяции жертвы и не зависит от плотности хищников, то уравнения можно несколько упростить и при этом обращаться с ними станет легче. Тогда уравнения (2.1) принимают вид:

В этом случае у=0, если кф(х)=е, т. е. если х имеет некое постоянное значение, как на рис. 2.1. При каких обстоятельствах уравнения (2.2) приводят к расходящимся колебаниям? Как показывает рис. 2.1, расходящиеся колебания возникнут при условии, что:

1) хищник является очень эффективным охотником и поэтому может выжить при плотности популяции жертвы, равной х и составляющей лишь небольшую долю емкости среды хе;

2) кривая в плоскости (х, у) имеет максимум или ≪горб≫, как его назвал Розенцвейг.

Рис. 2.1 Расходящиеся колебания в системе хищник—жертва.

Чтобы получить кривую х=0, мы должны действовать, как показано на рис. 2.2.

Рис. 2.2 Построение кривой х=0.

Сначала мы строим кривую f(x) и множество кривых уф(х) при разных значениях у как функции от х. Пересечение кривых f(x) с кривыми дает значения х, для которых х=0 при y=yi. Таким образом, мы получаем кривую х=0 в плоскоcти (x, у).

Очевидно, что кривая имеет максимум в том и только в том случае, если одна из кривых уф(х) пересекает кривую f(x) дважды. Это может произойти по одной из двух причин:

1. Кривая f(x) имеет форму, представленную на рис. 2.3.

Рис. 2.3 Кривая продуктивности вида жертвы при условии, что

плодовитость снижается как при низких, так и при высоких плотностях популяций.

Чтобы понять смысл этого, вспомннм, что уравнение дает нам скорость увеличения численности вида в отсутствие хищничества. Это уравнение можио записать в виде , где — скорость размножения жертвы. Обратите внимание, что кривая типа представленной на рис. 2.3 означает, что R убывает как при малых, так и при больших значениях х; это может быть связано с трудностью нахождения брачных партнеров при низкой плотности популяции или с какими-либо другими причинами. Таким образом, если скорость размножения жертвы при низких плотностях популяции снижается и если ее плотность поддерживается на низком уровне вследствие воздействия хищника, то взаимодействие хищник—жертва приведет к расходящимся колебаниям. Это вряд ли должно вызывать удивление.

Pис. 2.4 Возможные типы кривой для хищника.

2. На кривой может быть горб, даже если R непрерывно убывает с увеличением х, при условии что кривая для хищника относится к типу А (рис. 2.4), а не к типу Б (как это принимается в уравнениях Вольтерра) или В. Кривая типа А вполне вероятна.

Из этого анализа непрерывно воспроизводящихся систем можно сделать следующие выводы:

1) взаимодействие хищник — жертва, при котором численность хищника ограничивается наличием жертвы, может приводить к регулярным колебаниям численности;

2) если численность жертвы ограничивается количеством необходимых ресурсов, а не воздействием хищника, то это приводит к затуханию колебаний;

3) если численность хищника ограничивается не количеством жертвы, а каким-то другим фактором, то это приводит к затуханию колебаний;

4) если имеются убежища, благодаря чему некоторое постоянное число особей жертвы оказывается недосягаемым для хищников, то это приводит к затуханию колебаний;

5) амплитуда колебаний, по всей вероятности, будет возрастать, и эти колебания могут привести к вымиранию одного или обоих видов, если хищник может прокормиться при такой плотности популяции жертвы, которая гораздо ниже допускаемой емкостью среды, и если кривая интенсивности хищничества относится к типу А (рис. 2.4), а не Б или В.

3. Система типа хищник-жертва. Модель Холлинга-Тэннера.

Еще одна нелинейная модель системы «хищник-жертва» была предложена Холлингом и Тэннером.

Скорость роста популяции жертв в этой модели равна сумме трех величин:

* скорости размножения в отсутствие хищников – rx;

* влиянию межвидовой конкуренции за пищу при ограниченных ресурсах (для случая конкурирующих производителей это влияние ограниченных сырьевых ресурсов) –

* влиянию хищников, в предположении, что хищник перестает убивать, когда насыщается –

Таким образом, имеем модель Холлинга—Тэннера:

где r, s, K, D J – положительные постоянные величины.

Эта модель имеет две важные особенности. Ее нелинейность довольно сильна, что видно из вида фазового портрета, витки которого заметно отличны от эллипсов (рис. 3.1).

Рис.3.1 Фазовый портрет колебаний численности системы хищник-жертва согласно модели Холлинга-Тэннера.

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

4. Система типа хищник-жертва. Модель Вольтерра.

Пусть x(t) и y(t) — численность жертв и хищников соответственно. Предположим, что единственным лимитирующим фактором, ограничивающим размножение жертв, является давление на них со стороны хищников, а размножение хищников ограничивается количеством добытой ими пищи (количеством жертв). Тогда в отсутствие хищников численность жертв должна расти экспоненциально с относительной скоростью а, а хищники в отсутствие жертв — также экспоненциально вымирать с относительной скоростью т. Коэффициенты а и т — коэффициенты естественного прироста жертв и естественной смертности хищников соответственно.

Пусть V=V(x) — количество (или биомасса) жертв, потребляемых одним хищником за единицу времени, причем k-я часть полученной с этой биомассой энергии расходуется хищником на воспроизводство, а остальное тратится на поддержание основного обмена и охотничей активности. Тогда уравнение системы хищник — жертва можно записать в виде:

Функцию V(x) обычно называют трофической функцией хищника или функциональным откликом (functional resроnсе) хищника на плотность популяции жертвы. Именно эти функции обычно определяются в экспериментальных работах, посвященных изучению хищничества, и к настоящему времени считается установленным, что эти функции обычно принадлежат к одному из следующих трех типов (рис. 4.1). По-видимому, динамическое поведение системы в значительной степени зависит от вида трофической функции.

Рис. 4.1. Различные типы трофических функций в системе хищник — жертва: а) этот тип характерен для беспозвоночных и некоторых видов хищных рыб; б) трофическая функция с резко выраженным порогом насыщения характерна для хищников-фильтраторов (например, многих моллюсков); в) такой тип характерен для позвоночных — организмов, проявляющих достаточно сложное поведение (нанример, способных к обучению). Аналогичный вид будет иметь трофическая функция, если жертвы могут вырабатывать защитную стратегию (например, прятаться в убежище, недоступное хищникам)

При малых, значениях х, например, когда трофические отношения в системе напряжены и почти все жертвы становятся добычей хищника, который всегда голоден и насыщения не наступает (ситуация довольно обычная в природе), трофическую функцию V(x) можно считать линейной функцией численности жертв, т. е. V=βx. Кроме того, предположим, что k-const. Тогда:

Система (1.2) с точностью до обозначений совпадает с классической моделью хищник — жертва В. Вольтерра, который показал, что эта система имеет интеграл вида:

где X=x/x*, Y=y/y*, x*=m/kβ, y*=α/β. Если x0, y0 — начальные значения численностей жертв и хищников соответственно, то:

и уравнение (4.3) описывает семейство вложенных друг в друга замкнутых кривых, соответствующих фазовым траекториям периодических решений системы (4.2).

Рис. 4.2. Фазовый портрет классической вольтерровской системы хищник — жертва.

Заметим, что при увеличении С амплитуды колебаний x и у возрастают (рис. 4.2). При минимальном значении эти кривые стягиваются в точку с координатами (x*, y*). Легко видеть, что являются решениями системы (1.2) при т. е. ее нетривиальным равновесием. Для случая малых колебаний возле этого состояния уравнения модели можно записать в виде:

Корни соответствующего характеристического уравнения есть , т. е. точка (x,y) — центр. Период малых колебаний причем колебания численности одного вида сдвинуты по фазе относительно колебаний другого на π/2.

В системе (4.2) имеется еще одно положение равновесия — начало координат. Нетрудно видеть, что эта точка — седло. Оси кооординат являются сепаратрисами, причем ось Oy входит в седло, а ось Ox – выходит из него (рис. 4.2).

Несмотря на то, что модель Вольтерра смогла объяснить многие реально наблюдавшиеся явления, у нее есть большой недостаток — негрубость (в математическом смысле этого слова) вольтерровских циклов, так что при любых сколь угодно слабых возмущениях фазовых координат система переходит с одного цикла на другой. По-видимому, более адекватные модели должны обладать этим свойством «грубости».

С точки зрения теории устойчивости состояние равновесия системы (х*, у*) — это состояние безразличного равновесия, устойчивое по Ляпунову, но не асимптотически. Отсутствие асимптотической устойчивости равновесия указывает на то, что в вольтерровской системе отсутствуют механизмы, стремящиеся сохранить ее нетривиальное равновесное состояние.

Высказанное выше утверждение об устойчивости равновесия достаточно легко доказывается построением соответствующей функции Ляпунова '2'.

5. Практическая часть

5.1 Описание моделей

5.1.1 Задание 3. Система типа хищник-жертва. Модель Базыкина.

1. Постановка задачи '2,3'. Рассмотрим динамику популяции двух видов, взаимодействующих между собой по типу хищник-жертва, при наличии внутривидовой конкуренции жертв за ограниченные ресурсы и учете фактора нелинейности размножения жертв при малых плотностях популяции. Обозначим через x=x(t) и y=y(t) плотности популяций жертв и хищников в момент времени t. Уравнения имеют следующий вид '3':

где a, b, c, d, N, K — неотрицательные числа.

Структура уравнений следующая:

* Величина ax2/(N +x)(K−x)/K определяет скорость размножения жертв в отсутствии хищников. При малых x скорость определяется величиной ax/N и является малой (гиперболический закон). При больших плотностях (до величины K) популяция растет, при x>K – уменьшается в размерах (скорость отрицательна). Таким образом, это слагаемое описывает ограниченность ресурсов: окружающая среда может обеспечивать существование только популяции плотности меньшей K.

* Слагаемое bxy описывает влияние хищников на популяцию жертв. Функция bx характеризует количество жертв, убиваемых одним хищником в единицу времени (реакция хищника на плотность популяции жертвы). Как видим, в данной модели хищники чрезвычайно кровожадны.

* Второе уравнение определяет изменение популяции хищников. Постоянная определяется естественной нормой смертности хищников. Второе слагаемое (функция dx) характеризует прирост хищников в зависимости от плотности жертв (в данной модели хищники еще к тому же и чрезвычайно плодовиты).

2. Безразмерная форма уравнений. Вводя безразмерные величины:

преобразуем уравнения (5.1) к виду:

и дополним их начальными условиями:

5.1.2 Задание 8. Система типа хищник-жертва. Модель Maк-Артура.

1. Постановка задачи '2,3'. Рассмотрим динамику популяции двух видов, взаимодействующих между собой по типу хищник-жертва, при наличии внутривидовой конкуренции жертв за ограниченные ресурсы и при учете фактора насыщения хищника. Обозначим через x=x(t) и y=y(t) плотности популяций жертв и хищников в момент времени t. Уравнения имеют следующий вид:

где a, b, c, d, A, K — неотрицательные числа.

Структура уравнений следующая:

* Слагаемое a(1−x/K) x определяет скорость размножения жертв в отсутствии хищников. При малых x скорость определяется величиной a и, т.о., a характеризует норму рождаемости при малой плотности популяции. При большей плотности (до величины K) популяция растет, при x>K — уменьшается в размерах (скорость отрицательна). Т.о. это слагаемое описывает ограниченность ресурсов: окружающая среда может обеспечивать существование только популяции плотности меньшей K.

* Слагаемое bxy/(1+Ax) описывает влияние хищников на популяцию жертв. Функция bx/(1+Ax) характеризует количество жертв, убиваемых одним хищником в единицу времени (реакция хищника на плотность популяции жертвы). Здесь учтено, что при большой плотности жертв хищник убивает лишь b/A жертв единицу времени; т.е., перестает убивать, когда насыщается.

* Второе уравнение определяет изменение популяции хищников. Постоянная определяется естественной нормой смертности хищников. Второе слагаемое характеризует прирост хищников в зависимости от плотности жертв (функция dx/(1+Ax)). При большой плотности жертв скорость прироста хищников определяется величиной d/A и, т.о., d/A характеризует норму рождаемости при благоприятных для хищников условиях.

2. Безразмерная форма уравнений. Вводя безразмерные величины:

преобразуем уравнения (5.4) к виду '4':

и дополним их начальными условиями:

5.1.3 Задание 12. Система типа хищник-жертва. Модель Холлинга-Тернера.

1. Постановка задачи '3'. Рассмотрим динамику популяции двух видов, взаимодействующих между собой по типу хищник-жертва, при наличии внутривидовой конкуренции жертв за ограниченные ресурсы и учете фактора насыщения хищника. Обозначим через x=x(t) и y=y(t) плотности популяций жертв и хищников в момент времени t. Уравнения имеют следующий вид:

где r, s, ω, K, D, J — неотрицательные числа.

Структура уравнений следующая:

* Слагаемое r(1−x/K) x определяет скорость размножения жертв в отсутствии хищников. При малых x скорость определяется величиной r и, таким образом, r характеризует норму рождаемости при малой плотности популяции. При большей плотности (до величины K) популяция растет, при x>K — уменьшается в размерах (скорость отрицательна). Таким образом, это слагаемое описывает ограниченность ресурсов: окружающая среда может обеспечивать существование только популяции плотности меньшей K.

* Слагаемое ωxy/(D+x) описывает влияние хищников на популяцию жертв. Функция ωx/(D+x) характеризует количество жертв, убиваемых одним хищником в единицу времени (реакция хищника на плотность популяции жертвы). Здесь учтено, что при большой плотности жертв хищник убивает лишь ω жертв единицу времени; т.е., перестает убивать, когда насыщается.

* Второе уравнение определяет изменение популяции хищников. При большой плотности жертв и достаточно малой плотности хищников (величина Jy/x — мала) скорость прироста хищников определяется величиной s и, т.о., s характеризует норму рождаемости при благоприятных для хищников условиях. Величина J определяется количеством жертв, необходимых для поддержания жизни одного хищника. Поскольку популяция из x жертв может поддерживать не более чем x/J хищников, то при y>x/J численность хищников должна уменьшится. Модель, т.о., ограничивает численность хищников критической величиной x/J.

2. Безразмерная форма уравнений. Вводя безразмерные величины:

преобразуем уравнения (5.7) к виду:

и дополним их начальными условиями:

5.1.4 Задание 21. Система типа хищник-жертва. Модель Вольтерра.

1. Постановка задачи '1'. Рассмотрим динамику популяции двух видов, взаимодействующих между собой по типу хищникжертва. Обозначим через x=x(t) и y=y(t) плотности популяций жертв и хищников в момент времени t. Предположим, что: •'

жертва может найти достаточно пищи для пропитания; •'

при каждой встрече с хищником последний убивает жертву; •'

норма рождаемости жертв xb, нормы естественной смертности жертв xd и хищников c являются постоянными, a=xb−xd>0. •'

число случаев, когда хищник убивает жертву, зависит от вероятности их встречи и, следовательно, пропорционально произведению xy; в результате встреч с жертвами число хищников увеличивается.

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

известной как уравнения Лотки-Вольтерра (1925 г.).

2. Безразмерная форма уравнений. Вводя безразмерные величины;

преобразуем уравнения (5.10) к виду:

Дополним их начальными условиями:

Из вида уравнений (5.11) следует, что динамика популяции характеризуется только одним cущественным положительным параметром σ (равным отношению прироста жертв к естественной смертности хищников).

5.2 Программная реализация

5.2.1 Задание 3

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

на произвольном отрезке 'a, b', используя метод Рунге-Кутта 4-го порядка точности с постоянным шагом h '5,6,7,8':

2. Тестируем программу на примере системы уравнений:

на отрезке '0,4' с точным решением:

3. Для тестовой задачи построим графики зависимости максимальной погрешности решения e и e/h4 от выбранного шага h.

' Рис. 5.1

4. Найдем стационарные решения (состояния равновесия) системы (5.2). Как они зависят от параметров задачи?

При стационарным решении все производные должны быть равны нулю:

Отсюда находим:

5. Для ряда значений параметра m из интервала '0,1; 0,35' решим систему уравнений (5.2), (5.3) при помощи разработанной программы.

Расcчитаем динамику популяции при следующих исходных данных:

n=0,1; γ=1, X0=0,3; Y0=0,3.

При каких значениях параметра m в системе появляются и исчезают автоколебания? Приведите графики наиболее характерных решений в координатах (X, Y), (X(t),t) и (Y (t),t) и дайте их интерпретацию.

5.2.2. Задание 8

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

на произвольном отрезке 'a, b', используя метод Рунге-Кутта 3-го порядка точности с постоянным шагом h:

2. Тестируем программу на примере системы уравнений:

на отрезке '0,5' с точным решением:

В результате решения задачи методом Рунге-Кутта получаем результаты:

Точные решения равны:

3. Для тестовой задачи построим графики зависимости максимальной погрешности решения e и e/h3 от выбранного шага h. Какие выводы можно сделать из полученных графиков?

' 'Рис. 5.2

4. Найдем стационарные решения (состояния равновесия) системы (5.5). Как они зависят от параметров задачи?

При стационарным решении все производные должны быть равны нулю:

Отсюда находим:

5. Решим систему уравнений (5.5), (5.6) при помощи разработанной программы. Для значений параметра α из интервала '0,1; 0,9' расcчитаем динамику популяции при следующих исходных данных:

ε=0,1; γ=1; X0=3; Y0=1.

Определим те значения параметра α при которых в системе появляются и исчезают устойчивые автоколебания (предельный цикл). Приведем графики наиболее характерных решений в координатах (X,Y ), (X(t),t) и (Y (t),t).

5.2.3 Задание 12.

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

на произвольном отрезке 'a, b', используя метод Рунге-Кутта 3-го порядка точности с постоянным шагом h:

2. Тестируем программу на примере системы уравнений:

на отрезке '0,5' с точным решением:

'

Результат решения задачи методом Рунге-Кутта:

Точные решения равны:

3. Для тестовой задачи построим графики зависимости максимальной погрешности решения e и e/h3 от выбранного шага h. Какие выводы можно сделать из полученных графиков?

4. Найдем стационарные решения (состояния равновесия) системы (5.8). Как они зависят от параметров задачи?

При стационарным решении все производные должны быть равны нулю:

Отсюда находим:

Или если перейти к переменной Y, получим уравнение четивертой степени:

Решим его методом Феррари.

5. Решим систему уравнений (5.8), (5.9) при помощи разработанной программы. Для двух наборов начальных условий (5.9) в окрестности состояния равновесия и значений параметра α из интервала '1,30' расcчитаем динамику популяции при следующих исходных данных:

β=0,13; γ=5; X0=0,5; Y0=0,5.

Определим те значения параметра α при которых в системе появляются и исчезают устойчивые автоколебания (предельный цикл). Приведем графики наиболее характерных решений в координатах (X, Y), (X(t),t) и (Y (t),t).

5.2.4 Задание 21.

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

на произвольном отрезке 'a, b', используя метод Рунге-Кутта 4-го порядка точности с постоянным шагом h:

2. Тестируем программу на примере системы уравнений:

на отрезке '0,3' с точным решением:

Решение задачи методом Рунге-Кутта дает ответ:

Точное решения:

3. Для тестовой задачи построим графики зависимости максимальной погрешности решения e и e/h4 от выбранного шага h. Какие выводы можно сделать из полученных графиков? ' Рис. 5.3

4. Найдем стационарные решения (состояния равновесия) системы (5.11). Как они зависят от параметра σ?

При стационарным решении все производные должны быть равны нулю:

Отсюда находим:

5. Решим систему уравнений (5.11), (5.12) при помощи разработанной программы. Для двух наборов начальных условий (5.12) в окрестности состояния равновесия и нескольких значений параметра σ расcчитаем динамику популяции. Какие типы решений наблюдаются? Приведем графики наиболее характерных решений в координатах (X,Y), (X(t),t) и (Y (t),t).

ЗАКЛЮЧЕНИЕ

В данной работе рассмотрены модели взаимодействия сообществ типа «хищник-жертва». Приведены примеры численного исследования данного процесса с помощью метода сравнения Е. В. Воскресенского. Проведено моделирование в программном комплексе PTC MathCAD.

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

СПИСОК ЛИТЕРАТУРЫ

1. Колпак Е. П., Горыня Е. В., Селицкая Е. А. О моделях А. Д. Базыкина «хищник — жертва» // Молодой ученый. — 2016. — №2. — С. 12-20.

2. Эрроусмит Д., Плейс К. Обыкновенные дифференциальные уравнения. Качественная теория с приложениями. - М.: Мир, 1986.

3. Базыкин А.Д. Математическая биофизика взаимодействующих популяций. - М.: Наука, 1985.

4. Хайрер Э., Нерсетт С., Ваннер Г. Решение обыкновенных дифференциальных уравнений. - М.: Мир, 1990.

5. Самарский А.А., Гулин А.В. Численные методы. - М.: Наука, 1989.

6. Бахвалов Н.С., Жидков Н.П., Кобельков Г.М. Численные методы. - М.: Наука, 1987.

7. Дж.Ортега, У.Пул. Введение в численные методы решения дифференциальных уравнений. - М.: Наука, 1986.

8. Самарский А.А., Гулин А.В. Численные методы. - М.: Наука, 1989.


Способ заказа и контакты