import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
import math
%matplotlib inline
def osc(X, t):
x, y = X
return np.array([y, (-math.sin(2 * t) + 2) * x])
odeint(osc, np.array([1, 2]), np.array([0, 2]))
def phase_flow(system, t_0, t_1):
def g(x):
return odeint(system, x, np.array([t_0, t_1]))[-1]
return g
g = phase_flow(osc, 1, 4)
e1 = np.array([1, 0])
e2 = np.array([0, 1])
g(e1)
g(e2)
g(e1 + e2)
g(e1) + g(e2)
np.isclose(g(e1 + e2), g(e1) + g(e2)).all()
g(3 * e1)
np.isclose(3 * g(e1), g(3 * e1)).all()
np.linalg.det(np.array([g(e1), g(e2)]))
def linear_system(X, t):
x, y = X
return np.array([x * t + 2 * y,
x * math.sin(t) - 3 * y])
g = phase_flow(linear_system, 0, 3)
np.isclose(g(e1 + e2), g(e1) + g(e2)).all()
Найти $\det g_0^3$ аналитически.
math.exp(3 ** 2 / 2 - 3 * 3)
np.linalg.det(np.array([g(e1), g(e2)]))
g(e1)
g(e2)
def pendulum(X, t):
x, y = X
return np.array([y,
-math.sin(x)])
initial_points = [x for x in np.random.uniform(-3, 3, size=(100000, 2))
if np.linalg.norm(x - np.array([0.8, 0.7])) < 0.4]
T = np.linspace(0, 50)
trajectories = [odeint(pendulum, initial_point, T)
for initial_point in initial_points]
from ipywidgets import interact
@interact(t=(0, 50), manual=True)
def draw_phase_volume(t):
X = np.linspace(-8, 8, 500)
Y = np.linspace(-4, 4, 500)
x, y = np.meshgrid(X, Y)
plt.contour(x, y, y ** 2 / 2 - np.cos(x), levels=20)
plt.plot(*zip(*np.array(trajectories)[:, t, :]), '.')