Как рисовать фазовый портрет
Перейти к содержимому

Как рисовать фазовый портрет

  • автор:

Построение фазовых траекторий

Фазовая плоскость — координатная плоскость, в которой по осям координат откладываются какие-либо две переменные (фазовые координаты), однозначно определяющие состояние системы второго порядка.Каждая точка фазовой плоскости отражает одно состояние системы и называется фазовой, изображающей или представляющей точкой. Изменение состояния системы отображается на фазовой плоскости движением этой точки. След от движения изображающей точки называется фазовой траекторией.

Уравнение состояния системы имеет вид:

[math]m\ddot x = -kx — μ\dot x [/math]

корни характеристического уравнения:

[math]mλ^2 + μλ + k = 0[/math]

Фазовые портреты зависят от параметров m, k, μ.

Частные случаи [ править ]

  • При [math]μ^2 \gt 4mk\ [/math] (трение велико) оба корня характеристического уравнения действительные, различные и отрицательные. Фазовый портрет типа «Устойчивый узел»
  • При [math]μ^2 = 4mk\ [/math] корень действительный отрицательный. Фазовый портрет типа «Вырожденный устойчивый узел»
  • При [math] 0 \lt μ^2 \lt 4mk\ [/math] (трение мало) имеем два комплексно сопряженных корня с отрицательный вещественной частью. Фазовый портрет типа «Фокус»
  • При [math]μ = 0, k \ne 0 [/math] (трение отсутствует) оба корня характеристического уравнения действительные, различные и отрицательные. Фазовый портрет типа «Центр»
  • При [math]μ \ne 0, k = 0\ [/math] (Упругая сила отсутсвует) оба корня действительные, причем один отрицательный, а другой равен 0.

Реализация [ править ]

  • за основу была взята программа Суперэллипс и суперпарабола

Фазовые портреты «на пальцах» или что можно узнать о решениях диффура, не решая его

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

Встречайте: фазовые портреты (они же фазовые диаграммы). Простым языком, фазовый портрет — это то, как величины, описывающие состояние системы (a.k.a. динамические переменные), зависят друг от друга. В случае механического движения это координата и скорость, в электричестве это заряд и ток, в известной популяционной задаче это количество хищников и жертв и т.д.

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

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

Где x — угол отклонения стержня от вертикали, точка над x означает производную по времени, а коэффициент перед синусом зависит от размера и массы стержня.

Если амплитуда (размах) колебаний достаточно мала, синус можно приближенно заменить его аргументом (вы ведь помните первый замечательный предел, нет?). В таком случае, уравнение принимает следующий вид:

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

Получилось выражение, первый член которого выглядит как кинетическая энергия. Это не случайно — на самом деле мы получили именно закон сохранения энергии. Постоянная Е в правой части (полная энергия системы на единицу массы) может принимать различные значения, которые соответствуют разным начальным состояниям системы.

Полученный нами закон сохранения превратился в уравнение кривой на плоскости (x,u):

Для разных значений Е мы получим разные кривые. Нарисуем несколько таких линий для разных значений энергии:

По горизонтальной оси отложена величина x, по вертикальной — u

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

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

Следует однако помнить, что замена синуса его аргументом оправдана лишь для малых углов отклонения (от 10 градусов и меньше), поэтому мы не можем доверять тем траекториям, которые выходят за границы области, ограниченной жирными пунктирными линиями, то есть из четырех приведенных траекторий лишь оранжевая достоверно отображает реальность. Кроме того, поскольку x это угол, то его значения, соответствующие 180 и -180 градусам описывают одно и то же положение стержня, то есть правая и левая пунктирные линии (тонкие) на графике это на самом деле одна и та же линия.

Теперь, поскольку нам понятна суть, можно перейти к чему-то посложнее. Выше мы очень сильно упростили уравнение и при этом ограничили себя только малыми колебаниями. Математик бы сказал, что мы линеаризовали уравнение и пренебрегли нелинейными эффектами. Так давайте включим в рассмотрение нелинейность. Вернёмся к самому первому уравнению — с синусом. Если мы повторим с ним то, что проделали с линейным уравнением, мы получим следующий закон сохранения:

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

По горизонтальной оси отложена величина x, по вертикальной — u

Как видите, процессы происходящее в системе стали более разнообразными:

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

При больших энергиях (зеленая траектория) колебаний уже нет, вместо этого мы получаем вращательное движение с переменной скоростью. И действительно, если достаточно сильно «толкнуть» стержень, он будет вращаться, замедляясь при подъёме и ускоряясь при спуске.

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

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

По горизонтальной оси отложено время, по вертикальной — x

По горизонтальной оси отложено время, по вертикальной — x

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

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

  • математика
  • динамические системы
  • качественный анализ
  • диффуры

Фазовый портрет

Чтобы научиться рисовать и без труда понимать фазовые портреты, рассмотрим сначала совсем простые задачи. Пусть точка равномерно движется по прямой и в начальный момент t = 0 ее координата s равна нулю, так что s = v0t. График этого движения — прямая линия с наклоном, пропорциональным скорости (рис. 4.8, ?).

