Integrator Structures

group IntegratorStructs

Structures for the various integrators.

Variables related to time, current number of particles and simulation status/control

double t

Current simulation time.

double G

Gravitational constant. Default: 1.

double softening

Gravitational softening parameter. Default: 0.

double dt

Current timestep.

double dt_last_done

Last dt used by integrator.

int N

Current number of particles on this node.

int N_var

Total number of variational particles. Default: 0.

int var_config_N

Number of variational configuration structs. Default: 0.

struct reb_variational_configuration *var_config

These configuration structs contain details on variational particles.

int N_active

Number of massive particles included in force calculation (default: N). Particles with index >= N_active are considered testparticles.

int testparticle_type

Type of the particles with an index>=N_active. 0 means particle does not influence any other particle (default), 1 means particles with index < N_active feel testparticles (similar to MERCURY’s small particles). Testparticles never feel each other.

struct reb_hash_pointer_pair *particle_lookup_table

Array of pairs that map particles’ hashes to their index in the particles array.

int hash_ctr

Counter for number of assigned hashes to assign unique values.

int N_lookup

Number of entries in the particle lookup table.

int allocatedN_lookup

Number of lookup table entries allocated.

int allocatedN

Current maximum space allocated in the particles array on this node.

struct reb_particle *particles

Main particle array. This contains all particles on this node.

struct reb_vec3d *gravity_cs

Vector containing the information for compensated gravity summation.

int gravity_cs_allocatedN

Current number of allocated space for cs array.

struct reb_treecell **tree_root

Pointer to the roots of the trees.

int tree_needs_update

Flag to force a tree update (after boundary check)

double opening_angle2

Square of the cell opening angle \( \theta \).

REB_STATUS status

Set to 1 to exit the simulation at the end of the next timestep.

int exact_finish_time

Set to 1 to finish the integration exactly at tmax. Set to 0 to finish at the next dt. Default is 1.

unsigned int force_is_velocity_dependent

Set to 1 if integrator needs to consider velocity dependent forces.

unsigned int gravity_ignore_terms

Ignore the gravity form the central object (1 for WHFast, 2 for WHFastHelio, 0 otherwise)

double output_timing_last

Time when reb_output_timing() was called the last time.

unsigned long display_clock

Display clock, internal variable for timing refreshs.

int save_messages

Set to 1 to ignore messages (used in python interface).

char **messages

Array of strings containing last messages (only used if save_messages==1).

double exit_max_distance

Exit simulation if distance from origin larger than this value.

double exit_min_distance

Exit simulation if distance from another particle smaller than this value.

double usleep

Wait this number of microseconds after each timestep, useful for slowing down visualization.

struct reb_display_data *display_data
int track_energy_offset

< Datastructure stores visualization related data. Does not have to be modified by the user.

Track energy change during collisions and ejections (default: 0).

double energy_offset

Energy offset due to collisions and ejections (only calculated if track_energy_offset=1).

Variables related to ghost/root boxes

struct reb_vec3d boxsize

Size of the entire box, root_x*boxsize.

double boxsize_max

Maximum size of the entire box in any direction. Set in box_init().

double root_size

Size of a root box.

int root_n

Total number of root boxes in all directions, root_nx*root_ny*root_nz. Default: 1. Set in box_init().

int root_nx

Number of root boxes in x direction. Default: 1.

int root_ny

Number of root boxes in y direction. Default: 1.

int root_nz

Number of root boxes in z direction. Default: 1.

int nghostx

Number of ghostboxes in x direction.

int nghosty

Number of ghostboxes in y direction.

int nghostz

Number of ghostboxes in z direction.

Variables related to collision search and detection

int collision_resolve_keep_sorted

Keep particles sorted if collision_resolve removes particles during a collision.

struct reb_collision *collisions

Array of all collisions.

int collisions_allocatedN

Size allocated for collisions.

double minimum_collision_velocity

Used for hard sphere collision model.

double collisions_plog

