Программа mep2 для исследования динамических систем с клеточной структурой Текст научной статьи по специальности «Математика»

ИНФОРМАЦИОННЫЕ ТЕХНОЛОГИИ

Вестн. Ом. ун-та. 2011. № 2. С. 167-172.

УДК 004.42

В.В. Коробицын, Ю.В. Фролова

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

ПРОГРАММА МЕР2 ДЛЯ ИССЛЕДОВАНИЯ ДИНАМИЧЕСКИХ СИСТЕМ С КЛЕТОЧНОЙ СТРУКТУРОЙ

Представлен обзор функциональности программы МЕР2, разработанной авторами для исследования систем обыкновенных дифференциальных уравнений с клеточной структурой. На примере описана последовательность действий для решения задачи с помощью этой программы. Возможности графического интерфейса демонстрируются на графиках и фазовых портретах системы в 20 и 30 пространствах. Даны рекомендации по использованию программы.

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

Введение

Разрабатываемый авторами проект MEP2 (Modelling Evolution Processes 2) [1] направлен на создание программного обеспечения для анализа динамических систем с разрывной правой частью [2]. Модели процессов различной природы описываются подобными динамическими системами: перенос биомассы в цепи питания организмов (растительность - травоядные - хищники) [3]; функционирование электрических схем с нелинейными элементами [4]; сложные кинематические механизмы и др. Многие из подобных моделей описываются динамическими системами с клеточной структурой [5]. Это системы обыкновенных дифференциальных уравнений, в которых правая часть непрерывна на замкнутых подпространствах области определения системы. То есть на границах этих областей правая часть может претерпевать разрыв. Для таких систем алгоритмы для численного решения находятся в стадии разработки различными коллективами [6; 7]. Представленная программа является программным решением алгоритмов [8; 9], разрабатываемых авторами в этой области.

Динамическая система с клеточной структурой

Динамической системой с клеточной структурой будем считать систему вида

Г fj(x),x е Cj, dx fz(x),xе C2,

dt ~= •••

{fk(x),xе Ck,

где x(t) - вектор-функция, определенная в некоторой области QeRn. Разбиение области Q на закрытые подмножества Ci, С2,..., Ck с пересечением на множестве меньшей размерности имеет вид клеток с об© В.В. Коробицын, Ю.В. Фролова, 2011

общей границей между клетками. Функции fi, J2,..., fk определены и непрерывны в соответствующих областях. Вопросы существования решения таких систем соприкасаются с монографией А.Ф. Филиппова [2], поскольку данная системы является частным случаем систем Филиппова. Особенностями поведения таких систем являются: негладкость решения, реализация скользящего режима вдоль границы клетки, неединственность решения.

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

Проект программы MEP2

Программа была спроектирована при помощи объектно-ориентированного подхода. Поскольку количество классов в проекте достаточно велико, то они были объединены в пакеты. Диаграмма пакетов (рис. 1) дает представление о структуре программы. Перечислим основные функции пакетов:

- mep2 описывает пользовательский интерфейс программы;

- graph содержит классы построения графиков функций;

- difsys включает классы представления динамической системы;

- mepfile реализует функции работы с файлами;

- progra обеспечивает компиляцию и выполнение текста описания моделей;

- btree содержит класс обработки бинарных деревьев;

- unistack включает классы унифицированного стека;

- unicalc описывает универсальный вычислитель;

- node сохраняет информацию о сущностях в памяти;

- param реализует основные операции с сущностями.

Программа была реализована на языке C++ в среде программирования Microsoft Visual Studio 2005. Программа построена на основе базовых API функциях Windows и предназначена для работы в операционной системе семейства Microsoft Windows. В проекте использована библиотека OpenGL для обеспечения построения графиков функций в 2D и 3D пространстве. Общий объем исходного кода программы составил 415 КБ.

Рис. 1. Диаграмма пакетов программы MEP2

3. Основные функции программы

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

• ввод и редактирование системы уравнений;

• загрузку и сохранение файла описания системы;

• нахождение решения систем ОДУ с разрывами;

• вывод графиков решений системы;

• построение проекций фазового пространства на плоскость и на трехмерное подпространство.

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

Описание системы уравнений можно сохранить в файл. Файл имеет текстовый вид и строго заданную структуру. При необходимости файл может быть отредактирован без использования программы МЕР2. При загрузке файла системы программа интерпретирует содержание файла и выдает сообщения об ошибках при их возникновении. Кроме самой сис-

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

В программе реализованы алгоритмы численного решения систем обыкновенных дифференциальных уравнений, в том числе уравнений с разрывной правой частью. В частности, реализованы численные схемы: Эйлера, Рунге - Кутты 2-го и 4-го порядков, Фельберга, Дормана-Принса [10] и авторские алгоритмы решения систем с разрывной правой частью, описанные в работах [8; 9].

Вывод графиков решений системы осуществляется в одной системе координат для всех переменных системы одновременно. Численная схема, временной

интервал и параметры сетки задаются в настройках программы.

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

а) б)

