## Classical Mechanics - Week 9

### Last Week:
- We saw how a potential can be used to analyze a system
- Gained experience with plotting and integrating in Python 

### This Week:
- We will study harmonic oscillations using packages
- Further develope our analysis skills 
- Gain more experience wtih sympy


```python
# Let's import packages, as usual
import numpy as np
import matplotlib.pyplot as plt
import sympy as sym
sym.init_printing(use_unicode=True)
```

Let's analyze a spring using sympy. It will have mass $m$, spring constant $k$, angular frequency $\omega_0$, initial position $x_0$, and initial velocity $v_0$.

The motion of this harmonic oscillator is described by the equation:

eq 1.) $m\ddot{x} = -kx$

This can be solved as 

eq 2.) $x(t) = A\cos(\omega_0 t - \delta)$, $\qquad \omega_0 = \sqrt{\dfrac{k}{m}}$

Use SymPy below to plot this function. Set $A=2$, $\omega_0 = \pi/2$ and $\delta = \pi/4$. 

(Refer back to ***Notebook 7*** if you need to review plotting with SymPy.)


```python
# Plot for equation 2 here
```

## Q1.) Calculate analytically the initial conditions, $x_0$ and $v_0$, and the period of the motion for the given constants.  Is your plot consistent with these values?

&#9989; Double click this cell, erase its content, and put your answer to the above question here.

#### Now let's make plots for underdamped, critically-damped, and overdamped harmonic oscillators.
Below are the general equations for these oscillators:

- Underdamped, $\beta < \omega_0$ : 

eq 3.)    $x(t) = A e^{-\beta t}cos(\omega ' t) + B e^{-\beta t}sin(\omega ' t)$ , $\omega ' = \sqrt{\omega_0^2 - \beta^2}$
 
 ___________________________________
 
 
- Critically-damped, $\beta = \omega_0$:

eq 4.)    $x(t) = Ae^{-\beta t} + B t e^{-\beta t}$

 ___________________________________
 
- Overdamped, $\beta > \omega_0$:

eq 5.)    $x(t) = Ae^{-\left(\beta + \sqrt{\beta^2 - \omega_0^2}\right)t} + Be^{-\left(\beta - \sqrt{\beta^2 - \omega_0^2}\right)t}$
    
 _______________________
 
In the cells below use SymPy to create the Position vs Time plots for these three oscillators. 

Use $\omega_0=\pi/2$ as before, and then choose an appropriate value of $\beta$ for the three different damped oscillator solutions. Play around with the variables, $A$, $B$, and $\beta$, to see how different values affect the motion and if this agrees with your intuition.


```python
# Put your code for graphing Underdamped here
```


```python
# Put your code for graphing Critical here
```


```python
# Put your code for graphing Overdamped here
```

## Q2.) How would you compare the 3 different oscillators?

&#9989; Double click this cell, erase its content, and put your answer to the above question here.

# Here's another simple harmonic system, the pendulum. 

The equation of motion for the pendulum is:

eq 6.) $ml\dfrac{d^2\theta}{dt^2} + mg \sin(\theta) = 0$, where $v=l\dfrac{d\theta}{dt}$ and $a=l\dfrac{d^2\theta}{dt^2}$

In the small angle approximation $\sin\theta\approx\theta$, so this can be written:

eq 7.) $\dfrac{d^2\theta}{dt^2} = -\dfrac{g}{l}\theta$

We then find the period of the pendulum to be $T = \dfrac{2\pi}{\sqrt{l/g}}$ and the angle at any given time 
(if released from rest) is given by 

$\theta = \theta_0\cos{\left(\sqrt{\dfrac{g}{l}} t\right)}$.

Let's use Euler's Forward method to solve equation (7) for the motion of the pendulum in the small angle approximation, and compare to the analytic solution.

First, let's graph the analytic solution for $\theta$. Go ahead and graph using either sympy, or the other method we have used, utilizing these variables:

- $t:(0s,50s)$
- $\theta_0 = 0.5$ radians
- $l = 40$ meters


```python
# Plot the analytic solution here
```

Now, use Euler's Forward method to obtain a plot of $\theta$ as a function of time $t$ (in the small angle approximation).  Use eq (7) to calculate $\ddot{\theta}$ at each time step.
Try varying the time step size to see how it affects the Euler's method solution.