Keep track of momentum exchange (used to calculate collisional viscosity in ring systems.

double max_radius[2]

Two largest particle radii, set automatically, needed for collision search.

long collisions_Nlog

Keep track of number of collisions.

Variables related to the chaos indicator MEGNO

int calculate_megno

Internal flag that determines if megno is calculated (default=0, but megno_init() sets it to the index of variational particles used for megno)

double megno_Ys

Running megno sum (internal use)

double megno_Yss

Running megno sum (internal use)

double megno_cov_Yt

covariance of MEGNO Y and t

double megno_var_t

variance of t

double megno_mean_t

mean of t

double megno_mean_Y

mean of MEGNO Y

long megno_n

number of covariance updates

Variables related to SimulationArchive

long simulationarchive_size_first

Size of the initial binary file in a SA.

long simulationarchive_size_snapshot

Size of a snapshot in a SA (other than 1st), in bytes.

double simulationarchive_interval

Current sampling cadence, in code units.

double simulationarchive_interval_walltime

Current sampling cadence, in wall time.

double simulationarchive_next

Next output time.

char *simulationarchive_filename

Name of output file.

double simulationarchive_walltime

Current walltime since beginning of simulation.

struct timeval simulationarchive_time

Time of last output.

Variables describing the current module selection

enum [anonymous]

Available collision routines.

Values:

REB_VISUALIZATION_NONE = 0

No visualization (default if OPENGL compiler flag is turned off)

REB_VISUALIZATION_OPENGL = 1

OpenGL visualization (default if OPENGL compiler flag is turned on)

REB_VISUALIZATION_WEBGL = 2

WebGL visualization, only usable from Jupyter notebook widget.

enum [anonymous]

Available collision routines.

Values:

REB_COLLISION_NONE = 0

Do not search for collisions (default)

REB_COLLISION_DIRECT = 1

Direct collision search O(N^2)

REB_COLLISION_TREE = 2

Tree based collision search O(N log(N))

REB_COLLISION_MERCURIUS = 3

Direct collision search optimized for MERCURIUS.

enum [anonymous]

Available integrators.

Values:

REB_INTEGRATOR_IAS15 = 0

IAS15 integrator, 15th order, non-symplectic (default)

REB_INTEGRATOR_WHFAST = 1

WHFast integrator, symplectic, 2nd order, up to 11th order correctors.

REB_INTEGRATOR_SEI = 2

SEI integrator for shearing sheet simulations, symplectic, needs OMEGA variable.

REB_INTEGRATOR_LEAPFROG = 4

LEAPFROG integrator, simple, 2nd order, symplectic.

REB_INTEGRATOR_HERMES = 5

HERMES Integrator for close encounters (experimental)

REB_INTEGRATOR_NONE = 7

Do not integrate anything.

REB_INTEGRATOR_JANUS = 8

Bit-wise reversible JANUS integrator.

REB_INTEGRATOR_MERCURIUS = 9

MERCURIUS integrator.

enum [anonymous]

Available boundary conditions.

Values:

REB_BOUNDARY_NONE = 0

Do not check for anything (default)

REB_BOUNDARY_OPEN = 1

Open boundary conditions. Removes particles if they leave the box.

REB_BOUNDARY_PERIODIC = 2

Periodic boundary conditions.

REB_BOUNDARY_SHEAR = 3

Shear periodic boundary conditions, needs OMEGA variable.

enum [anonymous]

Available gravity routines.

Values:

REB_GRAVITY_NONE = 0

Do not calculate graviational forces.

REB_GRAVITY_BASIC = 1

Basic O(N^2) direct summation algorithm, choose this for shearing sheet and periodic boundary conditions.

REB_GRAVITY_COMPENSATED = 2

Direct summation algorithm O(N^2) but with compensated summation, slightly slower than BASIC but more accurate.

REB_GRAVITY_TREE = 3

Use the tree to calculate gravity, O(N log(N)), set opening_angle2 to adjust accuracy.

REB_GRAVITY_MERCURIUS = 4

Special gravity routine only for MERCURIUS.

reb_simulation::@1 reb_simulation::visualization

Available collision routines.

reb_simulation::@2 reb_simulation::collision

Available collision routines.

reb_simulation::@3 reb_simulation::integrator

Available integrators.

reb_simulation::@4 reb_simulation::boundary

Available boundary conditions.

reb_simulation::@5 reb_simulation::gravity

Available gravity routines.

Integrator structs (the contain integrator specific variables and temporary data structures)

struct reb_simulation_integrator_sei ri_sei

The SEI struct.

struct reb_simulation_integrator_whfast ri_whfast

The WHFast struct.

struct reb_simulation_integrator_ias15 ri_ias15

The IAS15 struct.

struct reb_simulation_integrator_hermes ri_hermes

The HERMES struct.

struct reb_simulation_integrator_mercurius ri_mercurius

The MERCURIUS struct.

struct reb_simulation_integrator_janus ri_janus

The JANUS struct.

Callback functions

void (*additional_forces)(struct reb_simulation *const r)

This function allows the user to add additional (non-gravitational) forces.

void (*pre_timestep_modifications)(struct reb_simulation *const r)

This function allows the user to make changes before each timestep.

void (*post_timestep_modifications)(struct reb_simulation *const r)

This function allows the user to make changes after each timestep.

void (*heartbeat)(struct reb_simulation *r)

This function is called at the beginning of the simulation and at the end of each timestep.

void (*display_heartbeat)(struct reb_simulation *r)

This function is called at the beginning of the simulation and at the end of each timestep.

double (*coefficient_of_restitution)(const struct reb_simulation *const r, double v)

Return the coefficient of restitution. By default it is NULL, assuming a coefficient of 1.

The velocity of the collision is given to allow for velocity dependent coefficients of restitution.

int (*collision_resolve)(struct reb_simulation *const r, struct reb_collision)

Resolve collision within this function. By default it is NULL, assuming hard sphere model.

A return value of 0 indicates that both particles remain in the simulation. A return value of 1 (2) indicates that particle 1 (2) should be removed from the simulation. A return value of 3 indicates that both particles should be removed from the simulation.

void (*free_particle_ap)(struct reb_particle *p)

Free particle’s ap pointer. Called in reb_remove function.

Hooks for external libraries

void *extras

Pointer to connect additional (optional) libraries, e.g., reboundx.

Enums

enum [anonymous]

Chooses the coordinate system for the WHFast algorithm. Default is Jacobi Coordinates.

Values:

REB_WHFAST_COORDINATES_JACOBI = 0

Jacobi coordinates (default)

REB_WHFAST_COORDINATES_DEMOCRATICHELIOCENTRIC = 1

Democratic Heliocentric coordinates.

REB_WHFAST_COORDINATES_WHDS = 2

WHDS coordinates (Hernandez and Dehnen, 2017)

Variables

double x

x coordinate

double y

y coordinate

double z

z coordinate

double* REBOUND_RESTRICT reb_dp7::p0

0 substep

double* REBOUND_RESTRICT reb_dp7::p1

1 substep

double* REBOUND_RESTRICT reb_dp7::p2

2 substep

double* REBOUND_RESTRICT reb_dp7::p3

3 substep

double* REBOUND_RESTRICT reb_dp7::p4

4 substep

double* REBOUND_RESTRICT reb_dp7::p5

5 substep

double* REBOUND_RESTRICT reb_dp7::p6

6 substep

double shiftx

Relative x position.

double shifty

Relative y position.

double shiftz

Relative z position.

double shiftvx

Relative x velocity.

double shiftvy

Relative y velocity.

double shiftvz

Relative z velocity.

double epsilon

This parameter controls the accuracy of the integrator.

Set to 0 to make IAS15 a non-adaptive integrator. The default value is: 1e-9.

double min_dt

The minimum allowed timestep.

The default value is 0 (no minimal timestep). Set a finite value to this variable if the IAS15 integrator has problems and the timestep becomes excessively small.

unsigned int epsilon_global

Flag that determines how relative acceleration error is estimated.

If set to 1, estimate the fractional error by max(acceleration_error)/max(acceleration), where max is take over all particles. If set to 0, estimate the fractional error by max(acceleration_error/acceleration).

double rcrit

Critical radius in units of Hill radii.

unsigned int recalculate_coordinates_this_timestep

Setting this flag to one will recalculate heliocentric coordinates from the particle structure at the beginning of the next timestep.

After one timestep, the flag gets set back to 0. If you want to change particles after every timestep, you also need to set this flag to 1 before every timestep. Default is 0.

unsigned int recalculate_rhill_this_timestep

Setting this flag to one will recalculate Hill radii at the beginning of the next timestep.

After one timestep, the flag gets set back to 0. If you want to recalculate hill radii at every every timestep, you also need to set this flag to 1 before every timestep. Default is 0.

unsigned int safe_mode

If this flag is set (the default), the integrator will recalculate heliocentric coordinates and synchronize after every timestep, to avoid problems with outputs or particle modifications between timesteps.

Setting it to 0 will result in a speedup, but care must be taken to synchronize and recalculate jacobi coordinates when needed.

unsigned int keep_unsynchronized

Generate inertial coordinates at the end of the integration, but do not change the Jacobi/heliocentric coordinates.

Danger zone! Only use this flag if you are absolutely sure what you are doing. This is intended for simulation which have to be reproducible on a bit by bit basis.

unsigned int is_synchronized

Flag to determine if current particle structure is synchronized.

unsigned int mode

Internal. 0 if WH is operating, 1 if IAS15 is operating.

unsigned int encounterN

Number of particles currently having an encounter.

unsigned int globalN
unsigned int globalNactive
unsigned int allocatedN
unsigned int rhillallocatedN
unsigned int encounterAllocatedN
double m0
double *rhill
double *encounterRhill
unsigned int *encounterIndicies
struct reb_particle *encounterParticles
struct reb_particle* REBOUND_RESTRICT reb_simulation_integrator_mercurius::p_hold
struct reb_simulation *mini

Mini simulation integrated using IAS15. See Silburt et al 2016.

struct reb_simulation *global

Global simulation integrated using WHFast. Only set in mini simulation. See Silburt et al 2016).

