Main REBOUND Functions

group MainRebFunctions

These are the functions that typically need to be called by the user.

Functions

struct reb_simulation *reb_create_simulation(void)

Creates and initialises a REBOUND simulation.

Allocate memory for one reb_simulation structure, initialise all variables and returni the pointer to the reb_simulation sructure. This function must be called before any particles are added.

struct reb_simulation *reb_copy_simulation(struct reb_simulation *r)

Creates a deep copy of a REBOUND simulation.

All simulation data, including all particle data will be copied. Function pointers need to be set manually afterwards. Working on the copy will not affect the original simulation.

int reb_diff_simulations(struct reb_simulation *r1, struct reb_simulation *r2, int output_option)

Compares to REBOUND simulations bit by bit.

Return

If r1 and r2 are exactly equal to each other then 0 is returned, otherwise 1. The walltime parameter in the REBOUND struct is ignored in this comparison.

Parameters
  • r1: The first REBOUND simulation to be compared.

  • r2: The second REBOUND simulation to be compared.

  • output_option: Is set to 1, then the output is printed on the screen. If set to 2, only the return value indicates if the simulations are different.

void reb_init_simulation(struct reb_simulation *r)

Initialize reb_simulation structure.

Same as reb_create_simulation() but does not allocate memory for structure itself.

Parameters
  • r: Structure to be initialized (needs to be allocated externally).

void reb_step(struct reb_simulation *const r)

Performon one integration step.

You rarely want to call this function yourself. Use reb_integrate instead.

Parameters
  • r: The rebound simulation to be integrated by one step.

enum REB_STATUS reb_integrate(struct reb_simulation *const r, double tmax)

Performs the actual integration.

This function performs an integration from the current time t until time tmax.

Return

This function returns an integer, indicating the success of the integration.

Parameters
  • r: The rebound simulation to be integrated.

  • tmax: The time to be integrated to. Set this to INFINITY to integrate forever.

void reb_integrator_synchronize(struct reb_simulation *r)

Synchronize particles manually at end of timestep.

This function should be called if the WHFAST integrator is used, safe_mode is set to zero and an output is needed. This advances the positions and velocities to be synchronized. If enabled, it also applies the symplectic corrector. If safe_mode is enabled, this function has no effect.

Parameters
  • r: The rebound simulation to be synchronized

void reb_integrator_reset(struct reb_simulation *r)

Cleanup all temporarily stored integrator values.

Parameters
  • r: The rebound simulation to be considered

void reb_configure_box(struct reb_simulation *const r, const double boxsize, const int root_nx, const int root_ny, const int root_nz)

Configure the boundary/root box.

This function helps to setup the variables for the simulation box. Call this function when using open, periodic or shear periodic boundary conditions.

Parameters
  • r: The rebound simulation to be considered

  • boxsize: The size of the root box

  • root_nx: The numbe rof root boxes in the x direction.

  • root_ny: The numbe rof root boxes in the y direction.

  • root_nz: The numbe rof root boxes in the z direction.

void reb_free_simulation(struct reb_simulation *const r)

Frees up all space used by a REBOUND simulation and the reb_simulation structure itself.

The REBOUND simulation is not usable anymore after being passed to this function.

Parameters
  • r: The rebound simulation to be freed

void reb_add(struct reb_simulation *const r, struct reb_particle pt)

Adds a particle to the simulation.

This function adds the particle pt to the simulation.

Parameters
  • r: The rebound simulation to which the particle will be added

  • pt: The particle to be added. Note that this is a structure, not a reference to a structure.

void reb_remove_all(struct reb_simulation *const r)

Remove all particles.

Parameters
  • r: The rebound simulation to be considered

int reb_remove(struct reb_simulation *const r, int index, int keepSorted)

Remove a particle by the position in particles array.

Return

Returns 1 if particle was successfully removed, 0 if index passed was out of range.

Parameters
  • r: The rebound simulation to be considered

  • index: The index in the particles array of the particle to be removed.

  • keepSorted: Set to 1, then particles with indices higher than index are all shifted down one position, ensuring the ordering remains.

int reb_remove_by_hash(struct reb_simulation *const r, uint32_t hash, int keepSorted)

Remove a particle by its hash.

see examples/removing_particles_from_simulation.

Return

Returns 1 if particle successfully removed, 0 if hash was not found in the particles array.

Parameters
  • r: The rebound simulation to be considered

  • id: The hash of the particle to be removed.

  • keepSorted: If set to 1 keep the particles with indices in the particles array higher than the one with the passed id are all shifted down one position, ensuring the ordering remains.

struct reb_particle *reb_get_particle_by_hash(struct reb_simulation *const r, uint32_t hash)

Get a pointer to a particle by its hash.

see examples/uniquely_identifying_particles_with_hashes.

Return

A pointer to the particle if found, NULL otherwise.

Parameters
  • r: The rebound simulation to be considered.

  • hash: The hash of the particle to search for.

void reb_run_heartbeat(struct reb_simulation *const r)

Run the heartbeat function and check for escaping/colliding particles.

You rarely want to call this function yourself. It is used internally to call the function you set to the heartbeat variable in reb_simulation.

Parameters
  • r: The rebound simulation to be considered

double reb_integrator_mercurius_L_mercury(const struct reb_simulation *const r, double d, double dcrit)

A force switching function for the MERCURIUS integrator. This function implements the same polynomial switching function used in MERCURY.

double reb_integrator_mercurius_L_infinite(const struct reb_simulation *const r, double d, double dcrit)

A force switching function for the MERCURIUS integrator. This function implements an infinitely differentiable switching function.

int reb_collision_resolve_halt(struct reb_simulation *const r, struct reb_collision c)

Resolve collision by simply halting the integration and setting r->status=REB_EXIT_COLLISION (Default)

int reb_collision_resolve_hardsphere(struct reb_simulation *const r, struct reb_collision c)

Hardsphere collision resolving routine (default).

int reb_collision_resolve_merge(struct reb_simulation *const r, struct reb_collision c)

Merging collision resolving routine.

Merges particle with higher index into particle of lower index. Conserves mass, momentum and volume.