Моделирование реальных процессов в пакете расширений «Dynamics» системы Maxima

ФЕДЕРАЛЬНОЕ ГОСУДАРСТВЕННОЕ БЮДЖЕТНОЕ

ОБРАЗОВАТЕЛЬНОЕ УЧРЕЖДЕНИЕ ВЫСШЕГО ОБРАЗОВАНИЯ

«Мордовский государственный

педагогический институт имени М. Е. Евсевьева»

Физико-математический факультет

Кафедра информатики и вычислительной техники

Реферат:

ИМИТАЦИОННОЕ МОДЕЛИРОВАНИЕ В MAXIMA. ВОЗМОЖНОСТИ МАТЕМАТИЧЕСКОГО ПАКЕТА MAXIMA

Направление подготовки 44.03.05 Педагогическое образование.

Профиль Информатика. Математика

Руководитель работы

канд. физ.-мат. наук, доцент___________________________ Т. В. Кормилицына

Саранск 2020

 Системы компьютерной математики

CAS были созданы в 70-ые годы и развивались в рамках проектов, связанных с искусственным интеллектом. Поэтому сфера применения их достаточно большая и разнообразная. Первыми популярными системами были Reduce, Derive, Macsyma. Некоторые из них до сих пор находятся в продаже. Свободно распространяемая версия Macsyma — Maxima. На данный момент лидерами продаж являются Maple и Mathematica. Оба этих пакета активно используются в математических, инженерных и других научных исследованиях. Существует множество коммерческих систем компьютерной алгебры: Maple, Mathematica, MathCad и другие. Свободно распространяемые программы: Axiom, Eigenmath, Maxima, Yacas и др.

Успех в современном использовании САВ лежит в интеграции всех машинных возможностей (символьный и численный интерфейс, встроенная графика, мультипликация, базы и банки данных и т. д.). Все современные коммерческие системы компьютерной математики (Mathematica, Maple, MatLab и Reduce) обладают стандартным набором возможностей:

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

        имеются основные символьные (математические) объекты: полиномы, ряды, рациональные функции, выражения общего вида, векторы, матрицы;

        системы используют целые, рациональные, вещественные, комплексные числа;

        имеется несколько дополняющих друг друга режимов работы: редактирование, диагностика, диалог, протокол работы;

        присутствует связь со средствами разработки программ: возможны подстановки, вычисления значений, генерация программ, использование стандартного математического обеспечения (библиотек);

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

Области математики, поддерживаемые в Maxima

        Операции с полиномами (манипуляция рациональными и степенными выражениями, вычисление корней и т.п.)

        Вычисления с элементарными функциями, в том числе с логарифмами, экспоненциальными функциями, тригонометрическими функциями

        Вычисления со специальными функциями, в т.ч. эллиптическими функциями и интегралами

        Вычисление пределов и производных

        Аналитическое вычисление определённых и неопределённых интегралов

        Решение интегральных уравнений

        Решение алгебраических уравнений и их систем

        Операции со степенными рядами и рядами Фурье

        Операции с матрицами и списками, большая библиотека функций для решения задач линейной алгебры

        Операции с тензорами

        Теория чисел, теория групп, абстрактная алгебра

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

Достоинства программы

Основными преимуществами программы Maxima являются:

        возможность свободного использования (Maxima относится к классу свободных программ и распространяется на основе лицензии GNU);

        возможность функционирования под управлением различных ОС (в частности Linux и Windows™ );

        размер программы (дистрибутив занимает порядка 23 мегабайт, в установленном виде со всеми расширениями потребуется около 80 мегабайт);

        широкий класс решаемых задач;

        возможность работы как в консольной версии программы, так и с использованием одного из графических интерфейсов (xMaxima, wxMaxima или как плагин (plug-in) к редактору TexMacs);

        расширение wxMaxima (входящее в комплект поставки) предоставляет пользователю удобный и понятный интерфейс, избавляет от необходимости изучать особенности ввода команд для решения типовых задач;

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

        наличие справки и инструкций по работе с программой (русскоязычной версии справки нет, но в сети Интернет присутствует большое количество статей с примерами использования Maxima);

