From b149738ca2c965ffcb0a56846c9468bfc42558e4 Mon Sep 17 00:00:00 2001
From: Erik Strand <erik.strand@cba.mit.edu>
Date: Thu, 20 Feb 2020 12:23:59 -0500
Subject: [PATCH] Add something for 3.1 (d)

---
 _code/pset_02_math.py |  48 ++++++++++++++++
 _psets/02.md          | 125 ++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 173 insertions(+)
 create mode 100644 _code/pset_02_math.py

diff --git a/_code/pset_02_math.py b/_code/pset_02_math.py
new file mode 100644
index 0000000..eff2a23
--- /dev/null
+++ b/_code/pset_02_math.py
@@ -0,0 +1,48 @@
+from sympy import *
+from sympy.assumptions.assume import global_assumptions
+#from sympy.abc import x, y
+
+#expr = spy.sin(x) / spy.cos(x)
+#spy.pprint(expr)
+#expr = spy.simplify(expr)
+#spy.pprint(expr)
+
+alpha, omega_1, t, A, B = symbols('alpha omega_1 t A B')
+x = exp(-alpha * t) * (A * exp(I * omega_1 * t) + B * exp(-I * omega_1 * t))
+dx = simplify(diff(x, t))
+ddx = simplify(diff(dx, t))
+pprint(x)
+print("")
+pprint(dx)
+print("")
+pprint(ddx)
+
+m, k = symbols('m k')
+global_assumptions.add(Q.real(m))
+global_assumptions.add(Q.positive(m))
+global_assumptions.add(Q.real(k))
+global_assumptions.add(Q.positive(k))
+omega_0, gamma = symbols('omega_0 gamma')
+E = simplify(Rational(1, 2) * (m * dx**2 + k * x**2))
+E_0 = simplify(E.subs(t, 0))
+#E_0 = simplify(E.subs([(t, 0), (omega_1, sqrt(omega_0**2 - alpha**2))]))
+#E_0 = simplify(E_0.subs([(omega_0, sqrt(k / m)), (alpha, gamma / (2 * m))]))
+print("")
+pprint(E)
+print("")
+pprint(E_0)
+
+dE = -gamma * dx
+delta_E = simplify(integrate(dE, (t, 0, 2 * pi / omega_1)))
+print("")
+pprint(delta_E)
+
+Q = 2 * pi * E_0 / delta_E
+Q = Q.subs(omega_1, sqrt(omega_0**2 - alpha**2))
+Q = Q.subs([
+    (omega_0, sqrt(k / m)),
+    (alpha, gamma / (2 * m))
+])
+Q = simplify(Q)
+print("")
+pprint(Q)
diff --git a/_psets/02.md b/_psets/02.md
index 7dad27f..a800ec1 100644
--- a/_psets/02.md
+++ b/_psets/02.md
@@ -134,6 +134,17 @@ $$
 A = \frac{1}{-m \omega^2 + i \gamma \omega + k}
 $$
 
+Thus the general solution to the inhomogeneous problem (assuming distinct roots) is
+
+$$
+x(t) = e^{-\frac{\gamma}{2m} t} \left(
+    A e^{\frac{\sqrt{\gamma^2 - 4mk}}{2m} t} +
+    B e^{-\frac{\sqrt{\gamma^2 - 4mk}}{2m} t}
+\right) + \frac{1}{-m \omega^2 + i \gamma \omega + k} e^{i \omega t}
+$$
+
+Here are plots of the magnitude and phase of the particular solution.
+
 ![amplitude](../assets/img/02_amplitude.png)
 
 ![phase](../assets/img/02_phase.png)
@@ -149,12 +160,126 @@ to the energy lost per radian (one cycle is $$2\pi$$ radians). Show that these t
 equal, assuming that the damping is small. How long does it take the amplitude of a 100 Hz
 oscillator with a Q of 109 to decay by $$1/e$$?
 
