• Introduction
• Scattering with Central Forces
• The Ellipse
• Introduction to Orbits
• Orbital Mechanics
• Lagrangian Formulation
• Real Two-body Dynamics
• Many body Dynamics

Introduction

When my kids were little, I'd tell them about gravity and how the force between two bodies was attractive. Then we'd look up at the moon and I'd ask "what keeps it up there"? Their conclusion was that gravity must not be what I said it was!

The motion of bodies under the influence of gravity are called orbits, and the dynamics of such motion is indeed interesting and complex. The lesson here is that as is much of physics, motion can sometimes be counter-intuitive. So, what is the core concept that allows us to understand how the moon stays up there? Best to start at the beginning, which means consider the forces at work, in this case gravity.

The beginning starts with Newton's Law of Universal Gravitation, which says the following:

1. Every object exerts an attractive force on every other object
2. The force is along the line connecting the two objects
3. The magnitude of the force is proportional to the product of the masses and inversely proportional to the square of the distance between them

In the figure below, there are 2 masses $m_1$ and $m_2$ separated by a distance $r$:

By Newton's laws, the gravitational force is given by: $$\vec{F} = \frac{Gm_1m_2}{r^2}\hat{r}\label{edm1m2}$$ Given that the force is along the direction between them, we can say that the motion will have an acceleration only along that direction, and not along any other orthogonal direction. According to Newton's 1st law of motion, a body in motion will stay in motion unless acted upon by a force. That means that if there is any motion in the direction orthogonal to the radial direction, then that motion will be characterized by a constant velocity along that direction. This is actually the origin of conservation of angular momentum (more on that below), which will always occur whenever the forces are "central", as is the case for gravity (and the Coulomb force as well).

For central $1/r^2$ forces such as gravity, we can solve for the two-body motion analytically. As we will see below, the motion depends on the total energy of the system, and the interplay between the kinetic and potential energy:

• Negative total energy: The system is bound, and the bodies will orbit each other in an ellipse (or a circle as a special case)
• Zero total energy: The orbit is parabolic, and at infinity the final velocity is zero. In other words, a satellite can get infinitely far away but never escapes
• Positive total energy: The system is not bound and the bodies diverge from each other along a hyperbolic path with a non-zero velocity

This chapter mostly discusses bound state orbits ($E\lt 0$), which are generally elliptical. So to begin, we first discuss the mathematics of ellipses, just to have it ready when we apply it to the physics.

---

The quantity $f$ is called the focus, and the ellipse is defined such that if you calculate the distance of any point on the ellipse to the two foci and add together, you get a constant value no matter what point you use. So the distance $2g$ is equal to the distance $c+d$, which is also equal to the distance from each focus to the point $[x,y]=[a,0]$ given by $2a$. So we have our first formula:

$g=a$

Since $a>b$, $x$ is called the major axis and $y$ the minor axis, and the quantity $a$ and $b$ are often called the "semi-major" and "semi-minor" axis lengths.

Since the two focal points are equidistant from the origin, then we can calculate the hypotenuse using $g^2 = f^2 + b^2$, and substituting for $g=a$ we have

$a^2 = b^2 + f^2$

The ratio between the focal point and the semi-major axis length, $f/a$ is called the "eccentricity".

$$e \equiv f/a$$

Combining these equations gives

$b^2 = a^2 - f^2 = a^2(1-(f/a)^2) = a^2(1-e^2)$

which gives the other useful equation

$$(b/a)^2 = 1-e^2\label{eeba}$$

A circle would be the special case of an ellipse with $b=a$, which implies $e=0$, so we can say that $e$, the eccentricity, measures the deviation from a circle ($e=0$ means a circle).

In the figure below, you can drag around the black circles at the endpoints of the $x$ and $y$ axes, and see how the eccentricity and foci change. (The program will not let you make the semi-major axis swap from $x$ to $y$ however.)