Пример командного файла для решения жёсткой системы

 ОДУ в Maxima (данная система нелинейна, поэтому используем метод Рунге-Кутта, однако расчёты затрудняются жёсткостью системы):

load("dynamics");

load("draw");

k1:0.1; k2:100; k3:10;

eq1:-k1*ca+k3*cb*cc;

eq2:k1*ca-k3*cb*cc-k2*cb;

eq3:k2*cb;

t_range:[t,0,100,0.01];

sol: rk([eq1,eq2,eq3],[ca,cb,cc],[1,0,0],t_range)$

len:length(sol);

t:makelist(sol[k][1],k,1,len)$

ca:makelist(sol[k][2],k,1,len)$

cb:makelist(sol[k][3],k,1,len)$

cc:makelist(sol[k][4],k,1,len)$

draw2d(title="Chemical system",xlabel="ca",ylabel="cb",

     grid=true,points_joined =true,points(t,ca),

     points(t,cb),points(t,cc),terminal=eps);

Данная система достаточно трудно решается при помощи функции  . Увеличение констант до  и  делает задачу практически неразрешимой средствами пакета dynamics.

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

Рис. 6.7. Изменение концентраций при моделировании автоколебательной химической реакции (брюсселятора)

Здесть брюсселятор рассматривается как пример автоколебательной системы.

Описание модели брюсселятора в Maxima приведено в следующем командном файле:

load("dynamics");

load("draw");

B:0.5;

eq1:-(B+1)*y0+y0^2*y1+1;

eq2:B*y0-y0^2+1;

t_range:[t,0,10,0.1];

sol: rk([eq1,eq2],[y0,y1],[1,1],t_range)$

len:length(sol);

t:makelist(sol[k][1],k,1,len)$

y0:makelist(sol[k][2],k,1,len)$

y1:makelist(sol[k][3],k,1,len)$

draw2d(title="Brusselator",xlabel="t",ylabel="y0,y1",

grid=true,points_joined = true,

points(t,y0),points(t,y1),terminal=eps);

Графическая иллюстрация (автоколебательный режим в системе) — на рис.67., а фазовые портреты — на рис.6.8 и рис.6.9.

Рис. 6.8. Фазовый портрет для брюсселятора (В=0.5)

Рис. 6.9. Фазовый портрет для брюсселятора (В=2.5)

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

Модель динамики популяций

Модель взаимодействия "хищник—жертва" независимо предложили в 1925–1927 гг. Лотка и Вольтерра. Два дифференциальных уравнения моделируют временную динамику численности двух биологических популяций жертв  и хищников  . Предполагается, что жертвы размножаются с постоянной скоростью, а их численность убывает вследствие поедания хищниками. Хищники же размножаются со скоростью, пропорциональной количеству пищи, и умирают естественным образом.

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

Эта модель Вольтерра-Лотка с логистической поправкой описывается системой уравнений

с условиями заданной численности "жертв" и "хищников" в начальный момент  .

Решая эту задачу при различных значениях, получаем различные фазовые портреты (обычный колебательный процесс и постепенная гибель популяций). Результаты приведены на рис.6.10, рис.6.11 и рис.6.12.

(%i1) a:4$ b:2.5$ c:2$ d:1$ alpha=0$

     eq1:'diff(y(t),t)=(a-b*z(t))*y(t)-alpha*y(t)^2;

     eq2:'diff(z(t),t)=(-c+d*y(t))*z(t)-alpha*y(t)^2;

     atvalue(y(t),t=0,3); atvalue(z(t),t=0,1);

(%i10)     desolve([eq1,eq2],[y(t),z(t)]);

rat: replaced -2.5 by -5/2 = -2.5

rat: replaced -2.5 by -5/2 = -2.5

rat: replaced -2.5 by -5/2 = -2.5

rat: replaced 2.5 by 5/2 = 2.5

Очевидная проблема — неразрешимость данной системы в явном виде методом преобразования Лапласа, т.к. она нелинейна.

Используем численный метод Рунге-Кутта из пакета dynamics.

Результаты решения для значений  и  представлены на рис.6.11 и рис.6.12.

Рис. 6.10. Фазовый портрет для системы Лотка-Вольтерра

