More on Barycenters

Picking up a bit where I left off. According to someone with a Phd on a forum that knows a great deal about astrophysics, there was absolutely nothing wrong with my rendition of a barycenter orbit. The shifting of the barycenter if no "counter velocity" is established for the more massive body will occur because momentum must be conserved. Makes sense.

Anyway, here is a representation of the same barycenter orbit, the white one leaving the visual representation (the frame of reference) as it is, the orange one "correcting" the positions so that it appears that the more massive body is "stationary". It is a visual trick, they are really the same orbit.



And here's the "correcting" code.

barycenter_length = radius * (charon_mass / (charon_mass + pluto_mass));
bary_x = sin(angle_to_charon) * barycenter_length;
bary_y = cos(angle_to_charon) * barycenter_length;

Then, this code keeps the more massive body centered...

pluto_x + bary_x
pluto_y + bary_y

And so on for Charon. That is, by the way, the formula for establishing a barycenter;

Bc = Radius * (lesser_mass / (lesser_mass + greater_mass))

I am sort of up to my eyes in work this month, so progress will be slow.

'Till the next!

Gravitational interaction 2D

Back to 2D for simplicity's sake, now. I thought I might play a bit with interactions between bodies of less drastically different masses than, say, a baseball orbiting Earth. That is, I wanted to see a barycenter in action. The classical example of this is the "mutual orbit" of Charon and Pluto.

Now, there's a Newtonian formula that goes into this, but I thought I might take a more grass roots approach to the problem. In my first post, I determined how I might make a body of indeterminate mass (that is; relatively insignificant) orbit around a large mass (like Earth's) with a simplified formula. In this one, I will apply that same formula to a pair of bodies of much similar masses. Charon is 0.116 of Pluto's mass, approximately.

To speed things up a bit, I maintained the masses but shortened the semi major axis. I also upped the time interrupt to 60 frames per second, with each frame counting as a second (ie; time factor X60).

// The variable declarations...
const double grav_const = 6.6742e-11;
const double charon_mass = 1.52e21;
const double pluto_mass = 1.31e22;
double radius;

double charon_x = 2.75e6; // A shortened, arbitrary "axis"
double charon_y = 0;
double charon_vel_x = 0;
double charon_vel_y = 530; // An estimated sample velocity

double pluto_x = 0;
double pluto_y = 0;
double pluto_vel_x = 0;
double pluto_vel_y = -35;

// There is now going to be two separate gravity calculations,
// based on each bodies' mass, each attracting the other.double grav_accel_by_charon = 0;
double grav_accel_by_pluto = 0;

// Each angle could just as easily be done by calculating one and subtracting PI
// from the other, but I am going ahead and calculating each separatelydouble angle_to_charon;
double angle_to_pluto;

// What follows goes inside the time interrupt loop...
angle_to_charon = atan2f((charon_x - pluto_x), (charon_y - pluto_y));
angle_to_pluto = atan2f((pluto_x - charon_x), (pluto_y - charon_y));
radius = sqrtf(pow((charon_x - pluto_x), 2) + pow((charon_y - pluto_y), 2));

grav_accel_by_pluto = (grav_const * (pluto_mass / pow(radius, 2)));
grav_accel_by_charon = (grav_const * (charon_mass / pow(radius, 2)));

charon_vel_x += (sin(angle_to_charon) * grav_accel_by_pluto);
charon_vel_y += (cos(angle_to_charon) * grav_accel_by_pluto);

pluto_vel_x += (sin(angle_to_pluto) * grav_accel_by_charon);
pluto_vel_y += (cos(angle_to_pluto) * grav_accel_by_charon);

charon_x -= charon_vel_x;
charon_y -= charon_vel_y;

pluto_x -= pluto_vel_x;
pluto_y -= pluto_vel_y;

Now that sort of had the desired effect, except that that "Pluto's" counter velocity (-35) was off by a bit. To further complicate the matter, I did not fully calculate the required velocities to make a true circular (non elliptical) orbit of Charon, I just plonked in a rough estimate (530). The result was that, yes, a mutual orbit did occur, but because of unbalanced velocities it crept up the screen as it progressed. In reality, Charon and Pluto are almost perfectly harmonized, Charon having an almost 0 eccentricity around the barycenter. I am not discounting yet that the phenomenon that occurred in my demo can actually happen (screen shots a to c in Figure 3). After all, Charon and Pluto have had, literally, ages to stabilize themselves. I am not going to sit in front of my computer for the next 15 years running this program without interruption to see if my model stabilizes, too. Now, come on! I can, however, dramatically speed up the execution and see what happens, but that some other time. What does seem apparent is that both example are maintaining constant orbits, just that one shifts in space as there is an unbalanced "counter velocity".

For the moment, I just used the "sledge hammer" approach and set the velocities by the respective masses. As Charon is 0.116 of Pluto's mass, I multiplied "my" Charon's velocity (530) and tried that as Pluto's "counter velocity". The results were an improvement (screen shots d to f in Figure 3). The small error is now probably due only to the slight eccentricity.

Fig 3: For description, see text above...

Not bad for starters, in any case.

Post Edit: Hmmm. I do now think there is something not right, however. Look what I just found. That is more or less what I am aiming at.