Автор: | Дмитрий Шачнев |
---|---|
Группа: | 411 |
Дата: | 2013-10-10 |
Дана система дифференциальных уравнений в \mathbb{R}^2:
\left\{ \eqalign{ \dot{x} &= P(x, y) &= y \\ \dot{y} &= Q(x, y) &= ky - x - y^3 } \right.
При значениях параметра k = 1, -1, 3 требуется:
- Найти собственные значения и собственные вектора матрицы Якоби.
- Найти особую точку и определить её тип.
- Найти длину периода, если существует периодическое решение.
- Нарисовать фазовый портрет системы.
Из уравнения y = ky - x - y^3 = 0 видно, что единственная особая точка — это (0, 0). Чтобы найти её тип, надо сначала вычислить собственные значения матрицы Якоби:
\det\left( \begin{matrix} -\lambda & 1 \\ -1 & k - \lambda \end{matrix} \right) = 0
\lambda^2 - k\lambda + 1 = 0
\lambda = \frac{k \pm \sqrt{k^2 - 4}}{2}
При k = 1 имеем \lambda = \frac{1 \pm \sqrt{3}i}{2}, численные значения 0.5 + 0.866i и 0.5 - 0.866i, собственные вектора — (1, \frac{1 \pm \sqrt{3}i}{2}).
Корни комплексные с положительной вещественной частью, то есть особая точка — неустойчивый фокус.
Периодическое решение существует, длина периода равна 7.7.
При k = -1 имеем \lambda = \frac{-1 \pm \sqrt{3}i}{2}, численные значения -0.5 + 0.866i и -0.5 - 0.866i, собственные вектора — (1, \frac{-1 \pm \sqrt{3}i}{2}).
Корни комплексные с отрицательной вещественной частью, то есть особая точка — устойчивый фокус.
Периодического решения не существует.
При k = 3 имеем \lambda = \frac{3 \pm \sqrt{5}}{2}, численные значения 2.618 и 0.382, собственные вектора — (1, \frac{3 \pm \sqrt{5}}{2}).
Корни действительные и положительные, то есть особая точка — неустойчивый узел.
Периодическое решение существует, длина периода равна 17.1.
В случаях k = 1, 3 имеет место предельный цикл, на который «накручиваются» траектории.
Граничные значения: k = 0 — центр, k = 2 — нестабильный вырожденный узел.
Для построения траектории из точки (x, y) сначала выбирается некоторый луч \mathcal{L}, выходящий из особой точки (например, ось Ox).
Далее, для построения каждого шага траектории используется метод Рунге–Кутты. Ниже описан алгоритм метода Рунге–Кутты четвёртого порядка точности:
\left\{ \eqalign{ P_1 &= P \left( x + \frac{\tau}{2} \cdot P(x, y),\; y + \frac{\tau}{2} \cdot Q(x, y) \right), \\ Q_1 &= Q \left( x + \frac{\tau}{2} \cdot P(x, y),\; y + \frac{\tau}{2} \cdot Q(x, y) \right) } \right.
\left\{ \eqalign{ P_2 &= P \left( x + \frac{\tau}{2} \cdot P_1,\; y + \frac{\tau}{2} \cdot Q_1 \right), \\ Q_2 &= Q \left( x + \frac{\tau}{2} \cdot P_1,\; y + \frac{\tau}{2} \cdot Q_1 \right) } \right.
\left\{ \eqalign{ P_3 &= P ( x + \tau \cdot P_2,\; y + \tau \cdot Q_2 ), \\ Q_3 &= Q ( x + \tau \cdot P_2,\; y + \tau \cdot Q_2 ) } \right.
\left\{ \eqalign{ x' &= x + \frac{\tau}{6} \left( P(x, y) + 2 \cdot P_1 + 2 \cdot P_2 + P_3 \right), \\ y' &= y + \frac{\tau}{6} \left( Q(x, y) + 2 \cdot P_1 + 2 \cdot P_2 + P_3 \right) } \right.
Соединив (x, y) и (x', y'), получим шаг траектории. Такие шаги будем продолжать делать до тех пор, пока не будет выполнено одно из трёх условий конца цикла:
- Траектория пересечёт некоторый луч два раза в точках, отличающихся не более чем на \frac{\tau}{2} (циклическое решение).
- Тректория подойдёт к особой точке на расстояние \frac{\tau}{2} (стремление к особой точке).
- Траектория отойдёт от начала координат на расстояние больше \frac{10}{\tau}, либо будет достигнуто максимальное число шагов \frac{50}{\tau} (стремление к бесконечности).
Пусть задан луч \mathcal{L} и точка (x_0, y_0) = (x_0(\rho), y_0(\rho)) \in \mathcal{L} (где \rho — локальная координата на луче \mathcal{L}).
Пустим из этой точки траекторию движения системы в положительном времени, и пусть (x_1, y_1) = (x_1(\chi), y_1(\chi)) — первое пересечение этой траектории с лучом \mathcal{L}. Отображение \rho \to \chi = \chi(\rho) называется отображением Пуанкаре этого луча.
Обозначим \psi(\rho) = \chi(\rho) - \rho. Для нахождения неподвижной точки отображения \chi (то есть нуля \psi) используется упрощённый метод Ньютона:
Выбираются две точки a и b, такие, что \psi(a) \cdot \psi(b) < 0.
Пока расстояние \left|a - b\right| не достигнет необходимой минимальной величины, выполняется замена
\eqalign{ a &\to b - (b - a) \cdot \frac{\psi(b)}{\psi(b) - \psi(a)}, \\ b &\to a - (a - b) \cdot \frac{\psi(a)}{\psi(a) - \psi(b)} }
Для поиска периода из неподвижной точки отображения Пуанкаре выпускается траектория с шагом \tau. Длиной периода будет N \cdot \tau, где N — число шагов, за которое алгоритм Рунге–Кутты вернётся в исходную точку.
Определим погрешность метода Рунге–Кутты с шагом \tau, запущенного из начальной точки (x_0, y_0). Для этого из той же точки запустим алгоритм с шагом \frac{\tau}{2}.
Пусть (x_i, y_i) и (\tilde{x}_i, \tilde{y}_i) — точки с номером i первой и второй траекторий соответственно. Тогда погрешность метода вычисляется по формуле
\varepsilon = \max_{0 \leqslant i \leqslant N} \sqrt{ \left( x_i - \tilde{x}_{2i} \right)^2 + \left( y_i - \tilde{y}_{2i} \right)^2 }.
Неподвижная точка существует, если существует решение, пересекающее некоторый луч два раза в точках, остоящих друг от друга не более, чем на \frac{\tau}{2}. Тогда расстояние между этими точками и будет погрешностью поиска неподвижной точки \delta.
\tau | \varepsilon | \delta |
---|---|---|
0.0001 | 1.09613e-07 | 1.76144e-11 |
0.0002 | 7.7602e-07 | 3.24754e-11 |
0.0005 | 9.25363e-06 | 1.72758e-10 |
0.001 | 5.59294e-05 | 1.46514e-10 |
0.002 | 0.000129327 | 8.42441e-10 |
0.005 | 0.000938488 | 5.29843e-06 |
0.01 | 0.00334614 | 3.76117e-05 |
0.02 | 0.0130162 | 0.000242305 |
0.05 | 0.0885789 | 0.000137598 |
0.1 | 0.376933 | 0.00013097 |
\tau | \varepsilon | \delta |
---|---|---|
0.0001 | 1.31033e-10 | 8.31218e-09 |
0.0002 | 3.36056e-10 | 3.32489e-08 |
0.0005 | 3.96062e-10 | 2.07774e-07 |
0.001 | 8.73237e-10 | 8.30426e-07 |
0.002 | 1.72283e-09 | 3.30989e-06 |
0.005 | 2.09922e-09 | 2.0184e-05 |
0.01 | 3.59741e-09 | 7.44739e-05 |
0.02 | 4.21477e-09 | 0.00021614 |
0.05 | 3.99769e-08 | 0.000112655 |
0.1 | 6.56572e-07 | 0.000268919 |
\tau | \varepsilon | \delta |
---|---|---|
0.0001 | 2.6328e-07 | 1.73772e-11 |
0.0002 | 1.59201e-06 | 7.64651e-11 |
0.0005 | 1.43449e-05 | 3.25251e-11 |
0.001 | 8.42552e-05 | 4.14602e-12 |
0.002 | 0.000578843 | 1.51516e-09 |
0.005 | 0.00145192 | 5.17844e-11 |
0.01 | 0.010968 | 1.04374e-08 |
0.02 | 0.0388027 | 8.47973e-07 |
0.05 | 0.233205 | 3.07063e-06 |
0.1 | 1.00551 | 9.50295e-05 |
Программа написана на языке C++ с использованием фреймворка Qt.
Система уравнений считывается из файла, в котором хранятся коэффициенты многочленов \dot{x} и \dot{y}.
Реализованы два варанта отрисовки:
- В первом пользователю предоставляется возможность задать начальную точку, из которой запускается траектория.
- Во втором фазовый портрет рисуется в виде «сетки» из верхней и нижней сторон области рисования.
Основным объектом является система типа PolynomSystem
, для
которой реализованы следующие методы:
/* Получить следующую точку траектории, с шагом eps */
Point getNextValue(Point point, double eps);
/* Отображение Пуанкаре в точке p, с шагом eps */
double poincareFunction(double p, double eps);
/* Найти неподвижную точку отображения Пуанкаре
* с начальными точками a и b и шагом eps */
double findPoincareStaticPoint(Point a, Point b, double eps);
/* Найти собственные значения матрицы Якоби */
complex[2] getEigenValues();
/* Определить тип особой точки */
PointType getPointType();