diff --git a/src/gravity/Default/gravity.h b/src/gravity/Default/gravity.h index e658a6d6a21e9474a1e930deacd6c8c49d6e7767..5bc516735d679f72206a8eb37e48626d3aa5a9de 100644 --- a/src/gravity/Default/gravity.h +++ b/src/gravity/Default/gravity.h @@ -20,19 +20,20 @@ #include <float.h> #include "potentials.h" /** - * @brief Computes the gravity time-step of a given particle - * - * @param gp Pointer to the g-particle data + * @brief Computes the gravity time-step of a given particle. * + * @param phys_cont The physical constants in internal units. + * @param gp Pointer to the g-particle data. */ +__attribute__((always_inline)) INLINE static float gravity_compute_timestep( + const struct phys_const* const phys_const, const struct gpart* const gp) { -__attribute__((always_inline)) -INLINE static float gravity_compute_timestep(const struct phys_const* phys_const, struct gpart* gp) { float dt = FLT_MAX; #ifdef EXTERNAL_POTENTIAL_POINTMASS dt = fminf(dt, external_gravity_pointmass_timestep(phys_const, gp)); #endif + return dt; } @@ -64,10 +65,27 @@ __attribute__((always_inline)) gp->a_grav[2] = 0.f; } -__attribute__((always_inline)) INLINE static void external_gravity(const struct phys_const* phys_const, struct gpart *g) -{ +/** + * @brief Finishes the gravity calculation. + * + * Multiplies the forces and accelerations by the appropiate constants + * + * @param gp The particle to act upon + */ +__attribute__((always_inline)) + INLINE static void gravity_end_force(struct gpart* gp) {} + +/** + * @brief Computes the gravitational acceleration induced by external potentials + * + * @param phys_const The physical constants in internal units. + * @param gp The particle to act upon. + */ +__attribute__((always_inline)) INLINE static void external_gravity( + const struct phys_const* const phys_const, struct gpart* gp) { + #ifdef EXTERNAL_POTENTIAL_POINTMASS - external_gravity_pointmass(phys_const, g); + external_gravity_pointmass(phys_const, gp); #endif } @@ -80,6 +98,3 @@ __attribute__((always_inline)) INLINE static void external_gravity(const struct */ __attribute__((always_inline)) INLINE static void gravity_kick_extra( struct gpart* gp, float dt, float half_dt) {} - -__attribute__((always_inline)) INLINE static void gravity_end_force( - struct gpart* gp) {} diff --git a/src/potentials.h b/src/potentials.h index 2af6cedbd2d62ea4b61d0dcadec295384fe8bdf3..7d34d1b14b46cb746420841625be7a63a99d9e36 100644 --- a/src/potentials.h +++ b/src/potentials.h @@ -15,18 +15,26 @@ struct external_potential { /* Properties of Point Mass */ #ifdef EXTERNAL_POTENTIAL_POINTMASS + #define External_Potential_X (50000 * PARSEC_IN_CGS / const_unit_length_in_cgs) #define External_Potential_Y (50000 * PARSEC_IN_CGS / const_unit_length_in_cgs) #define External_Potential_Z (50000 * PARSEC_IN_CGS / const_unit_length_in_cgs) #define External_Potential_Mass \ (1e10 * SOLAR_MASS_IN_CGS / const_unit_mass_in_cgs) + +/** + * @brief Computes the time-step due to the acceleration from a point mass + * + * @param phys_cont The physical constants in internal units. + * @param gp Pointer to the g-particle data. + */ __attribute__((always_inline)) INLINE static float external_gravity_pointmass_timestep( - const struct phys_const* phys_const, struct gpart* g) { + const struct phys_const* const phys_const, + const struct gpart* const g) { const double G_newton = phys_const->newton_gravity; - /* Currently no limit is imposed */ const float dx = g->x[0] - External_Potential_X; const float dy = g->x[1] - External_Potential_Y; const float dz = g->x[2] - External_Potential_Z; @@ -47,8 +55,17 @@ __attribute__((always_inline)) return 0.03f * sqrtf(a_2 / dota_2); } + +/** + * @brief Computes the gravitational acceleration of a particle due to a point + * mass + * + * @param phys_cont The physical constants in internal units. + * @param gp Pointer to the g-particle data. + */ __attribute__((always_inline)) INLINE static void external_gravity_pointmass( - const struct phys_const* phys_const, struct gpart* g) { + const struct phys_const* const phys_const, struct gpart* g) { + const double G_newton = phys_const->newton_gravity; const float dx = g->x[0] - External_Potential_X; const float dy = g->x[1] - External_Potential_Y; @@ -60,6 +77,14 @@ __attribute__((always_inline)) INLINE static void external_gravity_pointmass( g->a_grav[2] += -G_newton * External_Potential_Mass * dz * rinv * rinv * rinv; } #endif /* EXTERNAL_POTENTIAL_POINTMASS */ + +/** + * @brief Initialises the external potential properties in the internal system + * of units. + * + * @param us The current internal system of units + * @param potential The external potential properties to initialize + */ void initPotentialProperties(struct UnitSystem* us, struct external_potential* potential);