Линейный анализ устойчивости

Найдем стационарные состояния системы (1). Для этого приравняем к нулю правые части уравнений системы:

Получили алгебраическую систему уравнений. Решая ее, получаем единственное решение .

Это стационарное состояние системы, при котором концентрации веществ X и Y не будут зависеть от времени. Это состояние может установиться в системе, если оно будет устойчивым.

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

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

Теперь стационарное состояние (точка покоя) имеет координаты (0, 0), а система после замены переменных выглядит следующим образом:

(2)

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

Очевидно, . Величины и - очень маленькие, ими можем пренебречь. В итоге получаем систему первого приближения по отношению к системе (2):

(3)

Построим характеристическую матрицу для этой системы:

,

И выпишем характеристическое уравнение:

(4)

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

(5).

То есть стационарное состояние будет устойчивым при выполнении условия (5).

Теперь необходимо определить типы точки покоя в зависимости от параметров a и b.

Выпишем матрицу А, компонентами которой являются коэффициенты при неизвестных функциях в правых частях уравнений системы (3):

.

Заметим, что , a . Тогда уравнение (4) будет выглядеть так:

(4*)

Вычислим дискриминант и корни уравнения (4*):

Если , тогда собственные значения матрицы А – кратные вещественные. Сама матрица А – недиагональная. Получим тип точки покоя – вырожденный узел. Найдем, при каких a и b дискриминант равен нулю:

При выполнении условия (5) и точка покоя – устойчивый вырожденный узел. Если , то и точка покоя - неустойчивый вырожденный узел.

Рассмотрим случай .

Отсюда видно, что корни уравнения (4*) одного знака. Это значит, что точка покоя – узел.

При выполнении условия (5) , точка покоя – устойчивый узел. При и при - неустойчивый узел.

Если , то корни уравнения (4*) – комплексные:

Если вещественная часть корня равна нулю, то есть , то тип точки покоя – центр. Если не равна нулю, то – фокус, при - устойчивый, при - неустойчивый.

На рис.1 приведена графическая иллюстрация – диаграмма устойчивости.

Исследование модели с помощью численных методов

Зафиксируем значение параметра а=1. На рис.1 изображен фазовый портрет при b=0.5. Согласно условию (5) при таких параметрах точка покоя устойчива. Все траектории, вышедшие из разных точек, стремятся к одному и тому же аттрактору – (1, 0.5). Такой аттрактор, как видно из диаграммы устойчивости, - устойчивый фокус. Колебаний концентраций веществ X и Y во времени не происходит (рис.2). Устойчивое стационарное состояние не устанавливает колебательный режим в системе.

Рисунок 2. Фазовый портрет брюсселятора при b=0.5

Рисунок 3. Изменение концентраций X и Y во времени при b=0.5

Эволюцию фазового портрета можно наблюдать, проводя расчеты с различным параметром b. При его увеличении узел будет сначала постепенно смещаться в точку с координатами (1, b), пока не достигнет бифуркационного значения b=2 (рис.3). в этой точке происходит качественная перестройка портрета, выражающаяся в рождении предельного цикла.

Рисунок 4. Фазовый портрет брюсселятора при b=2

При дальнейшем увеличении b происходит лишь количественное изменение параметров этого цикла.

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

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

В качестве точки покоя возьмем точку на границе области устойчивости - (1, 2.000001). Тип этой точки покоя – неустойчивый фокус. Зададим начальные условия , шаг изменения t – . На рис.4 изображена траектория, соответствующая этим условиям.

Рисунок 5. Траектория при y(0) = 2.000002

Видно, что траектория вторично пересекает прямую незначительно выше начальной точки.

Теперь изменим начальное условие: (рис.5).

Рисунок 6. Траектория при y(0) = 2.0000017

Видно, что траектория вторично пересекает прямую ниже начальной точки. Согласно теореме Хопфа, между траекторией, изображенной на рис.4, и траекторией, изображенной на рис.5, находится предельный цикл (рис.6).

Рисунок 7. Схематическое изображение траекторий, близких к предельному циклу

В следующих таблицах приведены некоторые, вычисленные программой MATLAB r2007b, значения функций X и Y, по которым строились фазовые портреты на рисунках 4 и 5.

Таблица 1. Значения Х и Y при y(0) = 2.000002

Таблица 2. Значения X и Y при y(0) = 2.0000017

Жирным шрифтом выделены те значения переменной t, при которых концентрации X и Y совершают одно колебание. Теперь можно определить приблизительный период колебаний: