#HW1 Problem 1(b): Math and Plots # 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 from numpy.linalg import solve import matplotlib.pyplot as plt from matplotlib.animation import FFMpegWriter from matplotlib.patches import Circle # Parameters g = 9.8 m = 1 # Define Time Vector for Simulation dt = 0.001 # Time Step sim_time = 10 # Length of Simulation t_vec = np.arange(0, 10, dt) # Initialize State Vectors: vec_size = len(t_vec) # We use this value alot so lets store it x_vec = np.zeros(vec_size) dx_vec = np.zeros(vec_size) y_vec = np.zeros(vec_size) dy_vec = np.zeros(vec_size) lamb_vec = np.zeros(vec_size) # Set our initial condition x_vec[0] = -1 # initial x position dx_vec[0] = 0 # initial x velocity y_vec[0] = 1 # initial y position dy_vec[0] = 0 # initial y velocity for i in range(1, vec_size): A = np.array([[m, 0, -2*x_vec[i-1]], [0, m, 1], [-2*x_vec[i-1], 1, 0]]) B = np.array([0, -m*g, 2*dx_vec[i-1]**2]) [ddx, ddy, lamb] = np.linalg.solve(A, B) 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 lamb_vec[i] = lamb # plots plt.plot(t_vec,x_vec, label = "x(t)") plt.plot(t_vec, y_vec, label = "y(t)") plt.plot(t_vec, lamb_vec, label = "lambda") plt.axis([0, 10, -5, 5]) plt.legend() plt.show()

#HW1 Problem 1(c): Animation # Create Animation: # Setup Figure: fig, ax = plt.subplots() p, = ax.plot([], [], color='royalblue') min_lim = -10 max_lim = 10 ax.axis('equal') ax.set_xlim([min_lim, max_lim]) ax.set_ylim([min_lim, max_lim]) ax.set_xlabel('X') ax.set_ylabel('Y') ax.set_title('Bead on a Wire Simulation:') title = "HW1pr1animation" # Setup Animation Writer: FPS = 20 sample_rate = int(1 / (dt * FPS)) # Real Time Playback dpi = 300 writerObj = FFMpegWriter(fps=FPS) # Initialize Patch: xy_bead=(x_vec[0], y_vec[0]) c = Circle(xy_bead, radius=0.1, color='cornflowerblue') # Draws a circle of radius 0.1 at (X, Y) #location (0, 0) ax.add_patch(c) # Add the patch to the axes # We need to calculate the cartesian location of the bead: # Initialize Array: simulation_size = len(t_vec) # We use this alot, so lets store it. x_bead = np.zeros(simulation_size) y_bead = np.zeros(simulation_size) x_wire = np.linspace(-2, 2, num = vec_size) y_wire = np.square(x_wire) with writerObj.saving(fig, title + ".mp4", dpi): plt.plot(x_wire, y_wire) plt.show() for i in range(0, vec_size, sample_rate): # Update bead x_bead[i] = x_wire[i] y_bead[i] = x_wire[i]**2 # Update bead patch: patch_center = x_bead[i], y_bead[i] c.center = patch_center # Update Drawing: fig.canvas.draw() # Save Frame: writerObj.grab_frame()

#HW1 Problem 2(b): Math and Plots # 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 from numpy.linalg import solve import matplotlib.pyplot as plt from matplotlib.animation import FFMpegWriter from matplotlib.patches import Circle from matplotlib.patches import Rectangle # Parameters g = 9.8 m1 = 1. m2 = 2. L = 1. h = 1.5 k= 10. lo=1. # Define Time Vector for Simulation dt = 0.001 # Time Step sim_time = 10 # Length of Simulation t_vec = np.arange(0, 10, dt) #I know that it is supposed to be a simulation of 10 secs, but I was unable to get my #animation to run for that amount due to issues I could not determine with my plot. #Animation will run for a simulation time of 3 secs # Initialize State Vectors: vec_size = len(t_vec) # We use this value alot so lets store it x_vec = np.zeros(vec_size) dx_vec = np.zeros(vec_size) theta_vec = np.zeros(vec_size) dtheta_vec = np.zeros(vec_size) # Set our initial condition x_vec[0] = 0. # initial x position dx_vec[0] = 0. # initial x velocity theta_vec[0] = 1. # initial angle dtheta_vec[0] = 0. for i in range(1, vec_size): A = np.array([[m1*L**2, 0.], [0., m2]]) B = np.array([[-L*m1*g*np.sin(theta_vec[i-1]) - 0.5*k*L*(2*h*np.sin(theta_vec[i-1]) - 2*x_vec[i-1]*np.cos(theta_vec[i-1]))*(1 - lo/np.sqrt(h**2 + L*(-2*h*np.cos(theta_vec[i-1]) - 2*x_vec[i-1]*np.sin(theta_vec[i-1])) + L**2 + x_vec[i-1]**2))], [-0.5*k*(-2*L*np.sin(theta_vec[i-1]) + 2*x_vec[i-1] - (-2*L*np.sin(theta_vec[i-1]) +2*x_vec[i-1])*lo/np.sqrt(h**2 + L*(-2*h*np.cos(theta_vec[i-1]) - 2*x_vec[i-1]*np.sin(theta_vec[i-1])) + L**2 + x_vec[i-1]**2))]]) [ddtheta, ddx] = np.linalg.solve(A, B) x_vec[i] = x_vec[i - 1] + dx_vec[i - 1] * dt theta_vec[i] = theta_vec[i - 1] + dtheta_vec[i - 1] * dt dx_vec[i] = dx_vec[i - 1] + ddx * dt dtheta_vec[i] = dtheta_vec[i - 1] + ddtheta * dt # plots plt.plot(t_vec,x_vec, label = "x(t)") plt.plot(t_vec, theta_vec, label = "theta(t)") plt.axis([0, 10, -10, 10]) plt.legend() plt.show()