To see the equation for the ellipse in polar coordinates, we transform from $[x,y]$ to $[r,\theta]$ ($\theta$ measures the angle between the line from the origin to the point on the ellipse with the major axis as in figure 2) using: $$x = r\cos{\theta}\label{eqx}$$ $$y = r\sin{\theta}\label{eqy}$$ And plug into equation $\ref{eellxy}$ to get:

$r^2(b^2\cos^2\theta + a^2\sin^2\theta) = a^2b^2$, which after some algebra gives

where here the origin is between the two foci and $\theta$ is the angle between any origin-ellipse point and the major axis as in figure 1 above.

This equation makes is clear that the condition $e=0$ means $r=b$, a constant, and thus a circle.

Figure 2. Ellipse centered at origin.

In discussing orbital mechanics, it is useful to have the equation for an ellipse centered at one of the foci. To do this, we would only have to take equation \ref{eellxy} and translate by $x\to x-f$. After some algebra, we would get the equation: $$r = \frac{a(1-e^2)}{1\pm e\cos\theta}\label{eellrtf}$$ where now $\theta$ is the angle between the major axis and a point on the ellipse relative to one of the foci as in figure 3 below. The $\pm$ sign mean tells you which foci is at the origin: "$+$" means that $r$ is minimal at $\theta=0$, so the origin coincides with the focus on the right, and "$-$" sign in the denominator means that the origin coincides with the focus on the left. The figure below would correspond to having $1-e\cos\theta$ in the denominator for equation $\ref{eellrtf}$, where the origin coincides with the focus on the left.

Figure 3. Ellipse centered at one of the focal points.

$$v^2 = \dot{x}^2 + \dot{y}^2 = \dot{r}^2 + r^2\dot{\theta}^2$$ Now we can form the Lagrangian and crank through the equations of motion. What we are after is an understanding of the dynamics, which starts with an understanding of what the various elliptical parameters $e,a,b,f$ depend on. For the full treatment, see below.

For a quicker way, we first start with an ellipse centered at the left foci as in figure 3 and use equation $\ref{eellrtf}$ (reproduced below) for the motion as a function of $\theta$ and the eccentricity $e$: $$r = \frac{a(1-e^2)}{1 + e\cos\theta}$$ We can differentiate to get

The important thing to notice here is that $\dot r$ only depends on $\theta$ ($\ell,a,e$ are all constant), which means $\dot r=0$ at $\theta=0,\pi$, which means the radial velocity vanishes at apicenter and pericenter.

Next, using the above equations $\ref{eqpotential}$ and $\ref{eqkinetic}$ and going to polar coordinates, we form the total energy $E$

At apicenter and pericenter, knowing that $\dot r_a=\dot r_p=0$, we have $$\varepsilon = \frac{\ell^2}{2r^2} - \frac{g}{r}\label{eapper}$$ Rearranging to get a quadratic in $r$ gives $$2\varepsilon r^2 + 2gr - \ell^2 = 0\label{erel}$$ Solving for $r$ gives $$r = \frac{-g\pm\sqrt{g^2 + 2\ell^2\varepsilon}}{2\varepsilon}$$ Remembering that $\varepsilon\lt0$ for bound orbits, we can write $\varepsilon = -|\varepsilon|$ to get $$r = \frac{g\pm\sqrt{g^2 - 2\ell^2|\varepsilon|}}{2|\varepsilon|}\label{ere}$$

The $\pm$ sign tells you the difference between $r_a$ and $r_b$: $$r_a = \frac{g + \sqrt{g^2 - 2\ell^2|\varepsilon|}}{2|\varepsilon|}$$ $$r_p = \frac{g - \sqrt{g^2 - 2\ell^2|\varepsilon|}}{2|\varepsilon|}$$