double hill_switch_factor

Criteria for switching between IAS15 and WHFast in terms of Hill radii (default: 3.).

double solar_switch_factor

Criteria for switching between IAS15 and WHfast in terms of the central star’s physical radius (default: 15.).

int adaptive_hill_switch_factor

Flag (default: 1) for automatically calculating the appropriate HSF value each iteration.

int mini_active

Flag that is set to 1 by HERMES if the mini simulation is active in this timestep.

double OMEGA

Epicyclic/orbital frequency.

double OMEGAZ

Epicyclic frequency in vertical direction.

unsigned int corrector

This variable turns on/off different symplectic correctors for WHFast. See Rein & Tamayo 2015 and Wisdom 2006 for a discussion.

  • 0 (default): turns off all correctors
  • 3: uses third order (two-stage) corrector
  • 5: uses fifth order (four-stage) corrector
  • 7: uses seventh order (six-stage) corrector
  • 11: uses eleventh order (ten-stage) corrector

reb_simulation_integrator_whfast::@0 reb_simulation_integrator_whfast::coordinates

Chooses the coordinate system for the WHFast algorithm. Default is Jacobi Coordinates.

unsigned int recalculate_coordinates_this_timestep

Setting this flag to one will recalculate Jacobi/heliocentric coordinates from the particle structure in the next timestep.

