I'm using Rebound to simulate 1000 semi-active particles, and I've disabled collision detection by setting it to None. I'm wondering if I can assume that the particles at specific indices in the orbits file (that is generated using the function reb_output_orbits) don't mix with each other. This is considering:
The simulation runs for 10 days, with each day cut into 24-hour segments. I need to restart the simulation from the archive file at the start of each new day.
The initial conditions are generated using the Python version of Rebound and stored in a file named "initial_conditions.txt."
This is because I would like to track the final semi-major axis of each particle.
I would appreciate your help. Here is my code:
Hello,
I'm using Rebound to simulate 1000 semi-active particles, and I've disabled collision detection by setting it to None. I'm wondering if I can assume that the particles at specific indices in the orbits file (that is generated using the function reb_output_orbits) don't mix with each other. This is considering:
This is because I would like to track the final semi-major axis of each particle. I would appreciate your help. Here is my code:
#########################################
include
include
include
include
include
include "rebound.h"
double rand_range(double min, double max); void heartbeat(struct reb_simulation r); double E0; double ang_m0[5]; double ang_m[5]; double com_rel[6]; double L_x; double L_y; double L_z; double com_tot_vmag; double com_vmag; struct reb_particle com0; const int N_planetesimals = 1000; void calculate_com(struct reb_simulation const r, struct reb_particle com0, double com_rel[]); void calculate_angular_momentum(struct reb_simulation* const r, double L_tot[]);
int main(int argc, char* argv[]){
}
void heartbeat(struct reb_simulation* r){ if (reb_output_check(r, 10000)){ // every ~ 1600 years //relative energy error
}
void calculate_angular_momentum(struct reb_simulation const r, double L_tot[]) { //reb_move_to_com(r); L_x = 0; L_y = 0; L_z = 0; struct reb_particle p; for (int i = 0; i < r->N; i++) { p = r->particles[i]; L_x += p.m(p.yp.vz - p.zp.vy); L_y += -p.m(p.xp.vz - p.zp.vx); L_z += p.m(p.xp.vy - p.yp.vx); } double L_t = sqrt(L_zL_z + L_xL_x + L_y*L_y); L_tot[0] = L_x; L_tot[1] = L_y; L_tot[2] = L_z; L_tot[3] = L_t; }
void calculate_com(struct reb_simulation* const r, struct reb_particle com0, double com_rel[]) {
}