Lab 5 - FouadM0426/comp.physics-spring-2025 GitHub Wiki
Goal
-
The goal of this lab is the following to compare different types of calculations Euler,Euler-cromer,and analytic. We will be using a mass on a spring, for our testing for our methods of calculating accurate data.
-
Part 1 First, we will be using a simple graph to visualize the graph, this will show a standard sin wave, in which we will be able to better picture the amplitude and period of the graph. Next we will calculate the analytic solution for spring on the mass. The code below will show all of part 1 steps complete
import numpy as np
from matplotlib import pyplot as plt
!pip install wget
import wget
# will import data
# use astropy to read in your data file
from astropy.io import ascii
data = ascii.read("shm_mass_spring_good.csv",delimiter=',')
plt.plot(time, position)
plt.xlabel("time (s)")
plt.ylabel("position (m)")
plt.title("position vs time")
plt.show()
# amplitude and period
A = .0431
T = .75
omega = 2*np.pi/T
print(omega)
def calc_SHO(A,omega,phi,t):
return A*np.cos(omega*t + phi)
t= np.linspace(0,5,1000)
plt.plot(time, position,"bo",label="my data")
plt.xlabel("time (s)")
plt.ylabel("position (m)")
plt.plot(t, calc_SHO(A,omega,.8,t),label="analytic Data Points")
plt.title("position vs time")
plt.legend()
plt.show()
Part 2
The goal of part 2 we will be modeling motion,with Euler method and compare it to our analytic solution for the energy functions The following code below shows all steps to create the Euler for the energy calculated.
def acceleration(k,m,x):
"""
description: will take the following paramers and will return acc. using hooks law
parameters:k,m,x
return:accerlation of the mass, using hooke law
"""
acc = -k*x/m
return acc
def initialize(x0,v0,tmax,dt):
"""
description: a function to create and int position and velocity arrays
parameters: x0,v0,tmax,dt
return: time,position,velocity
"""
nsteps = int(tmax/dt)
time = np.zeros(nsteps)
position = np.zeros(nsteps)
velocity = np.zeros(nsteps)
position[0] = x0
velocity[0] = v0
return time,position,velocity
def calculate(time,position,velocity,k,m,dt):
"""
description: will calcuate the position and velocity using the euler method
parameters: time,position,velocity,k,m,dt
return: time,position,velocity
"""
# your code here
for i in range(len(time)-1):
time[i+1] = time[i] + dt
position[i+1] = position[i] + velocity[i]*dt
velocity[i+1] = velocity[i] + acceleration(k,m,position[i])*dt
return time,position,velocity
x0 = 0.03
v0 = -0.32
# set the amplitude, mass and spring constant to the values from Part 1
A = .0431
m = .150
k = w**2*m
print(k)
# Enter the period you determined from Part 1
T = .75
# set the time step to a small fraction of the period
dt = T/100
# set tmax to 5x the period
tmax = 5 * T
# call initialize
time,position,velocity = initialize(x0,v0,tmax,dt)
# call calculate
time,position,velocity = calculate(time,position,velocity,k,m,dt)
# plot position versus time
t= np.linspace(0,5,1000)
plt.figure(figsize=(8,6))
plt.plot(time,position,label="my data")
plt.xlabel("time (s)")
plt.ylabel("position (m)")
plt.title("position vs time")
# set the figure size to (8,6)
plt.plot(t, calc_SHO(A,omega,.8,t),'b',label="analytic")
plt.legend()
plt.show()
plt.figure(figsize=(8,6))
plt.plot(time,0.5*m*velocity**2 + 0.5*k*position**2)
plt.xlabel("time (s)")
plt.ylabel("energy (J)")
plt.show()
Part 3
In this part of the lab, we will implement the Euler-cromer method and compare to the analytic
def euler_cromer_calculate(time,position,velocity,k,m,dt):
"""
description:The function will calcualate the euler cromer method
parameters: time,position,velocity,k,m,dt
return: time,position,velocity
"""
# your code here
for i in range(1,len(time)):
velocity[i] = velocity[i-1] + acceleration(k,m,position[i-1])*dt
position[i] = position[i-1] + velocity[i]*dt
time[i] = time[i-1] + dt
return time,position,velocity
x0 = 0.03
v0 = -0.32
# set the amplitude, mass and spring constant to the values from Part 1
A = .0431
m = .150
k = w**2*m
print(k)
# Enter the period you determined from Part 1
T = .75
# set the time step to a small fraction of the period
dt = T/100
# set tmax to 5x the period
tmax = 5 * T
# call initialize
time,position,velocity = initialize(x0,v0,tmax,dt)
# call calculate
time,position,velocity = euler_cromer_calculate(time,position,velocity,k,m,dt)
# plot position versus time
t= np.linspace(0,5,1000)
plt.figure(figsize=(8,6))
plt.plot(time,position,label="std euler method")
plt.xlabel("time (s)")
plt.ylabel("position (m)")
# set the figure size to (8,6)
plt.plot(t, calc_SHO(A,omega,.8,t),'b',label="cromer")
plt.legend()
# overplot your analytic solution from Part 1
plt.plot(time,position)
plt.title("position vs time")
plt.show()
plt.figure(figsize=(8,6))
plt.plot(time,0.5*m*velocity**2 + 0.5*k*position**2,label="calc")
plt.xlabel("time (s)")
plt.ylabel("energy (J)")
plt.plot(time,position,label="energy act")
plt.title("energy vs time")
# but plot the output from the Euler-Cromer solution here
plt.scatter(position, velocity, c=time,label="plot pos.")
plt.legend()
plt.show()
""" The circular graph which is the velocity vs pos. graph will show us the acc. at any given position of the SHM
"""