After the timestep, the flag gets set back to 0. If you want to change particles after every timestep, you also need to set this flag to 1 before every timestep. Default is 0.

unsigned int safe_mode

If this flag is set (the default), whfast will recalculate jacobi/heliocentric coordinates and synchronize every timestep, to avoid problems with outputs or particle modifications between timesteps.

Setting it to 0 will result in a speedup, but care must be taken to synchronize and recalculate jacobi coordinates when needed. See AdvWHFast.ipynb in the python_tutorials folder (navigate to it on github if you don’t have ipython notebook installed). The explanation is general, and the python and C flags have the same names.

struct reb_particle* REBOUND_RESTRICT reb_simulation_integrator_whfast::p_jh

Jacobi/heliocentric coordinates.

This array contains the Jacobi/heliocentric coordinates of all particles. It is automatically filled and updated by WHfast. Access this array with caution.

unsigned int keep_unsynchronized

Generate inertial coordinates at the end of the integration, but do not change the Jacobi/heliocentric coordinates.

Danger zone! Only use this flag if you are absolutely sure what you are doing. This is intended for simulation which have to be reproducible on a bit by bit basis.

double scale_pos

Scale of the problem. Positions get divided by this number before the conversion to an integer.

double scale_vel

Scale of the problem. Velocities get divided by this number before the conversion to an integer.

