Peeter Joot's Blog.

Math, physics, perl, and programming obscurity.

Posts Tagged ‘kepler problem’

Runge-Lenz vector conservation

Posted by peeterjoot on February 13, 2012

[Click here for a PDF of this post with nicer formatting and figures if the post had any (especially if my latex to wordpress script has left FORMULA DOES NOT PARSE errors.)]

Motivation.

Notes from Prof. Poppitz’s phy354 classical mechanics lecture on the Runge-Lenz vector, a less well known conserved quantity for the 3D 1/r potentials that can be used to solve the Kepler problem.

Motivation: The Kepler problem.

We can plug away at the Lagrangian in cylindrical coordinates and find eventually

\begin{aligned}\int_{\phi_0}^\phi d\phi = \int_{r_0}^r \frac{M}{m r^2} \frac{dr}{\sqrt{\frac{2}{M} ( E - U + \frac{M^2}{2 m r^2}) }}\end{aligned} \hspace{\stretch{1}}(2.1)

but this can be messy to solve, where we get elliptic integrals or worse, depending on the potential.

For the special case of the 3D problem where the potential has a 1/r form, this is what Prof. Poppitz called “super-integrable”. With 2N - 1 = 5 conserved quantities to be found, we’ve got one more. Here the form of that last conserved quantity is given, called the Runge-Lenz vector, and we verify that it is conserved.

Runge-Lenz vector

Given a potential

\begin{aligned}U = -\frac{\alpha}{r}\end{aligned} \hspace{\stretch{1}}(3.2)

and a Lagrangian

\begin{aligned}\mathcal{L} &= \frac{m \dot{r}^2}{2} + \frac{1}{{2}} \frac{M_z^2}{m r^2} - U \\ M_z &= m r^2 \dot{\phi}^2\end{aligned} \hspace{\stretch{1}}(3.3)

and writing the angular momentum as

\begin{aligned}\mathbf{M} = m \mathbf{r} \times \mathbf{v} \end{aligned} \hspace{\stretch{1}}(3.5)

the Runge-Lenz vector

\begin{aligned}\mathbf{A} = \mathbf{v} \times \mathbf{M} - \alpha \hat{\mathbf{r}},\end{aligned} \hspace{\stretch{1}}(3.6)

is a conserved quantity.

Verify the conservation assumption.

Let’s show that the conservation assumption is correct

\begin{aligned}\frac{d}{dt} \left( \mathbf{v} \times \mathbf{M} \right)=\frac{d{{ \mathbf{v}}}}{dt} \times \mathbf{M} + \mathbf{v} \times \not{{\frac{d{{\mathbf{M} }}}{dt}}}\end{aligned} \hspace{\stretch{1}}(3.7)

Here, we note that angular momentum conservation is really d\mathbf{M}/dt = 0, so we are left with only the acceleration term, which we can rewrite in terms of the Euler-Lagrange equation

\begin{aligned}\frac{d}{dt} \left( \mathbf{v} \times \mathbf{M} \right)&=-\frac{1}{{m}} \boldsymbol{\nabla} U \times M \\ &=-\frac{1}{{m}} \frac{\partial {U}}{\partial {r}} \hat{\mathbf{r}} \times M \\ &=-\frac{1}{{m}} \frac{\partial {U}}{\partial {r}} \hat{\mathbf{r}} \times (m \mathbf{r} \times \mathbf{v}) \\ &=- \frac{\partial {U}}{\partial {r}} \hat{\mathbf{r}} \times (\mathbf{r} \times \mathbf{v}) \end{aligned}

We can compute the double cross product

\begin{aligned}(\mathbf{a} \times (\mathbf{b} \times \mathbf{c}) )_i&=a_m b_r c_s \epsilon_{r s t} \epsilon_{m t i} \\ &=a_m b_r c_s \delta^{[rs]}_{i m} \\ &=a_m b_i c_m -a_m b_m c_i \end{aligned}

For

\begin{aligned}\mathbf{a} \times (\mathbf{b} \times \mathbf{c}) = (\mathbf{a} \cdot \mathbf{c}) \mathbf{b} -(\mathbf{a} \cdot \mathbf{b}) \mathbf{c}\end{aligned} \hspace{\stretch{1}}(3.8)

Plugging this we have

\begin{aligned}\frac{d}{dt} \left( \mathbf{v} \times \mathbf{M} \right) &= \frac{\partial U}{\partial r} \\ left( (\hat{\mathbf{r}} \cdot \mathbf{r}) \mathbf{v}-(\hat{\mathbf{r}} \cdot \mathbf{v}) \mathbf{r} \right) \\ &= \left( \frac{\alpha}{r^2} \right)\left( r \mathbf{v}-\frac{1}{r}(\mathbf{r} \cdot \mathbf{v}) \mathbf{r} \right) \\ &= \alpha\left( \frac{\mathbf{v}}{r}-\frac{(\mathbf{r} \cdot \mathbf{v}) \mathbf{r} }{r^3}\right) \\ \end{aligned}

