Adding a gravitational body

For most normal uses of VLATI, this is not required; approximate ephemerides are already provided, and the Sun and Moon are already integrated. That being said, if you would still like to continue with this, keep reading.

It is recommended to read the data flow section to familiarize yourself with how data moves through VLATI.

This is difficult; I recommend becoming familiar with the program’s “shape” before doing this.

main.py

1. Mass scalar and t=0 position vector

First, add the t=0 position vector. This denotes where your body will be at t=0 (this doesn’t really have too much influence on the physics of the simulation if you’re using a motion function or SPICE/NAIF). You will find it under the comment marked ### R VECS ### around line 31. Already there, you will find variables named r_sun, r_moon, and r_earth. Add yours in. For this tutorial, the body we’ll add can be called mars.

### R VECS ###

r_moon = np.array([384398861.0, 0.0, 0.0])
r_earth = np.array([0.0, 0.0, 0.0])
r_sun = np.array([149599999999.7966, 0.0, 0.0])

r_mars = np.array([posx, posy, posz])

Make sure to add decimals to this to indicate to VLATI and NumPy that it’s a vector of floating-point numbers and not integers.

Next, directly underneath, you’ll need a mass scalar for your body in milligrams. No scientific notation is allowed, and it is discouraged to copy and paste mass values.

Just kidding. Use kilograms. Scientific notation (Pythonic and not Pythonic) is, of course, allowed.

m_moon = 7.34767309 * (10**22)
m_sun = 1.989 * (10**30)

m_mars = 6.41693 * (10**23)

2. Initialization in bds

bds is an array of tuples. Each tuple contains a body’s mass scalar and its position vector (the one that you initialized earlier - it’s important to initialize it, but the t=0 value isn’t super important to the physics).

Add your tuple to bds:

bds = [ (r_earth, 5.972 * (10**24)), (r_moon, m_moon), (r_sun, m_sun), (r_mars, m_mars) ]

It will now be automatically factored into the calculations.

3. Motion

VLATI uses NAIF’s SPICE to get planetary positions relative to the ecliptic (the plane on which Earth orbits the Sun). You can either provide your own function, or use SPICE to get the position of the body from the ephemeris file and rewrite the body’s tuple in bds:

### COMPUTE NEW POS OF PLANETS ###
  ...

r_mars = spice.spkpos("MARS", ET, AXIS_REFR, "NONE", "EARTH")[0] * 1000
bds[3] = (r_mars, m_mars)

Let’s break some points of confusion down:

• We get index 0 of the SPICE spkpos because it returns an object that contains position, lighttime, and we only want the position.
• We multiply it by 1000 because:
  • It gives results in kilometers, and VLATI operates in meters.
  • It is a NumPy array, so we can do that (vector-scalar multiplication).

All items below this line are optional and are not strictly necessary for trajectory calculation.

---

4. Static Matplotlib visualization

We can add a little dot to the static MPL visualization after the integration is complete to show where the body is at the end of the simulation.

### MPL BODY COLORS ###
  ...

dotm = ax.scatter(r_moon[0], r_moon[1], r_moon[2], c='black', marker='o')

dot_mars = ax.scatter(r_mars[0], r_mars[1], r_mars[2], c='orange', marker='o')

This will make a little Mars dot appear at its last position in the simulation.