Рис. 2. Примеры фазовых портретов системы, построенных в МЕР2: а) на плоскости; б) в трехмерном пространстве

Порядок решения задач в МЕР2

Рассмотрим пример модели резонансного конвертора с обратной связью [4] в случае постоянного напряжения нагрузки

с ис = х и }ь = Х2:

' бх± = ^

< СО. С’

| Сх2 х1 + Я • х2 - и}

. С ~ ь ’

где управляющий параметр щ принимает значения щ1 = Е - Ц), щ2 = Ц - Е, и3 = -Ц) и4 = Ц в соответствующих областях Оі, } = 1, 2, 3, 4. Области определяются следующими условиями:

01 : х2 > О’Х^ + х2 < Г,

02 : х2 < 0,х1 + х§ < }2Г,

03 :х2 > 0,х2 + х2 > }2Г,

04 : х2 < О.х2! + х2 > Г,

где ir > 0 - заданная величина тока.

В данной цепи напряжение нагрузки U0 не может превосходить напряжения питания E, поэтому E>U0>0. Система обязательно имеет стационарные точки (xi ,x2) = (w/,0), а также другие стационарные точки, количество и расположение которых зависит от величин E, U0 и ir. В рассматриваемом примере будут использоваться следующие значения параметров: R = 0,2 Ом, L = 3110-6 Гн, C = 2-10-6 Ф.

Описание модели в программе MEP2 можно записать в виде текстового файла. Файл состоит из двух частей. Первая описывает систему обыкновенных дифференциальных уравнений с клеточной структурой. Вторая задает параметры программы, используемые для этой модели. В первой части ключевые слова: Variables, Parameters, Chambers, Chamber

- определяют области описания соответ-

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

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

Листинг 1. Первая часть файла описания модели резонансного конвертера:

System_of_Differential_Equations_with_Chamber_Structure Variables: 2

x1 -150 150 50 x2 -100 100 10 Parameters: 6 C 0 1 2e-006 L 0 1 3.1e-005 R 0 1 0.2 ir 0 100 50 E 0 1000 500 U0 0 1000 100 Chambers: 4 Chamber: 1

x2>0&x1A2+x2A2<irA2

dx1/dt=x2/C

dx2/dt=-(x1+R*x2-(E-U0))/L Chamber: 2

x2<0&x1A2+x2A2<irA2

dx1/dt=x2/C

dx2/dt=-(x1+R*x2-(U0-E))/L Chamber: 3

x2>0&x1A2+x2A2>irA2 dx1/dt=x2/C dx2/dt=-(x1+R*x2+U0)/L Chamber: 4

x2<0&x1A2+x2A2>irA2

dx1/dt=x2/C

dx2/dt=-(x1+R*x2-U0)/L

End_of_System

Корректно созданный файл модели открывается программой и отображается в главном окне программы (рис. 3). Это окно разделено на четыре части, которые имеют следующие названия: Variables

(Переменные), Parameters (Параметры), Chambers (Клетки), Equations (Уравнения). Эти области используются для отображения одноименных сущностей модели. Для редактирования сущностей модели можно использовать кнопки: Add (Добавить), Del (Удалить), Edit (Редактировать). Список уравнений определяется для каждой клетки отдельно.

Для построения графиков решения системы необходимо задать параметры построения графиков (Graphs Settings), здесь можно выбрать доступный численный метод, задать начальный и конечный момент времени, способ выбора шага ин-

тегрирования (автоматический или постоянный), параметры сетки графика. После этого запустить расчет модели и построение графика (Graph^Draw Graph). Тогда мы получим изображение графика функций решения системы и информацию о процессе вычисления (рис. 4). Информация о результате вычислений показывает количество принятых и отброшенных шагов интегрирования, минимальный и максимальный шаг интегрирования, оценку погрешности, количество вычислений правых частей системы, затраченное расчетное время. Одновременно с выводом на экран траектория решения сохраняется в текстовый файл sys_out.txt, который можно использовать для дальнейшего анализа результатов вычислений.

• МЕР2 [Е12р.тер]

В“®

