python 双摆系统仿真与画图

双摆系统是一个经典的力学问题，它由两个相互连接的摆构成，可以用来模拟物理系统中的振动和运动。在 Python 中，我们可以使用 scipy 库来求解微分方程，并使用 matplotlib 库来绘制图表。

以下是一个简单的双摆系统的仿真和绘图代码示例：

import numpy as np
import scipy.integrate as spi
import matplotlib.pyplot as plt

# 定义常量和初始条件
g = 9.8  # 重力加速度
L1 = 1.0  # 第一个摆的长度
L2 = 0.5  # 第二个摆的长度
m1 = 1.0  # 第一个摆的质量
m2 = 0.5  # 第二个摆的质量
tmax = 10.0  # 最大仿真时间
dt = 0.01  # 时间步长

# 定义双摆系统的微分方程
def equations(state, t):
    theta1, omega1, theta2, omega2 = state
    
    dtheta1 = omega1
    domega1 = (-g*(2*m1+m2)*np.sin(theta1)-m2*g*np.sin(theta1-2*theta2)-2*np.sin(theta1-theta2)*m2*(omega2**2*L2+omega1**2*L1*np.cos(theta1-theta2)))/(L1*(2*m1+m2-m2*np.cos(2*theta1-2*theta2)))
    dtheta2 = omega2
    domega2 = (2*np.sin(theta1-theta2)*(omega1**2*L1*(m1+m2)+g*(m1+m2)*np.cos(theta1)+omega2**2*L2*m2*np.cos(theta1-theta2)))/(L2*(2*m1+m2-m2*np.cos(2*theta1-2*theta2)))
    
    return [dtheta1, domega1, dtheta2, domega2]

# 定义初始条件
theta1_init = 0.2
omega1_init = 0.0
theta2_init = 0.2
omega2_init = 0.0
state_init = [theta1_init, omega1_init, theta2_init, omega2_init]

# 求解微分方程
t = np.arange(0.0, tmax, dt)
state = spi.odeint(equations, state_init, t)

# 绘制摆的轨迹
x1 = L1*np.sin(state[:, 0])
y1 = -L1*np.cos(state[:, 0])
x2 = x1 + L2*np.sin(state[:, 2])
y2 = y1 - L2*np.cos(state[:, 2])

plt.plot(x1, y1, 'r', label='Pendulum 1')
plt.plot(x2, y2, 'b', label='Pendulum 2')
plt.legend(loc='lower right')
plt.xlabel('x')
plt.ylabel('y')
plt.title('Double Pendulum Simulation')
plt.show()

运行此代码，将会得到一个绘制出双摆系统轨迹的图表。

注意，这只是一个简单的示例，实际上在双摆系统的模拟中还需要考虑一些因素，例如阻尼、空气阻力等等。