```python
# Perform Euler's Method Here
```

You should have found that if you chose the time step size small enough, then the Euler's method solution was
indistinguishable from the analytic solution.  

We can now trivially modify this, to solve for the pendulum **exactly**, without using the small angle approximation.
The exact equation for the acceleration is

eq 8.) $\dfrac{d^2\theta}{dt^2} = -\dfrac{g}{l}\sin\theta$.

Modify your Euler's Forward method calculation to use eq (8) to calculate $\ddot{\theta}$ at each time step in the cell below.



```python

```

# Q3.) What time step size did you use to find agreement between Euler's method and the analytic solution (in the small angle approximation)? How did the exact solution differ from the small angle approximation? 

&#9989; Double click this cell, erase its content, and put your answer to the above question here.

### Now let's do something fun:

In class we found that the 2-dimensional anisotropic harmonic motion can be solved as

eq 8a.) $x(t) = A_x \cos(\omega_xt)$

eq 8b.) $y(t) = A_y \cos(\omega_yt - \delta)$

If $\dfrac{\omega_x}{\omega_y}$ is a rational number (*i.e,* a ratio of two integers), then the trajectory repeats itself after some amount of time. The plots of $x$ vs $y$ in this case are called Lissajous figures (after the French physicists Jules Lissajous).  If $\dfrac{\omega_x}{\omega_y}$ is not a rational number, then the trajectory does not repeat itself, but it still shows some very interesting behavior.

Let's make some x vs y plots below for the 2-d anisotropic oscillator. 

First, recreate the plots in Figure 5.9 of Taylor.  (Hint: Let $A_x=A_y$. For the left plot of Figure 5.9, let $\delta=\pi/4$ and for the right plot, let $\delta=0$.)

Next, try other rational values of $\dfrac{\omega_x}{\omega_y}$ such as 5/6, 19/15, etc, and using different phase angles $\delta$.

Finally, for non-rational $\dfrac{\omega_x}{\omega_y}$, what does the trajectory plot look like if you let the length of time to be arbitrarily long?

\[For these parametric plots, it is preferable to use our original plotting method, *i.e.* using `plt.plot()`, as introduced in ***Notebook 1***.\]


```python
# Plot the Lissajous curves here
```

# Q4.) What are some observations you make as you play with the variables? What happens for non-rational $\omega_x/\omega_y$ if you let the oscillator run for a long time?


&#9989; Double click this cell, erase its content, and put your answer to the above question here.

# Notebook Wrap-up. 
Run the cell below and copy-paste your answers into their corresponding cells.


```python
from IPython.display import HTML
HTML(
"""

"""
)
```

# Well that's that, another Notebook! It's now been 10 weeks of class

You've been given lots of computational and programing tools these past few months. These past two weeks have been practicing these tools and hopefully you are understanding how some of these pieces add up. Play around with the code and see how it affects our systems of equations. Solve the Schrodinger Equation for the Helium atom. Figure out the unifying theory. The future is limitless!

97e4ea8cb9864cd7b107c814f69f83269316fbc5
ipynb
Jupyter Notebook
doc/AdminBackground/PHY321/CM_Jupyter_Notebooks/Student_Work/CM_Notebook9.ipynb
Shield94/Physics321
9875a3bf840b0fa164b865a3cb13073aff9094ca
CC0-1.0
2020-01-09T17:41:16.000Z
2022-03-09T00:48:58.000Z
doc/AdminBackground/PHY321/CM_Jupyter_Notebooks/Student_Work/CM_Notebook9.ipynb
Shield94/Physics321
9875a3bf840b0fa164b865a3cb13073aff9094ca
CC0-1.0
2020-01-08T03:47:53.000Z
2020-12-15T15:02:57.000Z
doc/AdminBackground/PHY321/CM_Jupyter_Notebooks/Student_Work/CM_Notebook9.ipynb
Shield94/Physics321
9875a3bf840b0fa164b865a3cb13073aff9094ca
CC0-1.0
2020-01-10T20:40:55.000Z
2022-02-11T20:28:41.000Z
Qwen/Qwen-72B
1. YES
2. YES