unsigned int order

Order of the scheme. Default is 6.

unsigned int recalculate_integer_coordinates_this_timestep

If this flag is set, then janus will recalculate integer coordinates at the next timestep.

int p1

One of the colliding particles.

int p2

One of the colliding particles.

struct reb_ghostbox gb

Ghostbox (of particle p1, used for periodic and shearing sheet boundary conditions)

int ri

Index of rootcell (needed for MPI only).

struct reb_simulation *sim

Reference to the simulation.

int order

Order of the variational equation. 1 or 2.

int index

Index of the first variational particle in the particles array.

int testparticle

Is this variational configuration describe a test particle? -1 if not.

int index_1st_order_a

Used for 2nd order variational particles only: Index of the first first order variational particle in the particles array.

int index_1st_order_b

Used for 2nd order variational particles only: Index of the first first order variational particle in the particles array.

double x

x-position of the particle.

double y

y-position of the particle.

double z

z-position of the particle.

double vx

x-velocity of the particle.

double vy

y-velocity of the particle.

double vz

z-velocity of the particle.

double ax

x-acceleration of the particle.

double ay

y-acceleration of the particle.

double az

z-acceleration of the particle.

double m

Mass of the particle.

double r

Radius of the particle.

double lastcollision

Last time the particle had a physical collision.

struct reb_treecell *c

Pointer to the cell the particle is currently in.

uint32_t hash

hash to identify particle.

void *ap

Functionality for externally adding additional properties to particles.

struct reb_simulation *sim

Pointer to the parent simulation.

double d

Radial distance from central object.

double v

velocity relative to central object’s velocity

double h

Angular momentum.

double P

Orbital period.

double n

Mean motion.

double a

Semi-major axis.

double e

Eccentricity.

double inc

Inclination.

double Omega

Longitude of ascending node.

double omega

Argument of pericenter.

double pomega

Longitude of pericenter.

double f

True anomaly.

double M

Mean anomaly.

double l

Mean Longitude.

double theta

True Longitude.

double T

Time of pericenter passage.

struct reb_simulation_integrator_ias15
#include <rebound.h>

This structure contains variables and pointer used by the IAS15 integrator.

Public Members

double epsilon

This parameter controls the accuracy of the integrator.

Set to 0 to make IAS15 a non-adaptive integrator. The default value is: 1e-9.

double min_dt

The minimum allowed timestep.

The default value is 0 (no minimal timestep). Set a finite value to this variable if the IAS15 integrator has problems and the timestep becomes excessively small.

unsigned int epsilon_global

Flag that determines how relative acceleration error is estimated.

If set to 1, estimate the fractional error by max(acceleration_error)/max(acceleration), where max is take over all particles. If set to 0, estimate the fractional error by max(acceleration_error/acceleration).

struct reb_simulation_integrator_mercurius
#include <rebound.h>

This structure contains variables and pointer used by the MERCURIUS integrator.

Public Members

double rcrit

Critical radius in units of Hill radii.

unsigned int recalculate_coordinates_this_timestep

Setting this flag to one will recalculate heliocentric coordinates from the particle structure at the beginning of the next timestep.

After one timestep, the flag gets set back to 0. If you want to change particles after every timestep, you also need to set this flag to 1 before every timestep. Default is 0.

unsigned int recalculate_rhill_this_timestep

Setting this flag to one will recalculate Hill radii at the beginning of the next timestep.

After one timestep, the flag gets set back to 0. If you want to recalculate hill radii at every every timestep, you also need to set this flag to 1 before every timestep. Default is 0.

unsigned int safe_mode

If this flag is set (the default), the integrator will recalculate heliocentric coordinates and synchronize after every timestep, to avoid problems with outputs or particle modifications between timesteps.

Setting it to 0 will result in a speedup, but care must be taken to synchronize and recalculate jacobi coordinates when needed.

unsigned int keep_unsynchronized

Generate inertial coordinates at the end of the integration, but do not change the Jacobi/heliocentric coordinates.

Danger zone! Only use this flag if you are absolutely sure what you are doing. This is intended for simulation which have to be reproducible on a bit by bit basis.

unsigned int is_synchronized

Flag to determine if current particle structure is synchronized.

