!apt update -y
!apt install ffmpeg -y
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FFMpegWriter
from matplotlib.patches import Circle
# Defining our Parameters
g = 9.8 # Gravity (m/s^2)
L = 1 # Length (m)
m_1 = 1 # In Kg
m_2 = 2 # In Kg
k = 10 # Spring Constant (N/m)
h = 1.5 # Height of Hinge from ground (m)
# Creating the time array for 10 second run time
dt = 0.001 # Seconds (s)
t_vec = np.arange(0,10,dt) # Time Vector
# Define vectors that will hold our state variables over time
# x Vectors of Zero
x_vec = np.zeros(len(t_vec))
dx_vec = np.zeros(len(t_vec))
ddx_vec = np.zeros(len(t_vec))
# theta Vectors of Zero
theta_vec = np.zeros(len(t_vec))
dtheta_vec = np.zeros(len(t_vec))
ddtheta_vec = np.zeros(len(t_vec))
# Vector for hinge location
y = np.zeros(len(t_vec))
# Establishing the initial conditions for x & theta
x_vec[0] = 0 # x Initial Position
dx_vec[0] = 0 # dx Initial Velocity
theta_vec[0] = 1 # theta Initial Position
dtheta_vec[0] = 0 # dtheta Initial Velocity
# Loop through time
# Creating Loop to run through the theta direction
# The Code works once you remove "...", I could not figure out text wrapping
for i in range(1, len(t_vec)):
theta_vec[i] = theta_vec[i-1] + dtheta_vec[i-1]*dt
ddtheta_vec[i-1] = -((g*np.sin(theta_vec[i-1]))/(L)) - (k*(np.sqrt((L*np.sin(theta_vec[i-1]) - x_vec[i-1])**2 + (h - L*np.cos(theta_vec[i-1]))**2) - L)*((L*np.sin(theta_vec[i-1]) - x_vec[i-1])*(L*np.cos(theta_vec[i-1])) + (h - L*np.cos(theta_vec[i-1]))*(L*np.sin(theta_vec[i-1]))))/(m_1*L**2*np.sqrt((L*np.sin(theta_vec[i-1]) - x_vec[i-1])**2 + (h - L*np.cos(theta_vec[i-1]))**2))
dtheta_vec[i] = dtheta_vec[i-1] + ddtheta_vec[i-1]*dt
# Creating Loop to run through the x direction
for i in range(1, len(t_vec)):
x_vec[i] = x_vec[i-1] + dx_vec[i-1]*dt
# The Code works once you remove "...", I could not figure out text wrapping
ddx_vec[i-1] = (k*(np.sqrt((L*np.sin(theta_vec[i-1]) - x_vec[i-1])**2 + (h - L*np.cos(theta_vec[i-1]))**2) - L)*(L*np.sin(theta_vec[i-1]) - x_vec[i-1]))/(m_2*np.sqrt((L*np.sin(theta_vec[i-1]) - x_vec[i-1])**2 + (h - L*np.cos(theta_vec[i-1]))**2))
dx_vec[i] = dx_vec[i-1] + ddx_vec[i-1]*dt
y[i] = -1.5
# Setting up for the x component on graph
plt.plot(t_vec, x_vec)
plt.title('X - component')
plt.show()
# Setting up for the y component on graph
plt.plot(t_vec, theta_vec)
plt.title('Theta - component')
plt.show()
# SETTING UP THE SIMULATION
# Setup Figure: Initialize Figure / Axe Handles
fig, ax = plt.subplots() # Initialize Figure and Axis
p1, = ax.plot([], [], color='dodgerblue') # Initialize Empty Plot
p2, = ax.plot([], [], color='dodgerblue') # Initialize Empty Plot
ax.axis('equal')
# Setting Parameters for plot
ax.set_xlim([-3, 3])
ax.set_ylim([-3, 3])
# Setting up labels for plot
ax.set_xlabel('X') # X Label
ax.set_ylabel('Y') # Y Label
ax.set_title('Pendulum Simulation:')
video_title = "simulation 2C" # Naming Of Animation
# Initialize Patch: (Circle From Pendelum)
# Draws a circle of on the plot at the origin
c = Circle((0, -0.75), radius=0.1, color='dodgerblue')
ax.add_patch(c) # Adds Patch to PLot
# Initialize Patch: (Rectangle on Ground)
# Couldnt get right positioning so used a circle instead
#width = .5 # Width of rectangle
#height = .3 # Hieght of rectangle
#xy_rec = (x_vec[0]+ width, -h + height)
#b = Rectangle(xy_rec, width,height, color='mediumseagreen')
#ax.add_patch(b) # Adds Patch to PLot
# Implemented circle to show movements
b = Circle((0, -1.5), radius=0.2, color='dodgerblue')
ax.add_patch(b)
# Setup Animation Writer:
FPS = 20 # Frames per second for the video
# Rate at which we sample the simulation to visualize
sample_rate = int(1 / (dt * FPS))
dpi = 300 # Quality of the Video
writerObj = FFMpegWriter(fps=FPS) # Video Object
# We need to calculate the cartesian location of the pendulum:
# Initialize Array:
simulation_size = len(t_vec) # We use this alot, so lets store it.
x_pendulum_arm = np.zeros(simulation_size)
y_pendulum_arm = np.zeros(simulation_size)
origin = (0, 0) # Pin Joint (X, Y) Location
# Find the X and Y trajectories of the pendulum:
for i in range(0, simulation_size):
x_pendulum_arm[i] = L * np.sin(theta_vec[i]) + origin[0]
y_pendulum_arm[i] = origin[1] - L * np.cos(theta_vec[i])
# Plot and Create Animation:
with writerObj.saving(fig, video_title+".mp4", dpi):
# We want to create video that represents the simulation
# So we need to sample only a few of the frames:
for i in range(0, simulation_size, sample_rate):
# Update Pendulum Arm:
x_data_points = [origin[0], x_pendulum_arm[i]]
y_data_points = [origin[1], y_pendulum_arm[i]]
x2_data_points = [x_vec[i], x_pendulum_arm[i]]
y2_data_points = [y[i], y_pendulum_arm[i]]
# We want to avoid creating new plots to make an animation (Very Slow)
# Instead lets take the plot we made earlier and just update it with new data.
p1.set_data(x_data_points, y_data_points) # Update plot with set_data
p2.set_data(x2_data_points, y2_data_points)
# Update Pendulum Patch:
patch_center = x_pendulum_arm[i], y_pendulum_arm[i]
c.center = patch_center
m_center = x_vec[i], y[i]
b.center = m_center
# Update Drawing:
fig.canvas.draw() # Update the figure with the new changes
# Grab and Save Frame:
writerObj.grab_frame()
# We can also run it inline with the following block of code:
from IPython.display import Video
Video("/work/simulation 2C.mp4", embed=True, width=640, height=480)
# You can find the video's path by right clicking on it in Notebooks & Files
# and selecting "Copy path to clipboard"