#HW1 Problem 2: Animation t_vec = np.arange(0, 3, dt) # Setup Figure: Initialize Figure / Axe Handles fig, ax = plt.subplots() # Initialize Figure and Axes p, = ax.plot([], [], color='cornflowerblue') # Initialize Empty Plot ax.axis('equal') ax.set_xlim([-3, 3]) # X Lim ax.set_ylim([-3, 3]) # Y Lim ax.set_xlabel('X') # X Label ax.set_ylabel('Y') # Y Label ax.set_title('Pendulum and Spring Simulation:') # Plot Title video_title = "HW1pr2animation" # name of the mp4 y_offset = 0 #where spring and box are connected # Initialize Patch: c = Circle((0, 0), radius=0.1, color='cornflowerblue') # Draws a circle of radius 0.1 at #(X, Y) location (0, 0) ax.add_patch(c) # Add the patch to the axes # Initialize Patch: (box) width = 1 # Width of box height = width / 2 # Height of Cart xy_cart = (x_vec[0] - width / 2, y_offset - height / 2) # Bottom Left Corner of Cart r = Rectangle(xy_cart, width, height, color='cornflowerblue') # Rectangle Patch ax.add_patch(r) # Add Patch to Plot # Setup Animation Writer: FPS = 20 # Frames per second for the video sample_rate = int(1 / (dt * FPS)) # Rate at which we sample the simulation to visualize 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, h) # 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]) #x_cart = np.linspace(-2,2, num = simulation_size) y_cart = np.zeros(simulation_size) # Plot and Create Animation: with writerObj.saving(fig, video_title+".mp4", dpi): for i in range(0, simulation_size, sample_rate): # Update Pendulum Arm: x_data_points = [origin[0], x_pendulum_arm[i], x_vec[i]] y_data_points = [origin[1], y_pendulum_arm[i], y_cart[i]] p.set_data(x_data_points, y_data_points) # Update plot with set_data # Update Pendulum Patch: patch_center_circle = x_pendulum_arm[i], y_pendulum_arm[i] c.center = patch_center_circle # Same idea here, instead of drawing a new patch update its location # Update Rectangle Patch: r.set(xy=(x_vec[i] - width / 2, y_offset - height / 2)) # Update Drawing: fig.canvas.draw() # Update the figure with the new changes # Grab and Save Frame: writerObj.grab_frame()

