N Bob Pendlum with parallel computing
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation
from concurrent.futures import ThreadPoolExecutor
# Constants
L = 1.0 # Length of each pendulum segment
M = 1.0 # Mass of each pendulum segment
# Derivatives function
def derivs(state, num_pendulums, gravity):
dydx = np.zeros_like(state)
angles = state[:num_pendulums]
velocities = state[num_pendulums:]
for i in range(num_pendulums):
if i == 0:
dydx[i] = velocities[i]
dydx[i + num_pendulums] = -gravity / L * np.sin(angles[i])
else:
dydx[i] = velocities[i]
dydx[i + num_pendulums] = -gravity / L * (np.sin(angles[i]) - np.sin(angles[i - 1]))
return dydx
# Integration function
def integrate(state, dt, num_pendulums, gravity):
return state + dt * derivs(state, num_pendulums, gravity)
# Input from the user
num_pendulums = int(input("Enter the number of pendulums: "))
gravity = float(input("Enter the gravitational acceleration (m/s^2): "))
# Initial conditions
initial_angles = [np.pi / 2 for _ in range(num_pendulums)]
initial_velocities = [0 for _ in range(num_pendulums)]
state = np.array(initial_angles + initial_velocities, dtype='float64')
dt = 0.01
t = np.arange(0, 20, dt)
# Animation setup
fig, ax = plt.subplots(figsize=(8, 8))
max_length = num_pendulums * L
ax.set_xlim(-max_length, max_length)
ax.set_ylim(-max_length, max_length)
# Adjust line and bob thickness based on the number of pendulums
line_thickness = max(1, 6 - num_pendulums)
bob_size = max(1, 10 - num_pendulums)
lines, = ax.plot([], [], 'o-', lw=line_thickness, markersize=bob_size)
# Initialize the plot
def init():
lines.set_data([], [])
return lines,
# Update function for animation
def update(frame):
global state
with ThreadPoolExecutor() as executor:
futures = [executor.submit(integrate, state, dt, num_pendulums, gravity) for _ in range(1)]
for future in futures:
state = future.result()
x = np.cumsum([L * np.sin(angle) for angle in state[:num_pendulums]])
y = np.cumsum([-L * np.cos(angle) for angle in state[:num_pendulums]])
lines.set_data([0] + x.tolist(), [0] + y.tolist())
return lines,
ani = FuncAnimation(fig, update, frames=len(t), init_func=init, blit=True, interval=10)
plt.show()