Если на вертикальной оси графика 1 см соответствует 1 с времени, а 1 см по горизонтали соответствует 1 см пути, то скорость, очевидно, равна tg ? см/с. Дальше мы не будем упоминать об этом соглашении и с производными величинами будем обращаться точно так же (на оси скорости 1 см соответствует скорость 1 см/с и т. д.). Отрицательным значениям угла отвечают движения в отрицательном направлении по оси Os (рис. 4.8, б).

Нетрудно нарисовать любую фазовую траекторию. Это просто прямые, параллельные горизонтальной оси Os и пересекающие вертикальную ось в точке, соответствующей значению скорости, равному v0. Когда скорость положительна, изображающая точка А пробегает фазовую траекторию слева направо, при отрицательной скорости — в обратном направлении. Если s = s0 + v0t, то график движения не проходит через точку О, но фазовая траектория такого движения совпадает с фазовой траекторией движения s = v0t. Это, конечно, легко проверить, но на самом деле это должно быть очевидным, так как фазовые траектории не зависят от момента t0, в котором мы начинаем отсчет времени.

Если точка покоится, то на графике движения ей соответствует прямая, параллельная оси времени, т. е. s = s0 и ? = 0. На фазовой диаграмме этой прямой соответствует точка s = s0 на оси Os, т. е. точка (s, v) = (s0, 0). При разных значениях s0 эти точки заполняют всю ось Os. Каждую точку оси Os нужно рассматривать как отдельную фазовую траекторию.

Таким образом, фазовые траектории точки, движущейся равномерно по прямой, — это прямые, параллельные оси Os, а также точки оси Os. Через каждую точку фазовой плоскости (s, v) проходит только одна фазовая траектория, если договориться, что выбор начала отсчета времени t0 несуществен (т. е. важно лишь, какую кривую пробегает изображающая точка). Чтобы больше не возвращаться к этому, можно, как это делалось и раньше, условиться, что s = 0 при t = 0, а остальные движения получать сдвигом начала отсчета времени.

В качестве упражнения постройте фазовые диаграммы равномерно ускоренных движений грузика, падающего с высоты h или подбрасываемого вверх. Точка О на фазовой диаграмме представляет фазовую траекторию лежащего на земле грузика. Вообще, такие точки на фазовых диаграммах называются точками покоя. В самом нижнем положении наш грузик покоится устойчиво, иными словами, точка на фазовой диаграмме — устойчивая точка покоя. Если грузик слегка подбросить, он вернется назад. Дальнейшее движение грузика зависит от его устройства как реальной физической системы. Если грузик — модель упругого мячика, падающего на асфальт, то он будет отскакивать, пока вся его энергия не перейдет в тепло (попробуйте нарисовать фазовые траектории этих движений). Если же уронить на пол кусочек пластилина, то он останется в нижнем положении (какова фазовая траектория в этом случае?).

Точки покоя на фазовом портрете равномерно движущегося грузика, наоборот, неустойчивы. Если сообщить грузику небольшой импульс, то он начнет равномерно двигаться и в конце концов уйдет сколь угодно далеко от исходного положения. На фазовом портрете это будет выглядеть так, что точка (s0, 0) «перепрыгнет» на близкую фазовую траекторию и уйдет по ней сколь угодно далеко. В реальной системе (скажем, шайба на льду) этому помешает трение, но при очень малом трении шайба все равно улетит далеко, а при достаточном заметном трении нужно уже рисовать другой фазовый портрет, так как фазовые траектории не будут прямыми, параллельными оси Os (подумайте, как они могут выглядеть).

Digiratory

Лаборатория автоматизации и цифровой обработки сигналов

Построение фазовых портретов на языке Python

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

С помощью фазовых портретов можно синтезировать регуляторы (Метод фазовой плоскости) или проводить анализ положений устойчивости и характера движений системы.

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

В качестве примера воспользуемся моделью маятника с вязким трением:

Где \(\omega \) — скорость, \(\theta \) — угол отклонения, \(b \) — коэффициент вязкого трения, \(с \) — коэффициент, учитывающий массу, длину и силу тяжести.

Для работы будем использовать библиотеки numpy, scipy и matplotlib для языка Python.

Блок импорта выглядит следующим образом:

import numpy as np from scipy.integrate import odeint import matplotlib.pyplot as plt

Шаг 1. Реализация ОДУ в Python

Определим функцию, отвечающую за расчет ОДУ. Например, следующего вида:

def ode(y, t, b, c): theta, omega = y dydt = [omega, -b*omega - c*np.sin(theta)] return dydt

Аргументами функции являются:

  • y — вектор переменных состояния
  • t — время
  • b, c — параметры ДУ (может быть любое количество)

Функция возвращает вектор производных.

Шаг 2. Численное решение ОДУ

Далее необходимо реализовать функцию для получения решения ОДУ с заданными начальными условиями:

def calcODE(args, y0, dy0, ts = 10, nt = 101): y0 = [y0, dy0] t = np.linspace(0, ts, nt) sol = odeint(ode, y0, t, args) return sol