#HW1 Problem 3(a): Math and Plots # 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 from numpy.linalg import solve import matplotlib.pyplot as plt from matplotlib.animation import FFMpegWriter from matplotlib.patches import Circle from mpl_toolkits import mplot3d from mpl_toolkits.mplot3d import Axes3D %matplotlib inline # Parameters g = 9.8 m = 1 # Define Time Vector for Simulation dt = 0.001 # Time Step sim_time = 10 # Length of Simulation t_vec = np.arange(0, 10, dt) # Initialize State Vectors: vec_size = len(t_vec) # We use this value alot so lets store it 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) lamb_vec = np.zeros(vec_size) # Set our initial condition x_vec[0] = 1 # initial x position dx_vec[0] = 0 # initial x velocity y_vec[0] = 0 # initial y position dy_vec[0] = 1 # initial y velocity z_vec[0] = x_vec[0]**2 + y_vec[0]**2 #initial z position dz_vec[0] = 2*(x_vec[0])*(dx_vec[0]) + 2*(y_vec[0])*(dy_vec[0]) #initial z velocity for i in range(1, vec_size): A = np.array([[m, 0, 0, 2*x_vec[i-1]], [0, m, 0, 2*y_vec[i-1]], [0, 0, m, -1], [2*x_vec[i-1], 2*y_vec[i-1], -1, 0]]) B = np.array([0, 0, -m*g, -2*np.square(dx_vec[i-1]) - 2*np.square(dy_vec[i-1])]) [ddx, ddy, ddz, lamb] = np.linalg.solve(A, B) x_vec[i] = x_vec[i - 1] + dx_vec[i - 1] * dt y_vec[i] = y_vec[i - 1] + dy_vec[i - 1] * dt z_vec[i] = z_vec[i - 1] + dz_vec[i - 1] * dt dx_vec[i] = dx_vec[i - 1] + ddx * dt dy_vec[i] = dy_vec[i - 1] + ddy * dt dz_vec[i] = dz_vec[i - 1] + ddz * dt lamb_vec[i] = lamb # plots plt.plot(t_vec,x_vec, label = "x(t)") plt.plot(t_vec, y_vec, label = "y(t)") plt.plot(t_vec, lamb_vec, label = "lambda") #plt.axis([0, 10, -5, 5]) plt.legend() plt.show()

#HW1 Problem 3(b): Animation # Setup Figure: Initialize Figure / Axe Handles fig = plt.figure() # Makes the figure ax = fig.add_subplot(111,projection='3d') # Axes are 3d now #xyz, = ax.plot([], [], [], color='cornflowerblue') # Initialize Empty Plot ax.set_xlim([-3, 3]) ax.set_ylim([-3, 3]) ax.set_zlim([-3, 3]) ax.set_xlabel('x(t)') ax.set_ylabel('y(t)') ax.set_zlabel('z(t)') ax.set_title('Bead 3D Animation') video_title = "HW1pr3animationAttempt2" # Initialize Patch: #xyz_bead=(x_vec[0], y_vec[0], z_vec[0]) #c = Circle(xyz_bead, radius=0, color='cornflowerblue') # Draws a circle of radius 0 #at (X, Y, Z) #ax.add_patch(c) # Add the patch to the axes # Setup Animation Writer: FPS = 20 # Frames per second for the video sample_rate = int(1 / (dt * FPS)) # Rate at which we sample the simulation to visualize dpi = 300 # Quality of the Video writerObj = FFMpegWriter(fps=FPS) # Video Object ax = plt.axes(projection='3d') # Data for a three-dimensional line simulation_size = len(t_vec) # We use this alot, so lets store it. #x_bead = np.zeros(simulation_size) #y_bead = np.zeros(simulation_size) #z_bead = np.zeros(simulation_size) ax.plot(x_vec,y_vec,z_vec) with writerObj.saving(fig, video_title+".mp4", dpi): for i in range(0, vec_size, sample_rate): xyz = ax.scatter3D(x_vec[i],y_vec[i],z_vec[i]) # Update drawing fig.canvas.draw() # Grab and save the frame writerObj.grab_frame() xyz.remove() #with writerObj.saving(fig, title + ".mp4", dpi): #xline = np.linspace(-15, 15, num = vec_size) #yline = np.linspace(-15, 15, num = vec_size) #zline = np.square(xline) + np.square(yline) #ax.plot3D(xline, yline, zline, 'gray') ## tried to use the following code for a mesh grid so it would be 3D but gave #me an error saying I did not have enough memory #x = np.linspace(-6, 6, vec_size) #y = np.linspace(-6, 6, vec_size) #X, Y = np.meshgrid(x, y) #Z = np.square(X) + np.square(Y) #ax.plot_wireframe(x, y, Z, color='black') #for i in range(0, vec_size, sample_rate): # Update bead #x_bead[i] = xline[i] #y_bead[i] = yline[i] #z_bead[i] = np.square(xline[i]) + np.square(yline[i]) ## again, the following was for meshgrid #x_bead[i] = x[i] #y_bead[i] = y[i] #z_bead[i]= np.square(x[i]) + np.square(y[i]) # Update bead patch: #patch_center = x_bead[i], y_bead[i], z_bead[i] #c.center = patch_center # Update Drawing: #fig.canvas.draw() # Save Frame: #writerObj.grab_frame() #from IPython.display import Video #Video("/work/HW1pr3animationAttempt2.mp4", embed =True, width=640, height=480)