Parameters for Solar system simulation

1. Nov 16, 2013

CAF123

1. The problem statement, all variables and given/known data
I am constructing an N-body simulation in Java that will ultimately model the solar system. I am testing the reliability of it by slowly increasing the number of bodies in the simulation. So I am starting off with the Sun, Mercury and Venus. I am looking to find the initial position, velocity and masses of these three planets. Masses are easy enough to find, however I still need to find the initial position and velocity.

2. Relevant equations
I am told to find the initial velocity of the planets via, $$v = \sqrt{\frac{\mu}{r}},$$where $\mu = m_1 + m_2$ and $G = 1$ and I think $r$ is the distance of body 1 from body 2.

3. The attempt at a solution
If $\mu = m_1 + m_2$ then if one of the masses is the Sun, this will dominate the other planet. I think I can find the mean distance of each of Venus and Mercury from the Sun online. (Can I take this to be r?) So I can use the above formula to obtain the initial velocity. I am a bit confused though because this is a 3 body system, so how should I compute the resultant initial velocity using that formula?

Also, I can see that this eqn was derived redefining G=1, (that makes the equation dimensionally correct) but this does not alter the units of the quantities.

Many thanks

2. Nov 16, 2013

Staff: Mentor

Spend a bit of time tinkering with the JPL Horizons Web Interface.

The Horizons system allows you to generate ephemerides for any body in the solar system from just about any viewpoint you can name (use @sun to generate an ephemeris for an "observer" at the center of the Sun). Pick a particular date for your epoch, pick a target body (say Venus), and out pops an ephemeris with time-stamped coordinates as "seen" by the observer.
[STRIKE]
I think it's possible to change the outputs so that you can get rectangular coordinates, too.[/STRIKE] Nope. I checked.

Good luck!

Last edited: Nov 16, 2013
3. Nov 16, 2013

Staff: Mentor

That is good enough if you are not interested in details of the orbits (they are not circular, for example). If you want realistic orbits, use actual velocities and positions for the planets, there are databases that provide those values.
The forces between the planets are tiny compared to the forces between planets and sun.
It used a unit system where the units match. This can have weird consequences (like masses expressed as powers of length and time), but it works.

4. Nov 16, 2013

CAF123

Do you have reference to a link? I tried searching for this, but I could only find the current velocities/position. For what I am doing, I think I have to use the expression given, so I guess the mean distance of the planets from the Sun will do.

This makes sense, so then I could just apply the equation twice - for Mercury around the Sun and for Venus around the Sun.

So I would need to find the masses of the planets in terms of length and time. So [v] = m/s and [G] = [m3kg-1s-2]. Hence I should express my masses like [mass] = [length3 time-2] How would I obtain values like this so I can then use my formula?

5. Nov 16, 2013

D H

Staff Emeritus
You didn't check well enough. Change the "Ephemeris Type" (the first selection) from "Observer" to "Vector table".

6. Nov 16, 2013

Staff: Mentor

Oh ho! Very nifty. Hadn't seen that before. Good info! Thanks.

7. Nov 16, 2013

CAF123

Thank you, I presume what they mean by 'initial velocity' and 'initial position' are those quantities when the planets begun their orbits around the Sun. So what epoch is this?

My code is implemented in rectangular coordinates (that is I have seven parameters for each particle - label, 3 position coordinates, 3 velocity coordinates).

I am still trying to navigate around the Horizon interface and when I have more time, I'll look at the user guide. My instructor actually provided the equation I mentioned in the OP, so I think I should use it. Thanks.

8. Nov 16, 2013

Staff: Mentor

You can pick any instant in time as your epoch. Use the same instant for all the bodies and the positions and velocities will all be accurate to start your simulation from that instant.

If you assume circular orbits then that's fine. But your orbits will all be circular rather than elliptical, and you still have the problem of choosing suitable starting locations. Using actual ephemeris positions and velocities will give you the real thing.

Here's some sample output for an ephemeris for Mars:

Attached Files:

• Fig1.jpg
File size:
38.9 KB
Views:
144
9. Nov 16, 2013

Staff: Mentor

If you use F=Mm/r^2 as formula this is fine, if you use F=GMm/r^2 you can multiply all "mass" values with some power of G and whatever to get masses in kg.

10. Nov 16, 2013

D H

Staff Emeritus
That's exactly right. That's what your μ was expressed in in the opening post. Look at the units. You had $v=\sqrt{\mu/r}$. That means μ has to have dimensions of length3/time2. Solar system modelers don't use G*M. They use μ. The reason is that a planet's "standard gravitational parameter", or mu, can be inferred from planetary observations. Planetary mass is computed from μ/G. Here's the problem: G is one of the least well known physical constants. This uncertainty in G means that it's best if G never enters into the equations.