Now let’s look at the other term. We’ll need the derivative of \hat{\mathbf{r}}

\begin{aligned}\frac{d{{\hat{\mathbf{r}}}}}{dt} &=\frac{d}{dt} \frac{\mathbf{r}}{r} \\ &=\frac{\mathbf{v}}{r} + \mathbf{r} \frac{d{{\frac{1}{{r}}}}}{dt} \\ &=\frac{\mathbf{v}}{r} - \frac{\mathbf{r}}{r^2} \frac{d{{r}}}{dt} \\ &=\frac{\mathbf{v}}{r} - \frac{\mathbf{r}}{r^2} \frac{d{{ \sqrt{\mathbf{r} \cdot \mathbf{r}}}}}{dt} \\ &=\frac{\mathbf{v}}{r} - \frac{\mathbf{r}}{r^2} \frac{\mathbf{v} \cdot \mathbf{r}}{\sqrt{\mathbf{r}^2}} \\ &=\frac{\mathbf{v}}{r} - \frac{\mathbf{r}}{r^3} \mathbf{v} \cdot \mathbf{r}\end{aligned}

Putting all the bits together we’ve now verified the conservation statement

\begin{aligned}\frac{d}{dt} \left(\mathbf{v} \times \mathbf{M} - \alpha \hat{\mathbf{r}}\right)=\alpha\left( \frac{\mathbf{v}}{r}-\frac{(\mathbf{r} \cdot \mathbf{v}) \mathbf{r} }{r^3}\right) -\alpha \left( \frac{\mathbf{v}}{r} - \frac{\mathbf{r}}{r^3} \mathbf{v} \cdot \mathbf{r} \right)= 0.\end{aligned} \hspace{\stretch{1}}(3.9)

With

\begin{aligned}\frac{d}{dt} \left( \mathbf{v} \times \mathbf{M} - \alpha \hat{\mathbf{r}} \right) = 0,\end{aligned} \hspace{\stretch{1}}(3.10)

our vector must be some constant vector. Let’s write this

\begin{aligned}\mathbf{v} \times \mathbf{M} - \alpha \hat{\mathbf{r}} = \alpha \mathbf{e},\end{aligned} \hspace{\stretch{1}}(3.11)

so that

\begin{aligned}\boxed{\mathbf{v} \times \mathbf{M} = \alpha \left(\mathbf{e} + \hat{\mathbf{r}} \right).}\end{aligned} \hspace{\stretch{1}}(3.12)

Dotting 3.12 with \mathbf{M} we find

\begin{aligned}\alpha \mathbf{M} \cdot \left(\mathbf{e} + \hat{\mathbf{r}} \right)&=\mathbf{M} \cdot (\mathbf{v} \times \mathbf{M}) \\ &= 0\end{aligned}

With \hat{\mathbf{r}} lying in the plane of the trajectory (perpendicular to \mathbf{M}), we must also have \mathbf{e} lying in the plane of the trajectory.

Now we can dot 3.12 with \mathbf{r} to find

\begin{aligned}\mathbf{r} \cdot (\mathbf{v} \times \mathbf{M}) &= \alpha \mathbf{r} \cdot \left(\mathbf{e} + \hat{\mathbf{r}} \right) \\ &= \alpha \cdot \left( r e \cos(\phi - \phi_0) + r \right) \\ \mathbf{M} \cdot (\mathbf{r} \times \mathbf{v}) &= \\ \mathbf{M} \cdot \frac{\mathbf{M}}{m} &= \\ \frac{\mathbf{M}^2}{m} &=\end{aligned}

This is

\begin{aligned}\frac{\mathbf{M}^2}{m} = \alpha r \left( 1 + e \cos(\phi - \phi_0) \right).\end{aligned} \hspace{\stretch{1}}(3.13)

This is a kind of curious implicit relationship, since \phi is also a function of r. Recall that the kinetic portion of our Lagrangian was

\begin{aligned}\frac{1}{{2}} m (\dot{r}^2 + r^2 \dot{\phi}^2 )\end{aligned} \hspace{\stretch{1}}(3.14)

so that our angular momentum was

\begin{aligned}M_\phi = \frac{\partial }{\partial {\dot{\phi}}} \left( \frac{1}{{2}} m r^2 \dot{\phi}^2 \right) = m r^2 \dot{\phi},\end{aligned} \hspace{\stretch{1}}(3.15)