Рис. 6.11. Решения системы Лотка-Вольтерра в зависимости от времени (a = 0)

Рис. 6.12. Решения системы Лотка-Вольтерра в зависимости от времени (a = 0, 1)

Рассмотрим командный файл для задачи моделирования системы Лотка-Вольтерра в Maxima:

a:4; b:2.5; c:2; d:1; alpha1:0;

ode1:(a-b*x)*y-alpha1*x^2$ ode2:(-c+d*y)*x-alpha1*y^2$

alpha2:0.02;

ode3:(a-b*x)*y-alpha2*x^2$ ode4:(-c+d*y)*x-alpha2*y^2$

load("dynamics");

t1:[]$ xg1:[]$ yg1:[]$ t2:[]$ xg2:[]$ yg2:[]$

l1:rk([ode1,ode2],[y,x],[1,3],[t,0,9,0.01])$

l2:rk([ode3,ode4],[y,x],[1,3],[t,0,9,0.01])$

for i:1 thru length(l1) do(t1:append(t1,[l1[i][1]]),

     xg1:append(xg1,[l1[i][2]]),

     yg1:append(yg1,[l1[i][3]]));

for i:1 thru length(l2) do(t2:append(t2,[l2[i][1]]),

      xg2:append(xg2,[l2[i][2]]),

     yg2:append(yg2,[l2[i][3]]));

load("draw");

draw2d(terminal='eps, file_name="lotka1",

                            grid=true,xlabel = "x",

     ylabel = "y",

     title="Lotka-Volterra system, phaze portrait",

     key= "alpha=0",points_joined = true,

     point_type = none,

     line_width = 4,color = black, points(xg1,yg1),

     points_joined = true, color = black,

     point_type = none,

     line_width = 1,key="alpha=0.02", points(xg2,yg2))$

     draw2d(terminal='eps, file_name="lotka2",

                            grid=true,xlabel = "t",

     ylabel = "x,y", title="Lotka-Volterra system,

                                            alpha=0",

     key= "x(t)",points_joined = true, line_width = 1,

     color = black,point_type = none, points(t1,xg1),

     points_joined = true, line_width = 4,

                                      point_type = none,

     color = black, key= "y(t)", points(t1,yg1))$

draw2d(terminal='eps, file_name="lotka3",

                              grid=true,xlabel = "t",

     ylabel = "x,y", title="Lotka-Volterra system,

                                           alpha=0.02",

     key= "x(t)",points_joined = true,

                                   point_type = none,

     line_width = 1, color = black, points(t2,xg2),

     points_joined = true, line_width = 4,

                                    point_type = none,

     color = black, key= "y(t)", points(t2,yg2))$

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

Движение твердого тела

Рассмотрим пример построения трехмерного фазового портрета. Находим решение задачи Эйлера свободного движения твердого тела:

(%i1)    eq1:'diff(x1(t),t)=x2(t)*x3(t);

    eq2:'diff(x2(t),t)=-x1(t)*x3(t);

    eq3:'diff(x3(t),t)=-0.51*x1(t)*x2(t);

(%i4)    load("dynamics")$ l: rk([y*z, -x*z,0.51*x*y],

    [x,y,z],[1,2,3],[t,0,4,0.1])$

Фазовый портрет для данной динамической системы (трехмерная кривая) представлен на рис.6.13.

Рис. 6.13. Фазовый портрет трехмерной динамической системы

Аттрактор Лоренца

Одна из самых знаменитых динамических систем предложена в 1963 г. Лоренцем в качестве упрощенной модели конвективных турбулентных движений жидкости в нагреваемом сосуде тороидальной формы. Система состоит из трех ОДУ и имеет три параметра модели. Задаём правые части уравнений модели Лоренца:

(%i2) eq: [s*(y-x), x*(r-z) -y, x*y - b*z];

Задаём временные параметры решения

(%i3) t_range: [t,0,50,0.05];

Задаём параметры системы

(%i4) s: 10.0; r: 28.0; b: 2.6667;

Задаём начальные значения 

(%i7) init: [1.0,0,0];

Выполняемрешение

