Модель развития лесонасаждений

МОДЕЛЬ РАЗВИТИЯ ЛЕСОНАСАЖДЕНИЙ Гипотезы модели 1) Лесной массив развивается на замкнутой однородной территории, имеющей одинаковые климатические и географические условия. 2) Деревья, принадлежащие к разным возрастным группам, равномерно распределены по занимаемой территории. 3) Рассматривается однородный лес, либо смешанный лес без явного преобладания какого-либо вида, что обеспечивает отсутствие внутривидовой конкуренции. 4) Количество деревьев в каждой возрастной группе зависит только от времени (точеная модель, не учитывает пространственное распространение деревьев по занимаемой территории). 5) Появление молодых деревьев обусловлено количеством деревьев более старшего возраста, способных плодоносить. 6) По мере взросления деревья из одной возрастной группы переходят в более старшую возрастную группу. 7) Смертность деревьев младшего возраста зависит не только от своего текущего состояния, но и от количества деревьев старшего возраста. Обозначения модели Обозначим за xi (t ) - количество деревьев в i -ой возрастной группе в момент времени t , где i  1..n , n – количество классов разбиения. Будем полагать, что  ( ) - функция, характеризующая скорость появления молодых деревьев (деревьев первой возрастной группы). Показатель  отражает способность плодоношения деревьев, значение которого определяется численностями деревьев возрастных групп, способных давать семена. Пусть l , l  1,...,n номера возрастных групп деревьев способных плодоносить. Величина показателя  определяется суммой:   kl xl  ...  k n xn , где постоянные k j , j  l , n отражают степень влияния деревьев j -ой группы на появление молодых деревьев. Естественно потребовать выполнения условий: d  0 ,  (0)  0 . d Гибель деревьев Будем рассматривать две причины гибели деревьев, а именно, естественную смертность (смертность, связанную с влиянием природных факторов) и смертность, обусловленную подавлением деревьев старших возрастных групп. Пусть  i ( i ) - функция, характеризующая интенсивность гибели деревьев i -го возрастного класса обозначим за  i  rii1 xi 1  ...rni1 xn1  rni x n где r ji  0 , j  i - постоянные коэффициенты, которые отражают влияние деревьев j -го возрастного класса на деревья i -го класса. Полагаем, что для всех i  1, n  i (0)  0 - интенсивность естественной гибели. Интенсивность естественной гибели деревьев самой старшей группы обозначим h  const  0 . Пусть f i  const  0 - интенсивность перехода деревьев i -ой группы в (i  1) -ую, i  1, n  1. Модель динамики численности деревьев n-возрастных групп  dx1  dt   ( )   1 ( 1 ) x1  f1 x1  . . .   dx  i  f i 1 xi 1   i ( i ) xi  f i xi  dt  . . .   dxn 1  f n  2 xn  2   n 1 ( n 1 ) xn 1  f n 1 xn 1  dt   dxn  dt  f n 1 xn 1  hxn В качестве функции роста  ( ) может быть выбрана любая функция, в том числе и линейная функция:  ( )     , где   const  0 множитель, характеризующий вклад деревьев, способных плодоносить в процесс появления новых молодых деревьев Функция гибели Функция гибели  i ( i ) , i  1, n  1, так же может быть задана произвольным образом. Однако исследования в этой области позволили выделить следующие способы задания функции гибели  i ( i ) : 1.  i ( i ) , i  1, n  1 - монотонная неубывающая функция. Этот выбор объясняется фактом, что скорость гибели молодых деревьев растёт с увеличением числа деревьев более стершего возраста. Например, она может быть записана в одном из следующих видов: 1.1.  i ( i )  bi ; 1.2.  i ( i )  ai i  bi ; 1.3.  i ( i )  Aai i  bi ; где ai , bi  const  0 , A  const  1 - положительные параметры, bi коэффициент естественной гибели деревьев. Функция гибели 2.  i ( i ) , i  1, n  1 - немонотонная функция, имеющая один минимум, при ненулевом значении величины  i . В этом случае функция гибели  i ( i ) может быть записана в виде  i ( i )  ai ( i  ci ) n  bi , где n - чётное целое число; ai , bi , ci  const  0 - положительные параметры. Модель динамики численности деревьев n-возрастных групп, n=3 В частном случае, когда n  3 , будем исходить из предположений, что деревья старшего возраста затрудняют развитие только подроста (деревьев 1-ой возрастной группы) и не препятствуют развитию деревьев 2-ой возрастной группы, т.е.  1 (1 )   ( x3 ) . Обозначив за x  x1 , y  x2 , z  x3 , получим трёхвозрастную модель, описываемую следующей системой дифференциальных уравнений:  x   ( y, z )   ( z ) x  fx   y  fx  (q  d ) y  z  qy  hz  (2.3.2) где,  ( y, z ) - функция, характеризующая скорость рождения деревьев младшего возраста; функция  ( z) x - функция гибели подроста, то есть интенсивность его гибели под воздействием старшей возрастной группы и в результате естественной гибели подроста; d, h коэффициенты гибели деревьев второй и третьей возрастных групп. Анализ динамической системы, описывающей развитие лесных насаждений Проведём анализ динамической системы  x   ( y, z )   ( z ) x  fx   y  fx  (q  d ) y  z  qy  hz  (2.3.2) на устойчивость в зависимости от выбора функции  (z ) . Линейная система. динамической системы: Исследуем  x  (  f ) x   z   y  f x  (q  d ) y  z  qy  h z  поведение линейной (2.6.5) Анализ динамической системы, описывающей развитие лесных насаждений Эта система может быть записана в виде: X  AX , где x   X   y , z        (  f )   A f  (q  d ) 0  .  q  h   Практическое задание 1. Найти положение равновесия системы (2.6.5) Практическое задание 2. Самостоятельно получить характеристическое уравнение (найти определитель матрицы 3х3) характеристическое уравнение: 3  a1 2  a 2   a3  0 , (2.6.6) где a1  q  d  h    f ; a2  h(q  d )  (   f )(q  d  h); a3   fq  h(   f )(q  d ) . Практическое задание 3. Решить характеристическое уравнение и найти его корни Анализ динамической системы, описывающей развитие лесных насаждений Можно сделать вывод, что положение равновесия системы (2.6.5) будет асимптотически устойчиво в том и только том случае, когда q  d  h    f  0, (q  d  h    f )(h(q  d )  (  f )(q  d  h))  fq  h(  f )(q  d )  0, . h(  f )(q  d )  fq  0 Проанализируем возможные режимы системы (2.6.5), режим которой зависит от корней характеристического уравнения 3  a12  a2   a3  0 . В поле комплексных чисел это уравнение имеет три корня, для которых возможны следующие случаи: 1. если i  R , i  1,2,3 , тогда система имеет монотонный режим; 2. если 1  R , 2 , 3  bi , 3. если 1  R , 2 , 3  a  bi , b  0, b R - колебательный режим; a, b  R , a, b  0 возбуждающийся колебательный режим. - затухающий или Анализ динамической системы, описывающей развитие лесных насаждений Исследуем ситуации возникновения каждого из режимов. Для этого рассмотрим кубическую функцию f ( )  3  a12  a2   a3 . Она будет иметь только действительные корни в том и только том случае, когда уравнение f ' ( )  0 будет иметь два различных действительных корня, (дискриминант будет положительным) и для корней этого уравнения будет справедливо следующее соотношение f ( ' ) f ( ' ' )  0 . Рассмотрим ситуацию, когда корнями характеристического уравнения являются комплексные числа. Пусть 1  c , 2 , 3  a  bi , a, b  R , b  0 - корни характеристического уравнения. В этом случае уравнение будет выглядеть следующим образом: f ( )  (  c)(  a  bi)(  a  bi)  0 f ( )  3  (2a  c)2  (2ac  a 2  b 2 )  c(a 2  b 2 )  0 Анализ динамической системы, описывающей развитие лесных насаждений Найдём зависимость между коэффициентами полученного уравнения. Если 2 , 3  bi , то есть a  0 , тогда f ( )  3  c2  b 2   cb 2  0 . Выпишем зависимости коэффициентов полученного уравнения:  если c  0, то  cb 2  b2  0 ; c  если c  0 , то уравнение принимает вид f ( )  3  b 2   0 . Нетрудно проверить, что выполнение этих соотношений невозможно, в случае, когда 2 , 3  a  bi , a  0 . Делаем вывод, что уравнение (2.6.6) будет иметь два чисто комплексных корня, в том и только том случае, когда для коэффициентов этого уравнения справедливы соотношения:  a3   a2  0  a1 a  a  0, a  0 3 2  1 Анализ динамической системы, описывающей развитие лесных насаждений Анализ динамической системы, описывающей развитие лесных насаждений Итак, для динамической системы (2.6.5) характерны режимы: 1.) Монотонный режим, который имеет две возможные реализации  функции x(t ), y(t ), z(t ) являются убывающими и стремятся к положению равновесия, при этом положение равновесия устойчиво;  функции x(t ), y(t ), z(t ) являются возрастающими и удаляются от положения равновесия, при этом положение равновесия неустойчиво; 2.) Режим возбуждающихся (при a1  0) или затухающих (при a1  0) колебаний. Однако этот режим является ярко выраженным лишь на небольшом временном промежутке. При T   он вырождается в монотонный. 3.) Колебательный режим практически не реализуем. Анализ динамической системы, описывающей развитие лесных насаждений Нелинейная система. Проанализируем структуру динамической системы (2.3.2) в предположении, что функция гибели подроста является функцией четвёртой степени, т.е  ( x, z)   ( z) x ,  ( z )  az 4  b , где коэффициент b характеризует смертность подроста в отсутствии деревьев старшего возраста, а коэффициент a - степень влияния деревьев верхнего яруса на смертность подроста; функция рождения является линейной:  ( y, z)   z . Модель, описывающая процесс развития лесонасаждений примет вид: dx   z  (az 4  b) x  fx dt dy  fx  (q  d ) y dt dz  qy  hz dt (2.6.7) Анализ динамической системы, описывающей развитие лесных насаждений Система, описываемая полученной динамической системой, имеет одну точку равновесия: x*  qd * y f fh3 bh 4 fh4 y 4 3   aq (q  d ) aq 4 aq 4 * z*  (2.6.8) q * y h Из этих соотношений следует, что положение равновесия существует, если выполнены следующие условия: fh 3 bh 4 fh 4   0 aq 3 ( q  d ) aq 4 aq 4  , a, b, f , q, d , h  0 Практическое задание 4. Найти положение равновесия системы (2.6.7) Анализ динамической системы, описывающей развитие лесных насаждений Характер особой точки определяется характером поведения фазовых координат в её малой окрестности. Ниже представлены характерные режимы нелинейной динамической системы. Рассмотрим нулевое положение равновесия при следующих параметрах динамической системы:   0,26 , f  0,2, q  0,16, d  0,1, h  0,05 и  ( z)  az 4  b a  3,75 1014, b  0,44 . Положение равновесия в этом случае имеет координаты: X *  0,009; Y *  0,006; Z *  0,2. Система имеет одно нулевое положение равновесия, которое является устойчивым фокусом. С физической точки зрения данную ситуацию можно интерпретировать как вымирание системы. Анализ динамической системы, описывающей развитие лесных насаждений а) б) в) Рисунок 2.6.5. Нулевое положение равновесия Анализ динамической системы, описывающей развитие лесных насаждений Рассмотрим ещё один режим, а именно, монотонный при следующих параметрах динамической системы:   0,3, f  0,15, q  0,1, d  0,05, h  0,05 и  ( z)  az 4  b, a  5 1015, b  0,29 . Тогда координаты положения равновесия будут иметь вид: X *  668,7; Y *  668,7; Z *  1337,5. а) б) Анализ динамической системы, описывающей развитие лесных насаждений в) Рисунок 2.6.6. Монотонный режим Практическое задание 5. В MS Excel или другом программном продукте Реализовать монотонный режим системы (2.6.7) при заданных параметрах Анализ динамической системы, описывающей развитие лесных насаждений Рассмотрим колебательный режим нелинейной динамической системы, при следующих параметрах динамической системы:   0,3, f  0,31, q  0,18, d  0,08, h  0,06 и  ( z)  az 4  b, a  4 1014, b  0,29 . Координаты положения равновесия: X *  518,5; Y *  618,2; Z *  1854,5. Отметим, что колебания возникают в случае, когда начальные условия достаточно удалены от положения равновесия, в противном случае наблюдается монотонный режим. а) б) Анализ динамической системы, описывающей развитие лесных насаждений в) Рисунок 2.6.7. Колебательный режим