11. Nov 16, 2013

Staff: Mentor

This also means M is poorly known (well, relative to other measurements). Expressed in other words: the uncertainty on the mass of our sun is larger than the mass of earth.
We really need a space-based experiment to determine G...

12. Nov 17, 2013

CAF123

$\mu = G(m_1 + m_2) = m_1 + m_2 \approx m_1 = m_s$, where $G =1$ and $m_s$ is the mass of the Sun and this dominates all the other masses under consideration. How do I obtain this value of the mass of the Sun in terms of units of length^3/time^2?

]

13. Nov 17, 2013

Staff: Mentor

You know M, you want to know GM. Just multiply the mass with G (expressed in the SI system).

14. Nov 17, 2013

D H

Staff Emeritus
Google the term standard gravitational parameter. That's the name for the product of G and mass.

15. Nov 17, 2013

CAF123

I see, I was overcomplicating things. Now that I have the magnitude of the initial velocity for Mercury around the Sun (computed using the expression in OP) in order to put this into the program I need to supply the x,y,z components of initial velocity in an input file.

To find these, I was thinking I could use gneill's suggestion of the Horizon interface since I can see a key on the ephemerides about half way down 'JDCT X Y Z <newline> VX VY VZ'. If I understand this correctly, the first row lists the components of position and the second row lists the components of velocity. I can use this to get my initial position. Calculating gives a 0.14% difference in the value of the initial position using Horizon and the value of the mean distance.

Similarly, for the velocity I obtain a 6.9% difference. I think the only reason for this difference is the propagation of the error when I take the sqrt of $\mu/r$ since the differences in r are barely noticeable. Is this right? And thanks gneill for providing that link - it really helped.

Thanks everyone.

16. Nov 17, 2013

Staff: Mentor

Another (small) source of difference would be that the Horizons data has the solar system barycenter as its origin. This is close to, but not identical with, the center of the Sun (occasionally it lies outside of the Sun). The barycenter is determined by the instantaneous location of all the planets, not just the one you're looking at at the moment!

If you're modelling the Sun as free body, too, you'll want to give it initial location and velocity parameters.

17. Nov 17, 2013

CAF123

In the coordinate origin data entry, I changed this from the default 'Solar system barycenter' to the Sun. But see comment below; I may change this.

Since I set the Sun at the origin, I think I can let the initial position of the Sun to be at <0,0,0>. How would I obtain the velocity of the Sun? This will change depending on the point of observation. Since I set the Sun as my origin, the velocity of the Sun in this 'frame' is trivially zero. Do you think it is better, therefore, to assign values to everything with the Solar system barycenter as my origin?

Thank you.

18. Nov 17, 2013

Staff: Mentor

I would say yes. Otherwise you may find your entire simulation "creeping off the page" as time goes on. Why? Because your coordinate system origin won't be located at the center of momentum of the system, wherever it turns out to be once the 'physics' get going in the simulation. Unless the Sun happens to be at the true center of momentum for the ensemble, it's going to wander.

Two choices:
1) Set initial positions and velocities for all bodies, including the Sun, with the barycenter as the origin.
2) Set initial positions and velocities for all bodies and assume the Sun is motionless at the origin. THEN run through the data to determine the location and velocity of the center of mass of the ensemble. THEN tweak all the velocity values to make the center of mass stationary.

19. Nov 17, 2013

CAF123

So do you mean to say about the barycenter the total mechanical momentum of the system is conserved? Choice 1 seems simpler given I can use the Horizon interface.

20. Nov 17, 2013

Staff: Mentor

It's a closed system, so yes, momentum is conserved and the motion of the center of mass is inertial (constant velocity). You can choose your coordinate system to 'ride along with' this center of momentum, thus keeping your simulation "pinned" to the center of the coordinate system. Otherwise the center of mass of the system would drift linearly with time. It can be annoying to watch an interesting simulation gradually drift off the screen!

Yup, choice 1 is a good way to go. Even so, I would still recommend, as a matter of course, tweaking the velocities after assembling the system to zero any residual motion of the center of mass. This is because perhaps not all the existing solar system bodies will be included for every run, and the initial positions/velocities given by Horizons take all bodies into account. Imagine if Jupiter and the other gas giants are left out but the Horizon's location and velocity for the barycenter includes their influence...

Just find velocity vector for the center of mass (add all the momentum vectors, divide by total mass to yield velocity vector of center of mass). Subtract this vector from all the initial velocity vectors. Done.