REBOUND can add particles to simulations by obtaining ephemerides from NASA's powerful HORIZONS database. HORIZONS supports many different options, and we will certainly not try to cover everything here. This is meant to serve as an introduction to the basics, beyond what's in Churyumov-Gerasimenko.ipynb. If you catch any errors, or would either like to expand on this documentation or improve REBOUND's HORIZONS interface (rebound/horizons.py), please do fork the repository and send us a pull request.
Adding particles
When we add particles by passing a string, REBOUND queries the HORIZONS database and takes the first dataset HORIZONS offers. For the Sun, moons, and small bodies, this will typically return the body itself. For planets, it will return the barycenter of the system (for moonless planets like Venus it will say barycenter but there is no distinction). If you want the planet specifically, you have to use, e.g., "NAME=Pluto" rather than "Pluto". In all cases, REBOUND will print out the name of the HORIZONS entry it's using.
You can also add bodies using their integer NAIF IDs: NAIF IDs. Note that because of the number of small bodies (asteroids etc.) we have discovered, this convention only works for large objects. For small bodies, instead use "NAME=name" (see the SMALL BODIES section in the HORIZONS Documentation).
import rebound
sim = rebound.Simulation()
sim.add("Sun")
## Other examples:
# sim.add("Venus")
# sim.add("399")
# sim.add("Europa")
# sim.add("NAME=Ida")
# sim.add("Pluto")
# sim.add("NAME=Pluto")
sim.status()
Searching NASA Horizons for 'Sun'... Found: Sun (10) --------------------------------- REBOUND version: 4.0.3 REBOUND built on: Jan 12 2024 08:52:18 Number of particles: 1 Selected integrator: ias15 Simulation time: 0.0000000000000000e+00 Current timestep: 0.001000 --------------------------------- <rebound.particle.Particle object at 0x1120e7440, m=0.9999999999950272 x=-0.007850021941778935 y=-0.003071795945987032 z=0.00020917793069980578 vx=0.00029283533072769824 vy=-0.0003990033211746936 vz=-2.86982605624576e-06> --------------------------------- The following fields have non-default values: N: < 0 --- > 1 python_unit_l: < 0 --- > 4097248214 python_unit_m: < 0 --- > 2145773914 python_unit_t: < 0 --- > 1864791206 rand_seed: < 886170 --- > 676381 particles: > (128 bytes, values not printed)
Currently, HORIZONS does not have any mass information for solar system bodies. rebound/horizons.py has a hard-coded list provided by Jon Giorgini (10 May 2015) that includes the planets, their barycenters (total mass of planet plus moons), and the largest moons. If REBOUND doesn't find the corresponding mass for an object from this list (like for the asteroid Ida below), it will print a warning message. If you need the body's mass for your simulation, you can set it manually, e.g. (see Units.ipynb for an overview of using different units):
sim.add("NAME=Ida")
print(sim.particles[-1]) # Ida before setting the mass
sim.particles[-1].m = 2.1e-14 # Setting mass of Ida in Solar masses
print(sim.particles[-1]) # Ida after setting the mass
Searching NASA Horizons for 'NAME=Ida'... Found: 243 Ida (A884 SB) <rebound.particle.Particle object at 0x1120e7540, m=0.0 x=1.7283804595871177 y=-2.3772031668640636 z=-0.017147599198848423 vx=0.4522491120827006 vy=0.35453102610682247 vz=0.01092312049884117> <rebound.particle.Particle object at 0x1120e7540, m=2.1e-14 x=1.7283804595871177 y=-2.3772031668640636 z=-0.017147599198848423 vx=0.4522491120827006 vy=0.35453102610682247 vz=0.01092312049884117>
/Users/dtamayo/Documents/workspace/rebound/rebound/horizons.py:167: RuntimeWarning: Warning: Mass cannot be retrieved from NASA HORIZONS. Set to 0.
warnings.warn("Warning: Mass cannot be retrieved from NASA HORIZONS. Set to 0.", RuntimeWarning)
Time
By default, REBOUND queries HORIZONS for objects' current positions. Specifically, it caches the current time the first time you call rebound.add, and gets the corresponding ephemeris. All subsequent calls to rebound.add will then use that initial cached time to make sure you get a synchronized set of ephemerides.
You can also explicitly pass REBOUND the time at which you would like the particles ephemerides:
sim = rebound.Simulation()
date = "2005-06-30 15:24" # You can also use Julian Days. For example: date = "JD2458327.500000"
sim.add("Venus")
sim.add("Venus", date=date)
sim.status()
Searching NASA Horizons for 'Venus'... Found: Venus Barycenter (299) (chosen from query 'Venus') Searching NASA Horizons for 'Venus'... Found: Venus Barycenter (299) (chosen from query 'Venus') --------------------------------- REBOUND version: 4.0.3 REBOUND built on: Jan 12 2024 08:52:18 Number of particles: 2 Selected integrator: ias15 Simulation time: 0.0000000000000000e+00 Current timestep: 0.001000 --------------------------------- <rebound.particle.Particle object at 0x1120e7340, m=2.447838287784771e-06 x=-0.5245660051855766 y=-0.5089797719750951 z=0.023076250697312093 vx=0.8146329772804595 vy=-0.8461175012422552 vz=-0.05860561893434425> <rebound.particle.Particle object at 0x1120e75c0, m=2.447838287784771e-06 x=-0.6612400840226685 y=0.26970377066078294 z=0.04197305542931168 vx=-0.445669425641083 vy=-1.0953912971910784 vz=0.010725839790529726> --------------------------------- The following fields have non-default values: N: < 0 --- > 2 python_unit_l: < 0 --- > 4097248214 python_unit_m: < 0 --- > 2145773914 python_unit_t: < 0 --- > 1864791206 rand_seed: < 829806 --- > 210256 particles: > (256 bytes, values not printed)
We see that the two Venus positions are different. The first call cached the current time, but since the second call specified a date, it overrode the default. Any time you pass a date (with time in UTC), it will overwrite the default cached time, so:
sim = rebound.Simulation()
date = "2005-06-30 15:24" # You can also use Julian Days. For example: date = "JD2458327.500000"
sim.add("Venus", date=date)
sim.add("Earth")
sim.status()
Searching NASA Horizons for 'Venus'... Found: Venus Barycenter (299) (chosen from query 'Venus') Searching NASA Horizons for 'Earth'... Found: Earth-Moon Barycenter (3) (chosen from query 'Earth') --------------------------------- REBOUND version: 4.0.3 REBOUND built on: Jan 12 2024 08:52:18 Number of particles: 2 Selected integrator: ias15 Simulation time: 0.0000000000000000e+00 Current timestep: 0.001000 --------------------------------- <rebound.particle.Particle object at 0x1120e71c0, m=2.447838287784771e-06 x=-0.6612400840226685 y=0.26970377066078294 z=0.04197305542931168 vx=-0.445669425641083 vy=-1.0953912971910784 vz=0.010725839790529726> <rebound.particle.Particle object at 0x1120e72c0, m=3.0404326489511185e-06 x=-0.5575682509295975 y=0.813513379197337 z=0.00016705166681997484 vx=-0.8456273406927053 vy=-0.5626729140479645 vz=3.329436635861762e-05> --------------------------------- The following fields have non-default values: N: < 0 --- > 2 python_unit_l: < 0 --- > 4097248214 python_unit_m: < 0 --- > 2145773914 python_unit_t: < 0 --- > 1864791206 rand_seed: < 433115 --- > 842290 particles: > (256 bytes, values not printed)
would set up a simulation with Venus and Earth, all synchronized to 2005-06-30 15:24 UTC. All dates should either be passed in the format Year-Month-Day Hour:Minute or in JDxxxxxxx.xxxx for a date in Julian Days.
For reference HORIZONS interprets all times for ephemerides as Coordinate (or Barycentric Dynamical) Time.
Origin
REBOUND queries for particles' positions and velocities relative to the Solar System barycenter
sim = rebound.Simulation()
sim.add("Sun")
sim.status()
Searching NASA Horizons for 'Sun'... Found: Sun (10) --------------------------------- REBOUND version: 4.0.3 REBOUND built on: Jan 12 2024 08:52:18 Number of particles: 1 Selected integrator: ias15 Simulation time: 0.0000000000000000e+00 Current timestep: 0.001000 --------------------------------- <rebound.particle.Particle object at 0x1120e7e40, m=0.9999999999950272 x=-0.007850021941778935 y=-0.003071795945987032 z=0.00020917793069980578 vx=0.00029283533072769824 vy=-0.0003990033211746936 vz=-2.86982605624576e-06> --------------------------------- The following fields have non-default values: N: < 0 --- > 1 python_unit_l: < 0 --- > 4097248214 python_unit_m: < 0 --- > 2145773914 python_unit_t: < 0 --- > 1864791206 rand_seed: < 629898 --- > 447379 particles: > (128 bytes, values not printed)
The reference plane is the ecliptic (Earth's orbital plane) of J2000 (Jan. 1st 2000 12:00 GMT), with the x-axis along the ascending node of the ecliptic and the Earth's mean equator (also at J2000).