Астрозадача №6: 

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()