# Homework 1 Problem 3 # This code block installs ffmpeg to your Deepnote Machine # Allows Deepnote to write MP4 !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 from matplotlib.patches import Rectangle from matplotlib.pyplot import * from numpy import * from mpl_toolkits import mplot3d %matplotlib inline import numpy as np import matplotlib.pyplot as plt # Part b: ------------------------------------------------------------ # Simulate this system using Euler’s method (time step = 0.001 # seconds) in Python with initial conditions, x = -1m, dx/dt = 0m/s, # over the course of 10 seconds # Parameters g = 9.81 # m/s^2 m = 1 # kg # Time array dt = 0.001 t_vec = np.arange(0,10,dt) # Initialize a vector of zeros vec_size = len(t_vec) x_vec = np.zeros(vec_size) dx_vec = np.zeros(vec_size) y_vec = np.zeros(vec_size) dy_vec = np.zeros(vec_size) z_vec = np.zeros(vec_size) dz_vec = np.zeros(vec_size) lam_vec = np.zeros(vec_size) # Pole End effector Location for Animation: x_pole = np.zeros(vec_size) y_pole = np.zeros(vec_size) # Initial Conditions x_vec[0] = -1 y_vec[0] = 1 z_vec[0] = 1 lam_vec[0] = 0 # Initial Pole End effector Location: y_offset = 0 # The Y Location of where the Cart and Pole are connected x_pole[0] = x_vec[0] #+ L * np.sin(theta_vec[0]) y_pole[0] = y_vec[0] #y_offset #- L * np.cos(theta_vec[0]) # Euler Simulation: Using Matrix Form (A * x = B) # Initialize A and B: A = np.array([[m, 0, 0, 0], [0, m, 0, 0], [0, 0, m, 1], [0, 0, 1, 0]]) b = np.array([[0], [0], [-m*g], [0]]) for i in range (1, vec_size): # Replace non constant values A [0, 3] = -2 * x_vec[i-1] A [1, 3] = -2 * y_vec[i-1] A [3, 0] = -2 * x_vec [i-1] A [3, 1] = -2 * y_vec[i-1] b [3] = 2 * dx_vec[i-1]**2 + 2 * dy_vec[i-1]**2 # Solve for ddx, ddy, and lambda [ddx, ddy, ddz, lam] = np.linalg.solve(A, b) # Euler integration x_vec [i] = x_vec[i-1] + dx_vec[i-1] * dt dx_vec[i] = dx_vec[i-1] + ddx * dt y_vec [i] = y_vec[i-1] + dy_vec[i-1] * dt dy_vec[i] = dy_vec[i-1] + ddy * dt z_vec[i] = x_vec[i]**2 + y_vec[i]**2 #y_vec[i] = x_vec[i]**2 lam_vec[i] = lam # x position vs time plotting plt.plot(t_vec, x_vec) plt.title("x position vs time") plt.show() # y position vs time plotting plt.plot(t_vec, y_vec) plt.title("y position vs time") plt.show() # Part c: ------------------------------------------------------------ # Animate the motion from Problem 1b, showing a bead sliding on a # wire (at 20 frames per second). # Setup Figure: def f(x, y): return (x ** 2 + y ** 2) x = np.linspace(-6, 6, 30) y = np.linspace(-6, 6, 30) X, Y = np.meshgrid(x, y) Z = f(X, Y) fig = plt.figure() ax = plt.axes(projection='3d') ax.contour3D(X, Y, Z, 50, cmap='binary') ax.set_xlabel('x') ax.set_ylabel('y') ax.set_zlabel('z'); ax.set_title('Bead in a paraboloid:') title = "simulation" # Setup Animation Writer: FPS = 20 sample_rate = int(1 / (dt * FPS)) # Real Time Playback dpi = 300 writerObj = FFMpegWriter(fps=FPS) # Draw a point bead = ax.scatter([0], [0], [0], color="g", s=300) # Animate: with writerObj.saving(fig, title + ".mp4", dpi): for i in range(0, vec_size, sample_rate): bead = [x_vec[i], y_vec[i], z_vec[i]] #r.set_data(bead) # Update Bead Patch: # Plotting Circle bead #xy_bead = [x_vec[i], y_vec[i]] #c = Circle(xy_bead, radius = 0.1, color = 'cornflowerblue') #ax.add_patch(c) # Add Patch to Plot #p.set_data(xy_bead)#, y_bead) # Update Cart Patch: #r.set(xy =(x_vec[i]-0.2, y_vec[i]-0.1)) #c.set(xy =(x_vec[i], y_vec[i]), radius = 0.1) # Update Drawing: fig.canvas.draw() # Save Frame: writerObj.grab_frame() from IPython.display import Video Video("/work/simulation.mp4", embed=True, width=640, height=480)