Аргументами функции являются:

  • args — Параметры ОДУ (см. шаг 1)
  • y0— Начальные условия для первой переменной состояния
  • dy0 — Начальные условия для второй переменной состояния (или в нашем случае ее производной)
  • ts — длительность решения
  • nt — Количество шагов в решении (= время интегрирования * шаг времени)

В 3-й строке формируется вектор временных отсчетов. В 4-й строке вызывается функция решения ОДУ.

Шаг 3. Генерация и вывод фазового портрета

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

def drawPhasePortrait(args, deltaX = 1, deltaDX = 1, startX = 0, stopX = 5, startDX = 0, stopDX = 5, ts = 10, nt = 101): for y0 in range(startX, stopX, deltaX): for dy0 in range(startDX, stopDX, deltaDX): sol = calcODE(args, y0, dy0, ts, nt) plt.plot(sol[:, 1], sol[:, 0], 'b') plt.xlabel('x') plt.ylabel('dx/dt') plt.grid() plt.show()

Аргументами функции являются:

  • args — Параметры ОДУ (см. шаг 1)
  • deltaX — Шаг начальных условий по горизонтальной оси (переменной состояния)
  • deltaDX — Шаг начальных условий по вертикальной оси (производной переменной состояния)
  • startX — Начальное значение интервала начальных условий
  • stopX — Конечное значение интервала начальных условий
  • startDX — Начальное значение интервала начальных условий
  • stopDX — Конечное значение интервала начальных условий
  • ts — длительность решения
  • nt — Количество шагов в решении (= время интегрирования * шаг времени)

Во вложенных циклах (строки 3-4) происходит перебор начальных условий дифференциального уравнения. В теле этих циклов (строки 5-6) происходит вызов функции решения ОДУ с заданными НУ и вывод фазовая траектории полученного решения.

Далее производятся нехитрые действия:

  • Строка 7 — задается название оси X
  • Строка 9 — задается название оси Y
  • Строка 10 — выводится сетка на графике
  • Строка 11 — вывод графика (рендер)

Шаг 4. Запуск построения

Запустить построение можно следующим образом:

b = 0.25 c = 5.0 args=(b, c) drawPhasePortrait(args) drawPhasePortrait(args, 1, 1, -10, 10, -5, 5, ts = 30, nt = 301)

Строка 1-2 — задание значений параметрам ОДУ

Строка 3 — формирование вектора параметров

Строка 4 — вызов функции генерации фазового портрета с параметрами «по умолчанию»

Строка 5 — вызов функции генерации фазового портрета с настроенными параметрами

Итог

При запуске программы получаем следующий результат:

Полный текст программы под лицензией MIT (Использование при условии ссылки на источник):

# -*- coding: utf-8 -*- """ This function create a phase portrait of ode Created on Mon Jan 23 18:51:01 2017 Copyright 2017 Sinitsa AM site: digiratory.ru Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions: The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software. THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. @author: Sinitsa AM digiratory.ru """ import numpy as np from scipy.integrate import odeint import matplotlib.pyplot as plt ''' In this function you should implement your ode for example: d*theta/dt = omega d*omega/dt = -b*omega - c*np.sin(theta) @args: y - state variable t - time b, c - args of ode ''' def ode(y, t, b, c): theta, omega = y dydt = [omega, -b*omega - c*np.sin(theta)] return dydt ''' Calculate ode @args: args - arguments of ODE y0 - The initial state of y yd0 - The initial state of derivative y ts - the duration of the simulation nt - Number steps if simulation (=ts*deltat) ''' def calcODE(args, y0, dy0, ts = 10, nt = 101): y0 = [y0, dy0] t = np.linspace(0, ts, nt) sol = odeint(ode, y0, t, args) return sol ''' Drawing Phase portrait of ODE in function ode() @args: args - arguments of ODE deltaX - step x deltaDX - step derivative x startX - start value of x stopX - stop value of x startDX - start value of derivative x stopDX - stop value of derivative x ts - the duration of the simulation nt - Number steps if simulation (=ts*deltat) ''' def drawPhasePortrait(args, deltaX = 1, deltaDX = 1, startX = 0, stopX = 5, startDX = 0, stopDX = 5, ts = 10, nt = 101): for y0 in range(startX, stopX, deltaX): for dy0 in range(startDX, stopDX, deltaDX): sol = calcODE(args, y0, dy0, ts, nt) plt.plot(sol[:, 0], sol[:, 1], 'b') plt.xlabel('x') plt.ylabel('dx/dt') plt.grid() plt.show() b = 0.25 c = 5.0 args=(b, c) drawPhasePortrait(args) drawPhasePortrait(args, 1, 1, -10, 10, -5, 5, ts = 30, nt = 301)

Синица А.М.: Построение фазовых портретов на языке Python [Электронный ресурс] // Digiratory. 2017 г. URL: https://digiratory.ru/?p=435 (дата обращения: ДД.ММ.ГГГГ).

Построение фазовых портретов на языке Python : 2 комментария

  1. Уведомление: Устойчивость нелинейных систем — Digiratory
  2. Уведомление: Оценка параметров ДУ в Python — Digiratory

Добавить комментарий

Ваш адрес email не будет опубликован. Обязательные поля помечены *