dtamayo
commented
6 years ago Hola Fran,
The form of the forces (e.g. involving reduced masses) seems like the sort of thing you'd want when considering the effect on a 2 body orbit, in which case you also have to worry about indirect terms (which I think correspond to your second terms).
It's much easier in REBOUND/REBOUNDx. The interface with the user (and REBOUNDx) is always in an inertial frame, so if you have a force between two particles like tides, you just compute the acceleration on one using the actual (not reduced) masses, and apply the back reaction acceleration such that the forces are opposite and equal. You never have to worry about getting the right form for the indirect terms.
I think I've implemented some of the terms in the full tidal effect in tides_precession and synchronous_eccentricity_damping. They correspond to the conservative force in the radial direction and the force opposite the radial velocity, respectively. You might take a look at those. Let me know if you have any questions.
On Dec 31, 2017, at 3:54 PM, fran pozuelos romero notifications@github.com wrote:
Hi all,
I have resumed and old project about tides forces. I coded the equations for a reboundx module but I have a doubt in the last part, after compute the tides forces. The point is that the final particle acceleration is given by:
Where F^{T} are the tides forces and F^{GR} and F^{R} are not considered in this case. So, please ignore them. As one can notice, the main problem is that every planet acceleration depends on the forces of the others planets in the problem. So, when I coded the module (following the examples and the guide) I have something like that:
include
include
include
include "reboundx.h"
static void rebx_calculate_tides_forces(struct reb_simulation const sim, const double K2s, const double K2sDt, const double Ss_z, const double Rs, const int source_index){ struct reb_particle const particles = sim->particles; struct reb_particle* const source = &particles[source_index]; const double G = sim->G;
const int _N_real = sim->N - sim->N_var; for (int i=0;i<_N_real;i++){ const double* K2p = rebx_get_param_check(&particles[i], "K2p", REBX_TYPE_DOUBLE); if(K2p == NULL) continue; // only particles with K2p set feel tidal forces const double* K2pDt = rebx_get_param_check(&particles[i], "K2pDt", REBX_TYPE_DOUBLE); if(K2pDt == NULL) continue; // const double* Sp_z = rebx_get_param_check(&particles[i], "Sp_z", REBX_TYPE_DOUBLE); if(Sp_z == NULL) continue; // const double* Rp = rebx_get_param_check(&particles[i], "Rp", REBX_TYPE_DOUBLE); if(Rp == NULL) continue; // const struct reb_particle p = particles[i];Then, I computed the forces for every i-planet inside the for loop. I don't copy here all the equations involved to avoid the mesh. At the end, we have:
double F_x = (F_tr+(Ptop+Ptos)(v_dot_er/rps)) er_x + Ptop(Sp_cross_er_x - (1./rps_2)triple_product_x) + Ptos(Ss_cross_er_x-(1./rps_2)triple_product_x); double F_y = (F_tr+(Ptop+Ptos)(v_dot_er/rps)) er_y + Ptop(Sp_cross_er_y - (1./rps_2)triple_product_y) + Ptos(Ss_cross_er_y-(1./rps_2)triple_product_y); double F_z = (F_tr+(Ptop+Ptos)(v_dot_er/rps)) er_z + Ptop(Sp_cross_er_z - (1./rps_2)triple_product_z) + Ptos(Ss_cross_er_z-(1./rps_2)triple_product_z);
// Compute the accelaration
F_x = reduced_mass; F_y = reduced_mass; F_z *= reduced_mass;
particles[i].ax += F_x; particles[i].ay += F_y; particles[i].az += F_z;}
}
void rebx_tides_forces(struct reb_simulation const sim, struct rebx_effect const tides_forces){ double K2s = rebx_get_param_check(tides_forces, "K2s", REBX_TYPE_DOUBLE); if (K2s == NULL){ reb_error(sim, "Need to set potential love number of the star (K2s) in tides_forces effect.\n"); } double K2sDt = rebx_get_param_check(tides_forces, "K2sDt", REBX_TYPE_DOUBLE); if (K2sDt == NULL){ reb_error(sim, "Need to set potential love number times constant time lag of the star (K2sDt) in tides_forces effect.\n"); } double Ss_z = rebx_get_param_check(tides_forces, "Ss_z", REBX_TYPE_DOUBLE); if (Ss_z == NULL){ reb_error(sim, "Need to set the Spin in Z direction of the star (Ss_z) in tides_forces effect.\n"); } double Rs = rebx_get_param_check(tides_forces, "Rs", REBX_TYPE_DOUBLE); if (Rs == NULL){ reb_error(sim, "Need to set the radio of the star (Rs) in tides_forces effect.\n"); } const int N_real = sim->N - sim->N_var; struct reb_particle const particles = sim->particles; int source_found=0; for (int i=0; i<N_real; i++){ if (rebx_get_param_check(&particles[i], "Star", REBX_TYPE_INT) != NULL){ source_found = 1; rebx_calculate_tides_forces(sim, K2s, K2sDt, Ss_z, Rs, i); } } if (!source_found){ rebx_calculate_tides_forces(sim, K2s, K2sDt, Ss_z, *Rs, 0); // default source to index 0 if "radiation_source" not found on any particle } }
However, with this formulation, I computed for i-planet the tidal forces and its acceleration, but without taken into account the forces corresponding to others planets. So, basically I'm computing only the first term in the right hand of the equation showed at the beginning, and it is missing the second term.
Any suggestion about how to include this missing part?
Sorry in advance if it's not very clear the problem I have. I will be glad to re-explain if necessary.
Thank you very much =)
— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub, or mute the thread.
franpoz
Hi all,
I have resumed and old project about tides forces. I coded the equations for a reboundx module but I have a doubt in the last part, after compute the tides forces. The point is that the final particle acceleration is given by:
Where F^{T} are the tides forces and F^{GR} and F^{R} are not considered in this case. So, please ignore them. As one can notice, the main problem is that every planet acceleration depends on the forces of the others planets in the problem. So, when I coded the module (following the examples and the guide) I have something like that:
Then, I computed the forces for every i-planet inside the for loop. I don't copy here all the equations involved to avoid the mesh. At the end, we have:
However, with this formulation, I computed for i-planet the tidal forces and its acceleration, but without taken into account the forces corresponding to others planets. So, basically I'm computing only the first term in the right hand of the equation showed at the beginning, and it is missing the second term.
Any suggestion about how to include this missing part?
Sorry in advance if it's not very clear the problem I have. I will be glad to re-explain if necessary.
Thank you very much =)