unsigned int mode

Internal. 0 if WH is operating, 1 if IAS15 is operating.

unsigned int encounterN

Number of particles currently having an encounter.

struct reb_simulation_integrator_hermes
#include <rebound.h>

This structure contains variables and pointer used by the HERMES integrator.

Public Members

struct reb_simulation *mini

Mini simulation integrated using IAS15. See Silburt et al 2016.

struct reb_simulation *global

Global simulation integrated using WHFast. Only set in mini simulation. See Silburt et al 2016).

double hill_switch_factor

Criteria for switching between IAS15 and WHFast in terms of Hill radii (default: 3.).

double solar_switch_factor

Criteria for switching between IAS15 and WHfast in terms of the central star’s physical radius (default: 15.).

int adaptive_hill_switch_factor

Flag (default: 1) for automatically calculating the appropriate HSF value each iteration.

int mini_active

Flag that is set to 1 by HERMES if the mini simulation is active in this timestep.

struct reb_simulation_integrator_sei
#include <rebound.h>

This structure contains variables used by the SEI integrator.

This is where the user sets the orbital frequency OMEGA for shearing sheet simulations.

Public Members

double OMEGA

Epicyclic/orbital frequency.

double OMEGAZ

Epicyclic frequency in vertical direction.

struct reb_simulation_integrator_whfast
#include <rebound.h>

This structure contains variables used by the WHFast integrator.

Public Types

enum [anonymous]

Chooses the coordinate system for the WHFast algorithm. Default is Jacobi Coordinates.

Values:

REB_WHFAST_COORDINATES_JACOBI = 0

Jacobi coordinates (default)

REB_WHFAST_COORDINATES_DEMOCRATICHELIOCENTRIC = 1

Democratic Heliocentric coordinates.

REB_WHFAST_COORDINATES_WHDS = 2

WHDS coordinates (Hernandez and Dehnen, 2017)

Public Members

unsigned int corrector

This variable turns on/off different symplectic correctors for WHFast. See Rein & Tamayo 2015 and Wisdom 2006 for a discussion.

  • 0 (default): turns off all correctors
  • 3: uses third order (two-stage) corrector
  • 5: uses fifth order (four-stage) corrector
  • 7: uses seventh order (six-stage) corrector
  • 11: uses eleventh order (ten-stage) corrector

reb_simulation_integrator_whfast::@0 reb_simulation_integrator_whfast::coordinates

Chooses the coordinate system for the WHFast algorithm. Default is Jacobi Coordinates.

unsigned int recalculate_coordinates_this_timestep

Setting this flag to one will recalculate Jacobi/heliocentric coordinates from the particle structure in the next timestep.

After the timestep, the flag gets set back to 0. If you want to change particles after every timestep, you also need to set this flag to 1 before every timestep. Default is 0.

unsigned int safe_mode

If this flag is set (the default), whfast will recalculate jacobi/heliocentric coordinates and synchronize every timestep, to avoid problems with outputs or particle modifications between timesteps.

Setting it to 0 will result in a speedup, but care must be taken to synchronize and recalculate jacobi coordinates when needed. See AdvWHFast.ipynb in the python_tutorials folder (navigate to it on github if you don’t have ipython notebook installed). The explanation is general, and the python and C flags have the same names.

struct reb_particle* REBOUND_RESTRICT reb_simulation_integrator_whfast::p_jh

Jacobi/heliocentric coordinates.

This array contains the Jacobi/heliocentric coordinates of all particles. It is automatically filled and updated by WHfast. Access this array with caution.

unsigned int keep_unsynchronized

Generate inertial coordinates at the end of the integration, but do not change the Jacobi/heliocentric coordinates.

Danger zone! Only use this flag if you are absolutely sure what you are doing. This is intended for simulation which have to be reproducible on a bit by bit basis.

struct reb_simulation_integrator_janus

Public Members

double scale_pos

Scale of the problem. Positions get divided by this number before the conversion to an integer.

double scale_vel

Scale of the problem. Velocities get divided by this number before the conversion to an integer.

unsigned int order

Order of the scheme. Default is 6.

unsigned int recalculate_integer_coordinates_this_timestep

If this flag is set, then janus will recalculate integer coordinates at the next timestep.