Going back to equation $\ref{ere}$, if you remember that $2a=r_a+r_p$ you can easily write down that $2a = -g/\varepsilon$, or $$\varepsilon = -\frac{g}{2a}$$ This is a great equation - it says that for central force gravitational forces, knowing the energy of the satellite will determine the length of the major axis. Does it also tell you all you need to know the eccentricity? To answer this, go back to equation $\ref{erel}$ and substitute for $r$ using either $r_a = a(1+e) = -g(1+e)/2\varepsilon$ or the similar equation for $r_p$, and crank through. After a few lines of algebraic simplification, you should get the equation $$e = \sqrt{1+\frac{2\varepsilon\ell^2}{g^2}}\label{eeeLE}$$ From equation $\ref{eellxy}$, we see that we need 2 parameters to describe the ellipse: $a$ and $b$. So it should be no surprise that we need to know 2 constants of the motion: $\varepsilon$ and $\ell$. (And of course we can use equation $\ref{eeba}$ to find the eccentricity $e$ knowing $a$ and $b$.)

To complete the equations of motion, we can take equation $\ref{eeeLE}$ and substitute $\varepsilon = -\frac{g}{2a}$ to get $$\ell^2 = ag(1-e^2)\label{eell}$$ Now that we have $a(1-e^2)=\ell^2/g$ and $e$ as in equation $\ref{eeeLE}$, we can use equation $\ref{eellrtf}$ to get the equation of motion for a satellite relative to the origin at the foci: $$r = \frac{\ell^2/g}{1- \sqrt{1+\frac{2\varepsilon\ell^2}{g^2}}\cos\theta}$$ We can also use $\ref{eell}$ to tell us something about the velocities at apicenter and pericenter, $v_a$ and $v_p$. Since the angular momentum is conserved, we would have that $\ell = r_av_a = r_bv_b$. Plugging this into equation $\ref{eell}$ and using $r_a=a(1+e)$ and $r_b=a(1-e)$, we can get $$v_a = \sqrt{\frac{g}{a}\frac{1-e}{1+e}}$$ $$v_p = \sqrt{\frac{g}{a}\frac{1+e}{1-e}}$$ $$v_p = v_a\frac{1+e}{1-e}$$ Or, if we substitute $r_a=a(1+e)$ and $r_p=a(1-e)$, we get the following nice equations: $$e = 1 - \frac{r_av_a^2}{g} = \frac{r_pv_p^2}{g} - 1\label{eevavp}$$ Everything is therefore determined once you have the initial conditions, either in the form of $r_a,v_a$ (remember that $r_av_a=r_pv_p$) or $\varepsilon,\ell$.

Lagrangian Formulation

Back to top

Starting with equations in the above section, the Lagrangian $L=T-V$ will be given by:

$$L=\half m(\dot{r}^2 + r^2\dot{\theta}^2) + \frac{GMm}{r}\label{elagrange1}$$ Using the Euler-Lagrange equations we have the following equation for the components $r$: $$\frac{d}{dt}m\dot{r} - mr\dot{\theta}^2 + \frac{GMm}{r^2}$$ and $\theta$: $$\frac{d}{dt}mr^2\dot{\theta}=0\label{eLcons}$$ The last equation says that the quantity $mr^2\dot{\theta}$ is constant, and indeed this is the angular momentum ($\vec{L}=\vec{r} \times \vec{p}$) for motion in a plane: $$L\equiv mr^2\dot{\theta}\label{eL}$$ Using this we can subsitute into the radial equation to get: $$\ddot{r} - \frac{L^2}{mr^3} + \frac{GM}{r^2}=0\label{eLr}$$

Conservation of angular momentum played a big role in connecting Newton's laws of motion to Kepler's. Here, the motion in a plane looks like the following figure:

As the particle moves, it sweeps out an area under the curve subtended by the angle $\delta\theta$. That area is given by

$\delta A = \half r\cdot r\delta\theta + \half \delta r\cdot r\delta\theta$ + higher order...

The 2nd term is 2nd order in the differentials, so we can drop that when we take the limit as $\delta\theta \to 0$, leaving

$\delta A = \half r\cdot r\delta\theta$

which means that

$\frac{\delta A}{\delta t} = \half r^2\frac{\delta\theta}{\delta t}$

or equivalently:

$\dot A = \half r^2\dot{\theta}^2\label{eAdot}$

If we multiply equation $\ref{eLcons}$ by $1/2m$, we would get

