# Bead on a Wire
# Install ffmpeg (allows us to create videos)
# deepnote.com/@jeh15
!apt update -y
!apt install ffmpeg -y
# Import Libraries
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FFMpegWriter
from matplotlib.patches import Circle
# Parameters
g = 9.8
m = 1
# We need an array of time points from 0 to 10 in increments of 0.01 seconds
dt = 0.001
t_vec = np.arange(0,10,dt)
# Initialize a vector of zeros
x_vec = np.zeros(len(t_vec))
y_vec = np.zeros(len(t_vec))
dx_vec = np.zeros(len(t_vec))
dy_vec = np.zeros(len(t_vec))
lam_vec = np.zeros(len(t_vec))
# Set our initial condition
x_vec[0] = 1 # initial x position
dx_vec[0] = 0 # initial angular velocity
# Set the initial y values
y_vec[0] = x_vec[0] ** 2 # initial y position
dy_vec[0] = 0 # initial angular velocity
# Initialize Matrix and Constants Ax = B
A = np.array([[m, 0, 0],[0, m, -1],[0, -1, 0]], dtype = float)
B = np.array([[0],[-m*g],[0]], dtype = float)
# Loop through time
# Euler's Method (approximately integrates the differential equation)
for i in range(1, len(t_vec)):
# for i in range(1, 10):
# Update the Matrix
A[0,2] = 2 * x_vec[i-1]
A[2,0] = 2 * x_vec[i-1]
B[2,0] = -2 * dx_vec[i-1] ** 2
# Solve the equation Ax = B
[ddx, ddy, lam] = np.linalg.solve(A, B)
# Integrate Using Eulers Method
x_vec[i] = x_vec[i-1] + dx_vec[i-1]*dt
y_vec[i] = y_vec[i-1] + dy_vec[i-1]*dt
dx_vec[i] = dx_vec[i-1] + ddx*dt
dy_vec[i] = dy_vec[i-1] + ddy*dt
lam_vec[i] = lam
# Setup our figure for drawing
# fig, axs = plt.subplots(3)
# fig.suptitle('Positions and Forces')
# axs[0].plot(t_vec, x_vec)
# axs[0].set_ylabel('X Position')
# axs[1].plot(t_vec, y_vec)
# axs[1].set_ylabel('Y Position')
# axs[2].plot(t_vec, lam_vec)
# axs[2].set_ylabel('Lambda Value')
# Create a plot on those axes which is currently empty
fig, ax = plt.subplots()
p, = ax.plot([],[],color='black') # Initializes an empty plot
ax.set_xlim([-2, 2]) # x limits
ax.set_ylim([-1, 3]) # y limits
ax.axis('equal')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_title('Bead on Wire Simulation')
video_title = "simulation"
# Now we want add the circle to the plot
c = Circle((0,0), radius=0.04, color='cornflowerblue')
ax.add_patch(c) # puts the circle on that set of axes
# Now we want add the wire
x = np.arange(-2,2,0.1)
wire = x, x**2
p.set_data(x, x**2)
# Define Information for our animation
FPS = 20
sample_rate = int(1/(FPS*dt))
vec_size = len(t_vec)
dpi = 300 # Quality of the video
writerObj = FFMpegWriter(fps = FPS)
# Plot and Create Animation
with writerObj.saving(fig, video_title+".mp4", dpi):
# Loop through the frames we want to animate
for i in range(0,vec_size,sample_rate):
# Update the bob circle
patch_center = x_vec[i], y_vec[i]
# NOTE: commas w/o bracket create a list
# update the circle
c.center = patch_center
# Update the Drawing
fig.canvas.draw()
# Grab and save the frame
writerObj.grab_frame()
# Import the video you just made and run it on the notebook
from IPython.display import Video
Video("/work/simulation.mp4", embed=True, width=640, height=480)