with no \phi dependence in the Lagrangian we have

\begin{aligned}\frac{d}{dt} (m r^2 \dot{\phi}) = 0,\end{aligned} \hspace{\stretch{1}}(3.16)

or

\begin{aligned}\mathbf{M} = m r^2 \dot{\phi} \hat{\mathbf{z}} = \text{constant}\end{aligned} \hspace{\stretch{1}}(3.17)

Our dynamics are now fully specified, even if this not completely explicit

\begin{aligned}\boxed{\begin{aligned}r &= \frac{M^2}{m \alpha} \frac{1}{{1 + e \cos(\phi - \phi_0)}} \\ \frac{d\phi}{dt} &= \frac{M}{ m r^2}.\end{aligned}}\end{aligned} \hspace{\stretch{1}}(3.18)

What we can do is rearrange and separate variables

\begin{aligned}\frac{1}{{r^2}} = \frac{m^2 \alpha^2}{M^4} (1 + e \cos(\phi - \phi_0))^2 = \frac{m}{M} \frac{d\phi}{dt},\end{aligned} \hspace{\stretch{1}}(3.19)

to find

\begin{aligned}t - t_0 = \frac{M^3}{m \alpha^3} \int_{\phi_0}^\phi d\phi \frac{1}{{(1 + e \cos(\phi - \phi_0))^2}}=\frac{M^3}{m \alpha^3} \int_0^{\phi - \phi_0} du\frac{1}{{(1 + e \cos u)^2}}\end{aligned} \hspace{\stretch{1}}(3.20)

Now, at least \phi = \phi(t) is specified implicitly.

We can also use the first of these to determine the magnitude of the radial velocity

\begin{aligned}\frac{dr}{dt} &=-\frac{M^2}{m \alpha} \frac{1}{{(1 + e \cos(\phi - \phi_0))^2}} (-e \sin(\phi - \phi_0)) \frac{d\phi}{dt} \\ &=\frac{e M^2}{m \alpha} \frac{1}{{(1 + e \cos(\phi - \phi_0))^2}} \sin(\phi - \phi_0) \frac{M}{m r^2} \\ &=\frac{e M^3}{m^2 \alpha r^2} \frac{1}{{(1 + e \cos(\phi - \phi_0))^2}} \sin(\phi - \phi_0) \\ &=\frac{e M^3}{m^2 \alpha r^2} \left( \frac{ m r \alpha }{M^2} \right)^2 \sin(\phi - \phi_0) \\ &=\frac{e }{M } \sin(\phi - \phi_0),\end{aligned}

with this, we can also find the energy

\begin{aligned}E &= \dot{r}( m \dot{r}) + \dot{\phi} ( m r^2 \dot{\phi}) - \left( \frac{1}{{2}} m \dot{r}^2 + \frac{1}{{2}} m r^2 \dot{\phi}^2 - U \right) \\ &= \frac{1}{{2}} m \dot{r}^2 + \frac{1}{{2}} m r^2 \dot{\phi}^2 + U \\ &= \frac{1}{{2}} m \dot{r}^2 + \frac{1}{{2}} m r^2 \dot{\phi}^2 - \frac{\alpha}{r} \\ &= \frac{1}{{2}} m \frac{e^2}{M^2} \sin^2(\phi - \phi_0) + \frac{1}{{2 m r^2 }} M^2 - \frac{\alpha}{r}.\end{aligned}

Or

\begin{aligned}E= \frac{m}{2 M^2} (\mathbf{e} \times \hat{\mathbf{r}})^2 + \frac{1}{{2 m r^2 }} M^2 - \frac{\alpha}{r}.\end{aligned} \hspace{\stretch{1}}(3.21)

Is this what was used in class to state the relation

\begin{aligned}e = \sqrt{1 + \frac{2 E M^2}{m \alpha^2}}.\end{aligned} \hspace{\stretch{1}}(3.22)

It’s not obvious exactly how that is obtained, but we can go back to 3.18 to eliminate the e^2 \sin^2 \Delta \phi term

\begin{aligned}E = \frac{1}{{2}} m \frac{1}{M^2} \left( e^2 - \left( \frac{M^2}{r m \alpha} - 1\right)^2 \right) + \frac{1}{{2 m r^2 }} M^2 - \frac{\alpha}{r}.\end{aligned} \hspace{\stretch{1}}(3.23)

Presumably this simplifies to the desired result (or there’s other errors made in that prevent that).

Posted in Math and Physics Learning. | Tagged: , , , , , , , | Leave a Comment »

 
Follow

Get every new post delivered to your Inbox.