$\frac{d}{dt}\half r^2\dot{\theta}=0$

which means that $\dot A$ is conserved, which means that $\delta A = k\delta t$, where $k$ is some constant, which means that a body sweeps out equal areas in equal times, aka Kepler's 2nd Law.

We can use the Lagrange equations to find $x(t)$ and $y(t)$, but it is usaully more useful to find an equation for the "orbit", or in other words how the radius depends on angle: $r(\theta)$. To accomplish this, we will need to solve equation $\ref{eLr}$.

To do this, we use the chain rule:

$\frac{d}{dt} = \frac{d\theta}{dt}\frac{d}{d\theta}=\dot\theta\frac{d}{d\theta}$

Then, going back to equation $\ref{eL}$, we can substitute for $\dot\theta$ to get

$\frac{d}{dt} = \frac{L}{mr^2}\frac{d}{d\theta}$

We can then apply this to simplify:

$\ddot{r} = \frac{d}{dt}\frac{d}{dt}r = \frac{L}{mr^2}\frac{d}{d\theta}\frac{L}{mr^2}\frac{d}{d\theta}r$

This doesn't look like a simplification, until you make the following substitution: $$u\equiv 1/r\label{eur}$$ Then $\ddot r$ becomes

$-\frac{L^2}{m^2}u^2\frac{d^2u}{d\theta^2}$ and equation $\ref{eLr}$ becomes

$-\frac{L^2}{m^2}u^2\frac{d^2u}{d\theta^2} - \frac{L^2}{m^2}u^3 + GMu^2 = 0$

Multiplying each term by $m^2/(L^2u^2)$ and rearranging gives

$u'' + u = \frac{GMm^2}{L^2}$

where $u''$ means $d^2u/d\theta^2$. As we did above, we can collect symbols by defining $g\equiv GM$ and $\ell=L/m$ to give $$u'' + u = g/\ell^2$$ This equation is clearly sinosoidal, with solution

$u(\theta) = g/\ell^2 + B\cos(\theta - \theta_0)$ or going back to $r(\theta)$:

$$r(\theta) = \frac{1}{g/\ell^2 \pm B\cos(\theta - \theta_0)}\label{ertheta}$$

where $B$ and $\theta_0$ are initial conditions.

For the initial condition $\theta_0$, let's set $\theta_0=0$, which is pretty innocuous as a starting angle (see figure 2), and lets also set the origin at the left foci as in figure 3, which resolves the $\pm$ ambiguity in favor of $-$ (so that $r$ is a maximum at $r_0$). We then would have $r_0$, $v_0$, and $E_0$ evaluated at the same point (more on this later, see figure 4 below).

Figure 4. Ellipse centered at the left focal point with initial conditions $v_0$ and $r_0$ shown.

To solve for $B$, we can use the initial condition

$1/r_0 \equiv 1/r(0) = g/\ell^2 + B$

and go back to the other constant of the motion, the energy $E$, which is the sum $T+V$ and is given by

Since $E$ is constant, we can evaluate it anywhere, and the convenient place would be where $\dot{r}=0$. If you take the derivative of equation \ref{ertheta} and set it to 0, you will find that $\dot r=0$ happens when $\theta = 0$, and at that point $r=r_0$, which allows you to write

$$\varepsilon = \half r_0\frac{\ell^2}{r_0^4} - \frac{g}{r_0} = \frac{\ell^2}{2r_0^2} - \frac{g}{r_0}$$

Here is the trick: replace $1/r_0$ in the above equation for the energy with $B + g/\ell^2$, and solve for B. You should get

$$B = \frac{g}{\ell^2}\sqrt{ 1+\frac{2\varepsilon \ell^2}{g^2} }$$

Substituting that back into equation $\ref{ertheta}$ and rearranging gives:

$$r(\theta) = \frac{\ell^2/g}{ 1-\sqrt{ 1+\frac{2\varepsilon \ell^2}{g^2} }cos\theta }\label{eellipser}$$

