Numerical Integration - Simple Harmonic Oscillator

In this problem, we want to find the velocity and position as a function of time of a simple harmonic oscillator. A common example of a harmonic oscillator is a weight (with a mass $m$) that has been attached to the end of a spring with some resistance to it. If you pull the weight away from the resting position, there is a restoring force that tries to pull the spring back to where it was (we pull, the spring pulls back).

Mathematically, if we let the resting position be $x = 0$, then the restoring force is equal to

$F = -kx$

where $k$ is known as the spring constant, which simply describes how powerful the spring is (i.e. how stiff it is). If $k$ is large, the string is more stiff, if $k$ is small, it's easier to pull.

However, we know, thanks to Newton's Second Law,

$F = ma = m\frac{\Delta v}{\Delta t}$

So, this implies,

$a = \frac{F}{m} = \frac{-k}{m} x$

and by definition, velocity (v) is related to acceleration and position via

$a = \frac{\Delta v}{\Delta t}$ and $v = \frac{\Delta v}{\Delta t}$

So, let's say at some sample time $t$, the weight will be moving at some velocity $v$ and will be at position $x$. What will be it's position and velocity be at some small time later ($t+\Delta t$)?

$v_{\rm new} = v + \Delta v = v + (a \Delta t)$

$x_{\rm new} = x + \Delta x = x + (v \Delta t)$

This is the same idea we used in lecture to create a numerical integrator for a ball being thrown through the air.

Part 1

Using the discussion above, for a oscillator with some mass $m = 1$ and a spring constant $k = 1$, can you calculate the position and velocity of the weight at some time later? Use a while loop to do this, and let your $\Delta t$ = 0.01 between each step.

Let your initial $x$ = 3, initial $v$ = 0, and initial $t$ = 0. What are $v$ and $x$ at $t = 5$.


In [ ]:
m = 1
k = 1
x = 3
v = 0
t = 0
tFinal = 5
deltaT = 0.01
while (t < tFinal):
    deltaV = ((-k/m)*x)*deltaT
    deltaX = v*deltaT
    v = v + deltaV
    x = x + deltaX
    t = t + deltaT
print(x)
print(v)

Part 2

Use arrays to store the values of x and t for the loop, and then plot $t$ v.s. $x$. Can you do the same thing for $t$ v.s. $v$?

Remember, you can create an empty array using,

storage_array = np.array([])

And then add to the end of the array using np.append()

storage_array = np.append(storage_array,[new_value])

For this, let the final time be $t$ = 50.


In [ ]:
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np

m = 1
k = 1
x = 3
v = 0
t = 0
tFinal = 50
deltaT = 0.01

xArray = np.array([])
tArray = np.array([])

while (t < tFinal):
    deltaV = ((-k/m)*x)*deltaT
    deltaX = v*deltaT
    v = v + deltaV
    x = x + deltaX
    t = t + deltaT
    xArray = np.append(xArray,[x])
    tArray = np.append(tArray,[t])
print(x)
print(v)

plt.plot(tArray,xArray)
plt.show()

In [ ]:


In [ ]: