Брюсселятор
До сих пор в этой главе в качестве примеров расчета динамических систем мы приводили графики 1-3 траекторий на фазовой плоскости. Однако для надежного исследования фазового портрета необходимо решить систему ОДУ большое количество раз с самыми разными начальными условиями (и, возможно, с разным набором параметров модели), чтобы посмотреть, к каким аттракторам сходятся различные траектории. В Mathcad можно реализовать эту задачу, например, в форме алгоритма, приведенного в листинre 9.16 для решения системы уравнений автокаталитической химической реакции с диффузией. Эта модель, называемая моделью брюсселятора, предложена в 1968 г. Лефевром и Пригожиным. Неизвестные функции отражают динамику концентрации промежуточных продуктов некоторой реальной химической реакции. Параметр модели в равен исходной концентрации катализатора.
Листинг 9.16. Построение фазового портрета для модели брюсселятора:
Предложенный алгоритм формирует из отдельных матриц решений системы ОДУ с разными начальными условиями объединенную матрицу U. Пары начальных условий задаются в первой строке листинга в виде матрицы v размера 2х10. Последнее означает построение десяти траекторий. Для того чтобы поменять количество траекторий, измените соответствующим образом размер этой матрицы. Затем (рис. 9.24) элементы матрицы и выводятся на график в виде отдельных точек. Отсутствие соединения точек линиями является недостатком алгоритма, но это плата за возможность представить в Mathcad несложным образом сразу большое количество траекторий на фазовой плоскости.
Рис. 9.24. Фазовый портрет брюсселятора при В=0.5 (продолжение листинга 9.16)
Как видно из рис. 9.24, все траектории, вышедшие из разных точек, асимптотически стремятся к одному и тому же аттрактору (1.0.5). Из теории динамических систем нам известно, что такой аттрактор называется узлом (с узлом мы уже встречались в примере разд. 9.4.1). Конечно, в общем случае при анализе фазового портрета желательно "прощупать" большее число траекторий, задавая более широкий диапазон начальных условий. Не исключено, что в других областях фазовой плоскости траектории будут сходиться к другим аттракторам.
Эволюцию фазового портрета брюсселятора можно наблюдать, проводя расчеты с различным параметром в. При его увеличении узел будет сначала постепенно смещаться в точку с координатами (1, в), пока не достигнет бифуркационного значения B=2. В этой точке происходит качественная перестройка портрета, выражающаяся в рождении предельного цикла. При дальнейшем увеличении в происходит лишь количественное изменение параметров этого цикла. Решение, полученное при B=2.5, показано на рис. 9.25.
Примечание
Чтобы найти аттракторы динамической системы, как известно, нужно решить систему алгебраических уравнений, получающуюся из системы ОДУ заменой нулями их левых частей. Эти задачи также удобно решать средствами Mathcad (см. главу 5). В частности, исследование зависимости фазового портрета от параметров системы ОДУ и поиск бифуркаций можно проводить методами продолжения (см. разд. 5.3.3).
Рис. 9.25. Фазовый портрет брюсселятора при В=2.5
Читатели, сталкивающиеся с расчетом динамических систем, несомненно, оценят возможности Mathcad по построению фазовых портретов и исследованию бифуркаций. Возможно также, что они найдут лучшие программные решения этой задачи, чем алгоритм, предложенный в данном разделе автором.