This equation is clearly that of an ellipse: $-1\lt\cos\theta\lt +1$, therefore $(g/\ell^2+B)\lt\!1/r\lt(g/\ell^2-B)$, and $r$ is the same for $\theta = 90$ and $270$ degrees.

Equating equations $\ref{eellrtf}$ and $\ref{eellipser}$ allows you to solve for the eccentricity in terms of the initial conditions, giving $$e = \sqrt{ 1+\frac{2\varepsilon \ell^2}{g^2}}\label{eee}$$ (Remember that $\varepsilon\lt 0$ is constant as is $\ell\equiv L/m$.) Note that here, the eccentricity $e$ will also be given by the ratio of the distance between the two foci $f$ and the length of the major axis (the distance between there the ellipse intersects the $x$ axis).

Also by equating equations $\ref{eellrtf}$ and $\ref{eellipser}$, we can solve for the semi-major axis length $a$ to get a very interesting result: $$a = -\frac{g}{2\varepsilon}$$

This equation also says that $\varepsilon=-g/2a$, and of course $\varepsilon$ is a constant of the motion. To understand this better, let's first consider a circular orbit where $r$ and $v=r\dot\theta$ are both constants. The kinetic energy will be given by $T=\half mv^2$, and the potential by $V=-GMm/r$ so that total energy will be given by $E = \half mv^2 - GMm/r$. If we consider the centripetal acceleration, held into a circular orbit by a central gravitational force, we can write that

$\frac{mv^2}{r}=\frac{GMm}{r^2}$

or equivalently

$mv^2=\frac{GMm}{r}$

Substituting that into the above equation for $T$ gives $T=\half GMm/r$, or $T=-\half V$. This is called the "virial theorem" as applied to circular orbits. Plugging this new value for the kinetic energy into the equation for the energy gives:

$E = \half \frac{GMm}{r} - \frac{GMm}{r} = -\half\frac{GMm}{r}$

or equivalently

$\varepsilon = -\frac{g}{2r}$

for circular orbits, where here $r$ is the constant radius.

Same formula as for an elliptical orbit, except we use $a$ as the semi-major axis and write

Substituting $\varepsilon=-\frac{g}{2a}$ into equation $\ref{eee}$ gives $$e = \sqrt{ 1-\frac{l^2}{ga}}\label{eea}$$

Note that here we can see explicitly that $e\lt1$ for an ellipse.

This equation is interesting, but not as useful as equation $\ref{eee}$ when it comes to actually doing a calculation and simulation.

Real Two-body Dynamics

Back to top

So far, we've explored the orbits of one body of mass $m$ around another body of mass $M$, ignoring the motion of the latter mass. This is a good approximation for the condition $M\gt\gt m$, as in the case of satellites around the earth. But in general, we have to follow Newton's Third law on action/reaction.

The figure below sets up the generalized coordinate system, where $\vec{R_1}$ and $\vec{R_2}$ are the position vectors pointing to the centers of masses $m_1$ and $m_2$.

The kinetic and potential energy for the two masses will then be given by $$T = \half m_1 |\vec{\dot R_1}|^2 + \half m_2 |\vec{\dot R_2}|^2$$ and $$V = \frac{-Gm_1m_2}{|\vec{R_1}-\vec{R_2}|}$$ If we can write this Lagrangian in the same form as in equation $\ref{elagrange1}$ above, then we would be done. To do this, we make the following substitutions: $$\vec{r} = \vec{R_1}-\vec{R_2}\label{evecrr}$$ $$\vec{R} = \frac{m_1\vec{R_1} + m_2\vec{R_2}}{m_1+m_2}\label{evecRR}$$ These two new vectors are shown in red in the figure below.

Figure 5. Definition of relative vector $\vec{r}$ and CM vector $\vec{R}$.

where here $\alpha\equiv m_1/(m_1+m_2)$. The terms $\vec{R_1}\times\vec{R_1}$ and $\vec{R_2}\times\vec{R_2}$ vanish, leaving:

The only way that this can be true is if $\vec{R}$ ends at the connection between the two masses: QED.

