Collisions detected with unwanted particle

[sarahthiele]

I'm using REBOUND (in Python, on Linux) to integrate a system of satellites and orbital debris around Earth and detect close encounters amongst them (i.e. I want to detect satellite--debris close encounters, and debris--debris encounters, but NOT between satellites or debris with Earth), where I'd be integrating the system for up to maybe a year. I'm using WHFast to integrate, currently using a small timestep (2e-9) to make sure I'm catching everything, although I can maybe raise this later.

There's about 16000 particles that have perigees up to ~8300 km (altitudes <= 2000 km) with a wide range of eccentricities, mostly between 0. and 0.8.

To track these what I've been doing is using the collision detection of REBOUND using the "tree" method and setting the radius of all particles (I add satellites and debris as test particles with 0 mass) as half the distance that I want to count as a "close encounter". For example, right now I'm setting all of them to have a radius of 5km so a close encounter is when two particles come within 10km of each other.

However, I'm finding that sometimes rebound tracks collisions between Earth and the particles, even though they're not within 10km of each other. The orbits also seem to become messed up because when I print out the altitude of the test particle involved in a collision, it's less than that of the Earth (I checked to ensure all objects have perigees > R_earth + 100 km). Weirdly this only happens when I include both satellites and orbital debris in the simulation, and doesn't happen when I'm only integrating satellites. I'm not sure if it's because the debris has more extreme orbits (higher eccentricities, etc.) so is throwing things off? I've also written it such that Earth is not involved in the collision detection, nor are the test particles active. Here's a code snippet to illustrate what I mean:

sim = rebound.Simulation()
sim.add(m=Mearth, hash="Earth", r=RE_eq)   # RE_eq is Earth's radius at the equator

sim.integrator = "WHFAST"
sim.dt = 2e-9

sim.configure_box(10)
sim.collision = "tree"
sim.collision_resolve = collision_save  # collision_save is a function that writes the colliding particle's orbital info to a text file

 # I have a pandas dataframe called "catalog" with x,y,z,vx,vy,vz entries in it for each object, which I add to the sim:
for index in range(len(catalog)):
    entry = catalog.iloc[index]
    sim.add(x=entry.x, y=entry.y, z=entry.z, vx=entry.vx, vy=entry.vy, vz=entry.vz, r=5/aukm, hash='sat{}'.format(int(index+1)))

# note Earth as the only gravitationally active particle:
sim.N_active = 1

an example collision that results is:

collision.p1: 0 (aka Earth)
collision.p2: 4679, with altitude sqrt(x^2+y^2+z^2) - R_earth = -3231.241... km and eccentricity e = 0.00127...
time:  2.3252537185725912e-07 , distance between particles:  3146.7582047072597  km

Thanks,
Sarah Thiele (grad student at Princeton U.)

[hannorein]

Hi Sarah,

I can't answer your question without a minimal code example that I can run myself. Have you tried visualizing your simulation to see how the orbits look like? That might give you some clues as to what is going on.

Not directly related to your question, but note that the tree collision detection method only checks for overlapping particles after each timestep. You might want to use linetree instead. See https://rebound.readthedocs.io/en/latest/collisions/#linetree.

Hanno