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]))
array([[ 1. , 2. ], [14.64673206, 21.63805004]])
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)
array([36.77682016, 42.06341281])
g(e2)
array([29.09520864, 33.3047762 ])
g(e1 + e2)
array([65.87203428, 75.36819653])
g(e1) + g(e2)
array([65.8720288 , 75.36818901])
np.isclose(g(e1 + e2), g(e1) + g(e2)).all()
True
g(3 * e1)
array([110.33046451, 126.19024232])
np.isclose(3 * g(e1), g(3 * e1)).all()
True
np.linalg.det(np.array([g(e1), g(e2)]))
0.9999930497444511
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()
True
Найти $\det g_0^3$ аналитически.
math.exp(3 ** 2 / 2 - 3 * 3)
0.011108996538242306
np.linalg.det(np.array([g(e1), g(e2)]))
0.011194714438770365
g(e1)
array([196.84929067, 9.89932211])
g(e2)
array([115.89484579, 5.82827401])
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, :]), '.')
interactive(children=(IntSlider(value=25, description='t', max=50), Output()), _dom_classes=('widget-interact'…