We can now change from $\vec{R_1}$ and $\vec{R_2}$ to $\vec{R}$ and $\vec{r}$ in the kinetic and potential energies above to form the Lagrangian. The kinetic energy term becomes

where $M\equiv m_1+m_2$ and $\mu \equiv \frac{m_1m_2}{M}$.

The potential energy term will be

$$V = \frac{-Gm_1m_2}{|\vec{r}\!|}$$

so that the full Lagrangian is: $$L = \half M|\vec{\dot{R}}|^2 + \half\mu |\vec{\dot{r}}\!|^2 + \frac{GM\mu}{|\vec{r}\!|}$$ where we have used $M\mu = m_1m_2$.

What's going on here is clear: the vector $\vec R$ represents the center of mass, and $\vec r$ the relative position between the two masses, and so it's clear that 2 body motion can be described by the motion of the center of mass (CM) with total mass $M$ and position $\vec R$, and the motion within the center of mass with the mass $\mu$ and position $\vec r$.

In this case, let's assume that the CM is not changing (there are no forces on the system from outside). The motion is entirely internal and our final Lagrangian will be $$L = \half \mu |\vec{\dot{r}}\!|^2 + \frac{GM\mu}{|\vec{r}\!|}\label{elagrangeCM}$$ This reduces the problem of trying to find the motion of the 2 satellites to the motion of a single satellite around an infinitely larger satellite, given by equation $\ref{eellipser}$, and make the substitution $M\to m_1+m_2$ and $m\to\mu$. Once we have $r(\theta)$, we can use equations $\ref{er1cm}$ and $\ref{er2cm}$, set $\vec{R}=0$ (accomplished via initial conditions, see below), and get $$R_1(\theta)=\frac{\mu}{m_1}r(\theta)$$ $$R_2(\theta)=-\frac{\mu}{m_2}r(\theta)$$ In the first part of this chapter, we used the relative vector $\vec{r} =\vec{R_1}-\vec{R_2}$ because the potential goes like $1/|\vec{r}\!|$. However for the rest of this chapter, we are considering the dynamics of the satellites relative to the CM frame, and it's more convenient to define yet another set of vectors, $\vec{r_1}$ and $\vec{r_2}$, defined as directions relative to the CM position and shown in figure figure 6: $$\vec{R_1}=\vec{R}+\vec{r_1}\label{er1cm2}$$ $$\vec{R_2}=\vec{R}+\vec{r_2}\label{er2cm2}$$

Figure 6. Definition of vectors $\vec{r_1}$ and $\vec{r_2}$.

If we then plug these definitions into the definition for the vector $\vec{R}$ that points to the CM (equation $\ref{evecRR}$), we get the nice relation: $$M\vec{R}=m_1\vec{R_1}+m_2\vec{R_2}=m_1(\vec{R}+\vec{r_1})+m_2(\vec{R}+\vec{r_2}) = (m_1+m_2)\vec{R}+[m_1\vec{r_1}+m_2\vec{r_2}]=M\vec{R}+[m_1\vec{r_1}+m_2\vec{r_2}]$$ which can only be true if $$m_1\vec{r_1}+m_2\vec{r_2}=0\label{ercmcm}$$

Given that we now know how to describe the position of the CM, and the position of the satellites relative to the CM, we can now describe the velocities of the two satellites in a similar way. Starting with the diagram in figure 7 and the two velocities of the satellites relative to the origin, $\vec{V_1}$ and $\vec{V_2}$, we can form the analogous CM velocity analogous to equation $\ref{evecRR}$ via $$\vec{V} = \frac{m_1\vec{V_1} + m_2\vec{V_2}}{m_1+m_2}\label{evecVV}$$ and the relative velocities via: $$\vec{V_1}=\vec{V}+\vec{v_1}\label{ev1cm2}$$ $$\vec{V_2}=\vec{V}+\vec{v_2}\label{ev2cm2}$$ This is actually the same thing as considering the CM momentum, where we can define $P=MV=p_1+p_2$. Figure 7 shows how just like with the position vectors, the velocity vectors will add such that the CM velocity will lie on the vector $\vec{v}\equiv \vec{V_1} - \vec{V_2}$.

      Figure 7. Definition of vectors $\vec{V}$, $\vec{v_1}$ and $\vec{v_2}$. Note
