################## Homework 1 Problem 1 ###########################
# 992759-MOD-HW1
!apt update -y
!apt install ffmpeg -y
import numpy as np
import matplotlib.pyplot as plt
import sympy
from matplotlib.animation import FFMpegWriter
from matplotlib.patches import Circle
g = 9.8 # m/s^2
L = 1 # m
m = 1 # kg
# define our time vector for simulation
dt = .001
t_vec = np.arange(0,10,dt)
# define vectors that will hold our state variables over time (initializing vectors)
x_vec = np.zeros(len(t_vec))
dx_vec = np.zeros(len(t_vec))
ddx_vec = np.zeros(len(t_vec))
y_vec = np.zeros(len(t_vec))
dy_vec = np.zeros(len(t_vec))
ddy_vec = np.zeros(len(t_vec))
lambda_vec = np.zeros(len(t_vec)) # constraint
# inital conditions
x_vec[0] = -1 # position
dx_vec[0] = 0 # velocity
y_vec[0] = 1
dy_vec[0] = 0
# loop for x
for i in range(1, len(t_vec)):
x_vec[i] = x_vec[i-1] + dx_vec[i-1]*dt
ddx_vec[i-1] = -(2*x_vec[i-1]*(2*(dx_vec[i-1]**2) + g))/(4*(x_vec[i-1]**2) + 1)
dx_vec[i] = dx_vec[i-1] + ddx_vec[i-1]*dt
# loop for y
for i in range(1, len(t_vec)):
y_vec[i] = y_vec[i-1] + dy_vec[i-1]*dt
ddy_vec[i-1] = -(2*(2*g*x_vec[i-1]**2 - dx_vec[i-1]**2))/(4*x_vec[i-1]**2 + 1)
dy_vec[i] = dy_vec[i-1] + ddy_vec[i-1]*dt
# loop for constraint
for i in range(1, len(t_vec)):
lambda_vec[i] = (2*m*dx_vec[i-1]**2 + g*m)/(4*dx_vec[i-1]**2 + 1)
############ Plotting #############
# x plot
plt.plot(t_vec,x_vec)
plt.title('X positon vs time')
plt.xlabel('time (s)')
plt.ylabel('x (m)')
plt.show()
# y plot
plt.plot(t_vec,y_vec)
plt.title('Y postion vs time')
plt.xlabel('time (s)')
plt.ylabel('y (m)')
plt.show()
# constraint plot
plt.plot(t_vec,lambda_vec)
plt.title('Constraint force vs time')
plt.xlabel('time (s)')
plt.show()
######### Animation #############
fig, ax = plt.subplots()
p, = ax.plot([], [], color='navy')
ax.axis('equal')
ax.set_xlim([-3, 3])
ax.set_ylim([-3, 3])
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_title('Problem 1 Simulation')
video_title = "Problem 1 Simulation"
c = Circle((0, 0), radius=0.2, color='cornflowerblue')
ax.add_patch(c)
# Setup Animation Writer:
FPS = 20
sample_rate = int(1 / (dt * FPS))
dpi = 300
writerObj = FFMpegWriter(fps=FPS)
# calculating pendulum location in cartesian coords
simulation_size = len(t_vec)
x_var = np.linspace(-2,2,2000)
y_var = x_var**2
plt.plot(x_var, y_var, color='black')
plt.show()
# Plot and Create Animation:
with writerObj.saving(fig, video_title+".mp4", dpi):
for i in range(0, simulation_size, sample_rate):
x_data = x_vec[i]
y_data = y_vec[i]
p.set_data(x_data, y_data) # updating plot
# patching
patch_center = x_vec[i], y_vec[i]
c.center = patch_center
fig.canvas.draw() # update drawing
writerObj.grab_frame()
from IPython.display import Video
Video("/work/Problem 1 Simulation.mp4", embed=True, width=640, height=480)