365 Days of Astronomy: February 15th: Astronomy Word of the Week : Barycenter

I wrote a bit on this topic before (see Geocentrism vs. the Barycenter). Here's a follow-up, but with more mathematical graphic support.

Here we'll explore some more details about the two-body problem using Newton's laws and gravitation. An n-body code solves Newton's equations by computing the changing positions of multiple objects subjected to various forces. In this case, the forces used are the gravitational forces between the two bodies computed by Newton's Law of Gravitation:

Here M1 & M2 are the two masses and r12 is the distance between the. 'G' is the Newtonian gravitational constant. I provided movies of a simple 2-body problem in my first post on the n-body code (Doin' Real Science: Simulating Particles).

In this post, I'll examine the 2-body problem in more detail.

**Setup**

For these types of problems in science, we must first setup *Initial conditions* for the problem. The program will then compute the motions forward in time, updating velocities and positions based on the forces on the objects (in this case, gravity).

Let's start with two simple objects, which we call Primary and Secondary,

Primary, or mass1 = 10 solar masses

Secondary, or mass2 = 1 solar masses

I'll specify the object speeds in astronomical units per year (AU/yr) and the masses relative to the mass of the Sun. To convert the speed to more familiar units:

1 au/yr = (149.598e6 km)/(365.25*24*3600 s) = 4.74 km/s

= (92.955e6 miles/(365.25*24 hr) = 10604 miles per hour

So to convert my speeds to more familiar values, multiply the AU/year value by 4.74 to get kilometers/second or 10604 to get miles/hour or 2.9455 to get miles per second.

In the computer program, I give the Primary and Secondary objects initial positions and velocities consistent with this setup. We call them INITIAL positions and velocities because they will change with time as the force acts on them. For simplicity, I place the Primary object at the origin (0, 0, 0), the *center* of our universe as we have defined it, and give it an initial velocity of zero (0, 0, 0).

I place the Secondary object at (1AU, 0, 0) and give it a velocity in the y-direction (perpendicular to the initial line between the Primary & Secondary, which defines the x-axis). For this case, I set the magnitude of the velocity to 60% of the speed of a circular orbit at the distance of 1AU.

To make the problem more general, I add an additional velocity of 1AU/yr in the y-direction for both of these objects. Here's a printout of the table of initial parameters:

Q Object Mass ------------------ position -------------------

-------------------- velocity ---------------------

1) 0 10 M 1.00000e+00 Msun x=( 0.000000e+00, 0.000000e+00, 0.000000e+00) AU

v=( 0.000000e+00, 1.000000e+00, 0.000000e+00) AU/yr

2) 0 1 M 1.00000e-01 Msun x=( 1.000000e+00, 0.000000e+00, 0.000000e+00) AU

v=( 0.000000e+00, 4.769911e+00, 0.000000e+00) AU/yr

The position vector is 'x' and the velocity vector is 'v'. The 'Q' column specifies that the planets have no significant electric charge.

**Analytic Solution**

Using Kepler's laws, we can compute a few specifics about our orbit before running the simulation. This will help provide a check on the accuracy of the simulation run. N-body problems accumulate small errors at each step of the computation so as they run for longer and longer times, their accuracy declines. Here's some Keplerian parameters based on the input values.

Velocity of Center-of-Mass: (0.000000, 1.342719, 0.000000) AU/yr

Orbital Elements

Period: 0.440723 yr

Semi-major Axis: 0.597826 AU

Semi-minor Axis: 0.442326 AU

Eccentricity: 0.672727

Periapsis: 0.195652 AU

Apoapsis: 1.000000 AU

Because we chose a speed for the Secondary object that is less than the circular orbit velocity at 1AU, we see that our initial position will correspond to the apoapsis (or farthest distance) of our orbit from the Primary object. When we start running the simulation, we expect the Secondary object to fall *towards* the Primary object.

**Numerical Solution**

So we run the simulation in the coordinate system as we have defined. Here's a plot of the output dataset.

(download the entire movie here.)

The red curve tracks the motion of the Primary and the green curve tracks the motion of the Secondary. The black dashed curve corresponds to the motion of the center of mass of the two bodies, computed at each time we've plotted.

We note several things about the motion

- the primary object, while it started out at the origin with zero velocity, did not stay at the origin. This is a problem for Geocentrists since it will happen for all objects moving under mutual forces.
- both masses move around the center of mass point, the CM lying on the line between them. These are the motions measured when searching for planets by the radial velocity method (wikipedia).
- Remember that I gave both objects an initial motion of 1AU/yr in the y-direction. Yet the CM is moving slightly faster, 1.34 AU/yr. This is again due to the fact that gravitation is a mutual force and Newton's laws require an action & reaction.

Kepler's laws are actually defined in the reference frame were one of the objects is at rest, as I noted in the earlier post. So we must compute the RELATIVE positions of the objects. Here we do the case of the motion relative to the Primary (more massive) body.

x(1,ref), y(1,ref), z(1,ref) = position of object 1, in the coordinate system where the origin is object

*ref*.

For the Primary body, this becomes:

x(1,ref) = x(1) - x(1) = 0.0

y(1,ref) = y(1) - y(1) = 0.0

So the Primary resides at the origin in this new coordinate system. Similarly, we compute the new position of the Secondary object at each time step:

x(2,ref) = x(2) - x(1)

y(2,ref) = y(2) - y(1)

Here's a plot of the result:

(download the entire movie here.)

With this plot, we can more easily compare the plot to the computed orbital period, semi-minor axis, and periapsis in the table above with the aid of the grid overlay. The secondary (green) completes the orbit between t=0.440 yr and t=0.441 yr, consistent with P=0.4407 yr computed from the input values. We also observe the speed-up of the Secondary as it approaches the Primary mass, as predicted by Kepler's Laws. I also plot the center-of-mass (black dashed line) in this reference frame. In this example, the CM moves in an ellipse around the primary mass.

**Other Frames of Reference**

Because mathematically, the frames are equivalent, we can also convert to the frame of the less massive object as well!

(download the entire movie here.)

This time, the CM appears to move around the Secondary. We can also plot the motion of both objects around the center of mass.

(download the entire movie here.)

Mr. Martin claimed that Newton's laws and the barycenter were inconsistent (see Invalidations of Newton's and Kepler's orbital mechanics):

Here’s some simple invalidations of Newton’s and Kepler’s orbital mechanics –But now we can see these motions are perfectly compatible. The only evidence Mr. Martin presents for his claim is his 'say so'.

1. The orbital mechanics of Newton dictates the earth orbits the sun’s center of mass in an ellipse, yet Newtonian mechanics states the earth also orbits the solar system barycenter. As the solar system barycenter is almost never at the center of mass of the sun, then the earth simply cannot be orbiting the sun in an ellipse. Therefore Newton’s principle of barycentric motion invalidates Kepler’s laws of elliptical motion.

In fact, I can also plot the motion in the original input frame (green & red) and overlay the motion of the secondary mass relative to the primary mass (magenta).

(download the entire movie here.)

Notice that from the input frame, the Keplerian orbit is not fixed in space, but appears to be carried by the primary mass as the primary moves around the center of mass.

N-body codes are a regular tool for astronomy and celestial navigation. While I don't include many of the smaller forces (planetary oblateness, relativisitic effects, etc.) my tool is sufficiently accurate to demonstrate many principles in the application of Newton's laws with gravitation.

I've got a number of simulations of Geocentrists' claimed "problems" that I'm writing up for future posts.