that $\alpha\equiv m_1/(m_1+m_2)$ as discussed above.

Now we can discuss the 2-body dynamics in terms of the angular momentum. First form the total angular momentum $\vec{L_{tot}}$ about the origin, and expand using equations $\ref{er1cm2}$, $\ref{er2cm2}$, $\ref{ev1cm2}$, and $\ref{ev2cm2}$:

The first term in the last line can be written as: $\vec{R}\times(m_1\vec{V_1}+m_2\vec{V_2})=\vec{R}\times M\vec{V}$ using equation $\ref{evecVV}$ and $M\equiv m_1+m_2$. This is clearly the angular momentum $\vec{L_{CM}}$ of the CM relative to the coordinate origin: $$\vec{L_{CM}} = \vec{R}\times M\vec{V}\label{elcm}$$ The last 2 terms in that line are not the angular momentum of the two satellites relative to the CM, which would be given by $L_1 = \vec{r_1}\times m_1\vec{v_1}$ and etc for satellite 2, however we can show that the sum of these last 2 terms is the sum of the angular momenta about the CM. To see this,d we can expand those last two terms using equations $\ref{ev1cm2}$ and $\ref{ev2cm2}$. For satellite 1, we have: $$m_1\vec{r_1}\times\vec{V_1}=m_1\vec{r_1}\times(\vec{V}+\vec{v_1})= m_1\vec{r_1}\times\vec{V} + m_1\vec{r_1}\times\vec{v_1}$$ When we add this to a similar term for satellite 2, we can easily see that the terms $m_1\vec{r_1}\times\vec{V}+m_2\vec{r_2}\times\vec{V}=(m_1\vec{r_1}+m_2\vec{r_2})\times\vec{V}=0$ using equation $\ref{ercmcm}$ above. This leaves the 2 terms $m_1\vec{r_1}\times\vec{v_1}$ and $m_2\vec{r_2}\times\vec{v_2}$, which are in fact the angular momenta $\vec{L_1}$ and $\vec{L_2}$ relative to the CM: $$\vec{L_1} = \vec{r_1}\times m_1\vec{v_1}$$ $$\vec{L_2} = \vec{r_2}\times m_2\vec{v_2}$$ This means we can write the total angular momentum $\vec{L_{tot}}$ as: $$\vec{L_{tot}} = \vec{L_{CM}} + \vec{L_1}+ \vec{L_2}$$ There are several ways to eliminate the angular momentum of the CM via choosing the right initial conditions, but the easiest is to just place the origin at the CM. Then $\vec{R}=$, which will make $\vec{L_{CM}}=0$

Back to initial conditions. We can give one or both of the satellites an initial velocity and start the simulation. However, if the net momentum of the satellites does not sum to zero, then the CM frame will have a momentum and a velocity, and you then see that the whole simulation would travel up the page. Not convenient for simulation purposes! To fix this, we will insert a correction to the initial velocities something that will ensure that the overall CM velocity remains 0. In otherwords, a boost into the CM frame. This is easily done.

First, we choose an initial velocity for the 2 satellites, call it $U_{i0}$ and $U_{i1}$ (we leave off the vector notation to make things easier to follow). Then, we make a boost (add velocity) $\Delta U$ to each one: $$V_{i0} = U_{i0} + \Delta U$$ $$V_{i1} = U_{i1} + \Delta U$$ such that we would have $m_1V_{i0}+m_2V_{i1}=0$. Following this through we will come up with $$\Delta U = -\frac{m_1U_{i0}+m_2U_{i1}}{M}$$

In the simulation below, you can vary the masses of the red and blue satellites, drag them both horizontally to change their locations, and drag the end of the arrow up and down to change the initial velocity of the red satellite. The green ball shows the CM location, and the cross shows the coordinate system center.

Back to top