+#### Undamped
+
+As we found previously, the general solution is
+
+$$
+x(t) = e^{-\frac{\gamma}{2m} t} \left( A e^{\frac{\sqrt{\gamma^2 - 4mk}}{2m} t} + B e^{-\frac{\sqrt{\gamma^2 - 4mk}}{2m} t} \right)
+$$
+
+To get rid of a symbol, let's define $$\omega_0 = \sqrt{k/m}$$ and $$\alpha = \gamma / (2 m)$$.
+Then
+
+$$
+\frac{\sqrt{\gamma^2 - 4mk}}{2m} = \sqrt{\alpha^2 - \omega_0^2}
+$$
+
+Since we're assuming low damping, $$\alpha^2 - \omega_0^2$$ is negative. So we'll use
+
+$$
+\sqrt{\alpha^2 - \omega_0^2} = i \sqrt{\omega_0^2 - \alpha^2}
+$$
+
+Finally, to be as lazy as possible, we'll define $$\omega_1^2 = \omega_0^2 - \alpha^2$$ (which
+is a positive real number). All together this compresses our general solution to
+
+$$
+\begin{aligned}
+x(t)
+&= e^{-\alpha t} \left( A e^{\sqrt{\alpha^2 - \omega_0^2} t} + B e^{-\sqrt{\alpha^2 - \omega_0^2} t} \right) \\
+&= e^{-\alpha t} \left( A e^{i \sqrt{\omega_0^2 - \alpha^2} t} + B e^{-i \sqrt{\omega_0^2 - \alpha^2} t} \right) \\
+&= e^{-\alpha t} \left( A e^{i \omega_1 t} + B e^{-i \omega_1 t} \right) \\
+\end{aligned}
+$$
+
+The energy in the system is
+
+$$
+E = \frac{1}{2} m \dot{x}^2 + \frac{1}{2} k x^2
+$$
+
+It's easy enough to see that
+
+$$
+x^2(t) = e^{-2 \alpha t} \left( A e^{i \omega_1 t} + B e^{-i \omega_1 t} \right)^2
+$$
+
+Applying the product rule to the general form as written above, and doing some tedious bookkeeping,
+we can compute the $$\dot{x}$$ terms we need for the energy of the system.
+
+$$
+\begin{aligned}
+\dot{x}(t) &=
+-\alpha e^{-\alpha t} \left( A e^{i \omega_1 t} + B e^{-i \omega_1 t} \right)
++ i \omega_1 e^{-\alpha t} \left( A e^{i \omega_1 t} - B e^{-i \omega_1 t} \right) \\
+\dot{x}(t)^2 &=
+e^{-2 \alpha t} \big[ \alpha^2 \left( A e^{i \omega_1 t} + B e^{-i \omega_1 t} \right)^2 \\
+&- i \alpha \omega_1 \left( A e^{i \omega_1 t} + B e^{-i \omega_1 t} \right) \left( A e^{i \omega_1 t} - B e^{-i \omega_1 t} \right) \\
+&- \omega_1^2 \left( A e^{i \omega_1 t} - B e^{-i \omega_1 t} \right)^2 \big]
+\end{aligned}
+$$
+
+So at time $$t = 0$$ the energy in the system is
+
+$$
+\begin{aligned}
+E &= \frac{1}{2} m \big[ \alpha^2 (A + B)^2 - i \alpha \omega_1 (A + B)(A - B) - \omega_1^2 (A - B)^2 \big] \\
+&+ \frac{1}{2} k (A + B)^2 \\
+&= \frac{1}{2} \big[ (m \alpha^2 + k) (A + B)^2 - i m \alpha \omega_1 (A + B)(A - B) - m \omega_1^2 (A - B)^2 \big]
+\end{aligned}
+$$
+
+To compute the change in the energy of the system, we'd like to integrate $$dE/dt$$. Note that
+
+$$
+\begin{aligned}
+E &= \frac{1}{2} m \dot{x}^2 + \frac{1}{2} k x^2 \\
+\frac{dE}{dt} &= m \dot{x} \ddot{x} + k x \dot{x} \\
+&= \dot{x} (m \ddot{x} + k x) \\
+&= -\gamma \dot{x}^2
+\end{aligned}
+$$
+
+The last line follows because the ODE governing the system says $$m \ddot{x} + k x = -\gamma
+\dot{x}$$.
+
+...and I got caught in a sympy rabbit hole. Think it'll be valuable for future problem sets but
+cost me too much time on this one.
+
+
 ### (e)
 
 {:.question}
 Now find the solution to equation (3.58) by using Laplace transforms. Take the initial condition as
 $$x(0) = \dot{x}(0) = 0$$.
 
+Using the table in the book, and keeping in mind the initial conditions, I found the following
+transforms.
+
+$$
+\begin{aligned}
+\mathcal{L}[m \ddot{x}(t)]
+&= m \left( s^2 X(s) - s x(0) - \dot{x}(0) \right) \\
+&= m s^2 X(s) \\
+\mathcal{L}[\gamma \dot{x}(t)]
+&= \gamma \left( s X(s) - x(0) \right) \\
+&= \gamma s X(s) \\
+\mathcal{L}[k x(t)]
+&= k X(s) \\
+\mathcal{L}[e^{i \omega t}] &= \frac{1}{s - i \omega}
+\end{aligned}
+$$
+
+So our transformed ODE is
+
+$$
+m s^2 X(s) + \gamma s X(s) + k X(x) = \frac{1}{s - i \omega} \\
+X(s) = \frac{1}{(s - i \omega) (m s^2 + \gamma s + k)}
+$$
+
+Now it's just a bunch of algebra and an inverse transform...
+
 ### (f)
 
 {:.question}
-- 
GitLab