(%i8) sol: rk(eq, [x,y,z],init,t_range)$

(%i9) len:length(sol);

Выделяем компоненты решения и строим графические иллюстрации:

(%i10)     t:makelist(sol[k][1],k,1,len)$

     x:makelist(sol[k][2],k,1,len)$

     y:makelist(sol[k][3],k,1,len)$

     z:makelist(sol[k][4],k,1,len)$

     plot2d([discrete,t,x])$ plot2d([discrete,t,y])$

Результаты решения (хаотические колебания  ) представлен на рис.6.14 и рис.6.15 (фазовый портрет системы). На рисунках объединены в одних осях кривые  .

Рис. 6.14. Пример формирования динамического хаоса (аттрактор Лоренца)

Рис. 6.15. Трехмерный фазовый портрет (аттрактор Лоренца)

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

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

6.3.6 Модель автоколебательной системы: уравнение Ван дер Поля

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

Рис. 6.16. Решение уравнения Ван дер Поля

Решением уравнения Ван дер Поля являются колебания, вид которых для  показан на рис.6.16. Они называются автоколебаниями и принципиально отличаются от рассмотренных ранее (например, численности популяций в модели Вольтерpa) тем, что их характеристики (амплитуда, частота, спектр) не зависят от начальных условий, а определяются исключительно свойствами самой динамической системы. Через некоторое время расчетов после выхода из начальной точки решение выходит на один и тот же цикл колебаний, называемый предельным циклом. Аттрактор типа предельного цикла является замкнутой кривой на фазовой плоскости. К нему асимптотически притягиваются все окрестные траектории, выходящие из различных начальных точек, как изнутри (рис.6.17), так и снаружи предельного цикла.

Рис. 6.17. Фазовый портрет уравнения Ван дер Поля

Использованный командный файл Maxima (для построения графической иллюстрации использован пакет draw):

load("dynamics")$ load("draw")$

mu:1$ s:rk ([v,mu*(1-y^2)*v-y],[y,v],[1,0],[t,0,40,0.2])$

time:makelist(s[k][1],k,1,length(s))$

y:makelist(s[k][2],k,1,length(s))$

v:makelist(s[k][3],k,1,length(s))$

draw2d(points_joined = true, point_type=6, key= "y-v",

     xlabel="y",ylabel="v",points(y,v), terminal=eps)$

Список использованных источников

Документация по текущей версии пакета - 
 http://maxima. sourceforge.net/docs/manual/en/maxima.html

В. А. Ильина, П. К. Силаев Система аналитических вычислений Maxima для физиков-теоретиков.  М.:МГУ им. М. В. Ломоносова, 2007. — 113 с. 

Методическое пособие по изучению математического пакета Maxima. Математический практикум с применением пакета Maxima. 
http://www.pmtf.msiu.ru/chair31/students/spichkov/maxima2.pdf

Книги по Maxima (электронное руководство) 
http://maxima.sourceforge.net/ru/documentation.html

Gilberto E. Urroz  
http://www.neng.usu.edu/cee/faculty/gurro/Maxima.html

Аладьев В. З. Системы компьютерной алгебры: Maple: искусство программирования  М.: Лаборатория базовых знаний, 2006. — 792 с.

Васильев А. Н. Mathcad 13 на примерах  СПб.: БХВ-Петербург, 2006. — 528 с.  

Воробъев Е. М.Введение в систему символьных, графических и численных вычислений Mathematica 5  М.: ДИАЛОГ-МИФИ, 2005. — 368 с.

Цибулин В. Г.,  Говорухин В. Н. Введение в Maple. Математический пакет для всех  М.: Мир, 1997. — 208 с

Гурский  Вычисления в MathCAD  Мн.: Новое знание, 2003. — 814 с.

Гурский Д. А., Турбина Е. Mathcad для студентов и школьников. Популярный самоучитель СПб.: Питер, 2005. — 400 с.

Дьяконов В. П. Maple 9 в математике, физике и образовании  М.: СОЛОН-Пресс, 2004. — 688 с.

Адрес публикации: https://www.prodlenka.org/metodicheskie-razrabotki/413590-modelirovanie-realnyh-processov-v-pakete-rass