File Graph Phase Phase3D Help

н1: [-150., 150.], х1 (0)=50. С: [0.,1.]= 2.Є-006

х2: [-100., 100.], х2(0)=10. L: [0.,1.]= 3.1 е-005

R: [0.,1.]= 0.2

ir: [0., 100.]= 50.

Е: [0., 1000.]= 500.

U0: [0., 1000.]= 100.

Chamber Structure

Chambers: Add Del Edit For chamber: |

else Equations: Edit

х2<08сх1 л2+к2л2<ігл2 dx1/dt=x2/C

х2>0&х1 л2+х2л2>ігл2 dx2/dt=-(x1+R*x2-(E-U0]]/L

х2<08сх1 л2+х2л2> ігл2

Рис. 3. Основное окно программы MEP2

2 System Dynamics Graph □є®

*10л1 11. -'ІП - <1 Given paiameters: 2. Solving method: Runge-Kutta 4 Modified Local accuracy: 1.e-008 Maximum amount of steps: 10000 Final time: 4.e-005 P a-M ■ Ih r h ati

-

I u. Q - —_ __

У. й - / neSUH StausuCS. Finished time: 4 e-005 Count of accepted steps: 99 M inimum integration step: 1.547830981980828e-009 Maximum integration step: 6.760792472453236e-006 Right-hand calculation count: 1252 Main ratio error: 4.797360027115114e-009 Count of rejected steps: 6 Solving time (sec): 0.841

о 7 _

R -

0. tr _

J. 1

4. - о _ f\

■5. 2 h A _ 4 \

Ш \

I . 0. H 1 \

Л _

0. 0.5 1. 1.5 2. 2.5 3. 3.5 4. 4.5 *10A-5

Рис. 4. Окно отображения графиков функций решения системы и информации о процессе вычисления

Фазовое пространство системы можно построить в отдельном окне программы (рис. 5). Для этого необходимо задать соответствующие параметры расчета

(Phase^Settings) и запустить расчет (Phase^Draw Phase). Время, необходимое для построения портрета, зависит от модели компьютера, сложности модели и параметров фазового портрета. Синие точки на фазовой плоскости указывают на начальные точки траектории, красная и зеленая линии - это траектории в поло-

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

□С

Рис. б. Окно отображения фазового портрета системы

Рекомендации по использованию и заключение

Представленная программа МЕР2 позволяет исследовать класс систем с разрывной правой частью (динамических систем с клеточной структурой). Такими системами могут описываться модели различной природы. Представленный инструмент можно использовать для исследования моделей такого класса.

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

ЛИТЕРАТУРА

[1] Коробицын В. В. MEP2 - программа для моделирования эволюционных процессов. - URL: http://users.univer.omsk.su/~korobits/mep2 (дата обращения: 01.06.2010).

[2] Филиппов А. Ф. Дифференциальные уравнения с разрывной правой частью. - М. : Наука, 1985.

[3] Krivan V. Optimal Intraguild Foraging and Population Stability // Theoretical Population Biology. 2000. Vol. 58. P. 79-94.

[4] Melin J, Hultgren A. A Limit Cycle of a Resonant Converter // Engell H. G. S., Zaytoon J. Analysis and Design of Hybrid Systems 2003. IFAC, Elsevier, 2003. P. 169-174.

[5] Коробицын В. В., Маренич В. Б., Фролова Ю. В. Численный алгоритм склеивания ре-

шения на поверхности разрыва для систем с клеточной структурой // Вестн. Ом. ун-та. 2010. № 4.

[6] Piiroinen P. T., Kuznetsov Yu. A. An event-driven method to simulate Filippov systems with accurate computing of sliding motions // ACM Trans. Math. Software. 2008. Vol. 34. № 13. P. 1-24.

[7] Shampine L. F, Thompson S. Event location for

ordinary differential equations // Computer and Mathematics with Application. 2000. Vol. 39.

P. 43-54.

[8] Коробицын В. В., Фролова Ю. В., Маре-

нич В. Б. Алгоритм численного решения кусочно-сшитых систем // Вычисл. технологии. 2008. Т. 13. № 2. С. 70-81.

[9] Коробицын В. В., Фролова Ю. В. Алгоритм вычисления скользящего режима для системы с гладкой границей разрыва // Вычисл. технологии. 2010. Т. 15. № 2. С. 56-72.

[10] Хайрер Э, Нёрсетт С., Ваннер Г. Решение обыкновенных дифференциальных уравнений. Нежесткие задачи. М. : Мир, 1990.

• System Dynamics Phase Portrait