Решение
import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
from matplotlib.animation import ArtistAnimation
from sympy import *
# Движение плоского маятника заданной длины и массы, точка подвеса которого
# движется по вертикальной окрудности заданного радиуса и частоты,
# в плоскотси движения матяника
# Определяем функции
L = Function('L')
phi = Function('phi')
v_phi = Function('v_phi')
# Определяем переменные
t = Symbol('t')
# Определяем уравнение
omega = Symbol('omega')
l = Symbol('l')
R = Symbol('R')
g = Symbol('g')
m = Symbol('m')
L = m * l**2 / 2 * (v_phi(t))**2 + m * R * l * omega**2 * cos(phi(t) - omega * t) \
+ m * g * l * cos(phi(t))
# Составляем уравнение
Eq = diff(diff(L, v_phi(t)), t) - diff(L, phi(t))
print(Eq)
# Определяем переменную величину
N = 1000
t = np.linspace(0, 31, N)
# Определяем функцию для системы диф. уравнений
def move_func(s, t):
phi, v_phi = s
dphidt = v_phi
dv_phidt = ((R * omega**2 * np.sin(-phi+ omega * t)
- g * np.sin(phi)) / l)
return dphidt, dv_phidt
# Определяем начальные значения и параметры, входящие в систему диф. уравнений
phi0 = np.pi / 180 * 45
v_phi0 = 1
s0 = phi0, v_phi0
g = 9.8
R = 0.5
l = 0.6
omega = 0.8
# Решаем систему диф. уравнений
sol = odeint(move_func, s0, t)
# Возвращаемся в декартовы координаты
x1 = R * np.sin(omega * t[:])
y1 = - R * np.cos(omega * t[:])
x2 = x1 + l * np.sin(sol[:, 0])
y2 = y1 - l * np.cos(sol[:, 0])
# Строим решение в виде графика и анимируем
fig = plt.figure()
plt.style.use(['dark_background'])
fig = plt.figure(figsize=(12,12))
bodys = []
for i in range(1, len(t), 2):
thisx = [0, x1[i], x2[i]]
thisy = [0, y1[i], y2[i]]
body_line, = plt.plot(thisx, thisy, '-o', color='white', ms = 10)
bodys.append([body_line])
ani = ArtistAnimation(fig, bodys, interval=50)
plt.axis('off')
# ani.save('Движение маятника.gif')
plt.show()