diff --git a/src/cell.c b/src/cell.c index f000b824ee134ad02c8dcf41b1ee12acf5ff7914..37cd46b7dee38f55a961e11af56801a7562a6ce3 100644 --- a/src/cell.c +++ b/src/cell.c @@ -1733,7 +1733,8 @@ void cell_activate_subcell_grav_tasks(struct cell *ci, struct cell *cj, const double r2 = dx * dx + dy * dy + dz * dz; /* Can we use multipoles ? */ - if (gravity_multipole_accept(multi_i, multi_j, theta_crit2, r2)) { + if (gravity_multipole_accept(multi_i->r_max, multi_j->r_max, theta_crit2, + r2)) { /* Ok, no need to drift anything */ return; diff --git a/src/engine.c b/src/engine.c index 6f1864918875cec97214ffea60cb3bcfd5c64c46..99c4baee1e6856fa08789b7b81103970ded3d57e 100644 --- a/src/engine.c +++ b/src/engine.c @@ -1798,8 +1798,9 @@ void engine_make_self_gravity_tasks_mapper(void *map_data, int num_elements, const double r2 = dx * dx + dy * dy + dz * dz; /* Are the cells too close for a MM interaction ? */ - if (!gravity_multipole_accept_rebuild(multi_i, multi_j, theta_crit2, - r2)) { + if (!gravity_multipole_accept(multi_i->r_max_rebuild, + multi_j->r_max_rebuild, theta_crit2, + r2)) { /* Ok, we need to add a direct pair calculation */ scheduler_addtask(sched, task_type_pair, task_subtype_grav, 0, 0, diff --git a/src/multipole.h b/src/multipole.h index e37d2ea40566ab2506063d9ef6b6eeb73bdda0c2..b894fefc233991cfb2a80ca55171eccbbd08616a 100644 --- a/src/multipole.h +++ b/src/multipole.h @@ -2390,52 +2390,20 @@ INLINE static void gravity_M2P(const struct multipole *ma, /** * @brief Checks whether a cell-cell interaction can be appromixated by a M-M - * interaction using the CoM and cell radius at rebuild. + * interaction using the distance and cell radius. * * We use the multipole acceptance criterion of Dehnen, 2002, JCoPh, Volume 179, * Issue 1, pp.27-42, equation 10. * - * @param ma The #multipole of the first #cell. - * @param mb The #multipole of the second #cell. - * @param theta_crit2 The square of the critical opening angle. - * @param r2 Square of the distance (periodically wrapped) between the - * multipoles. - */ -__attribute__((always_inline)) INLINE static int -gravity_multipole_accept_rebuild(const struct gravity_tensors *const ma, - const struct gravity_tensors *const mb, - double theta_crit2, double r2) { - - const double r_crit_a = ma->r_max_rebuild; - const double r_crit_b = mb->r_max_rebuild; - const double size = r_crit_a + r_crit_b; - const double size2 = size * size; - - // MATTHIEU: Make this mass-dependent ? - - /* Multipole acceptance criterion (Dehnen 2002, eq.10) */ - return (r2 * theta_crit2 > size2); -} - -/** - * @brief Checks whether a cell-cell interaction can be appromixated by a M-M - * interaction using the CoM and cell radius at the current time. - * - * We use the multipole acceptance criterion of Dehnen, 2002, JCoPh, Volume 179, - * Issue 1, pp.27-42, equation 10. - * - * @param ma The #multipole of the first #cell. - * @param mb The #multipole of the second #cell. + * @param r_crit_a The size of the multipole A. + * @param r_crit_b The size of the multipole B. * @param theta_crit2 The square of the critical opening angle. * @param r2 Square of the distance (periodically wrapped) between the * multipoles. */ __attribute__((always_inline)) INLINE static int gravity_multipole_accept( - const struct gravity_tensors *const ma, - const struct gravity_tensors *const mb, double theta_crit2, double r2) { + double r_crit_a, double r_crit_b, double theta_crit2, double r2) { - const double r_crit_a = ma->r_max; - const double r_crit_b = mb->r_max; const double size = r_crit_a + r_crit_b; const double size2 = size * size; diff --git a/src/runner_doiact_grav.h b/src/runner_doiact_grav.h index bdb85a1fab2c85a17208eeee4c4af78eec114cb5..70f29041282a1f6aa475d7eb1932b962d30978bd 100644 --- a/src/runner_doiact_grav.h +++ b/src/runner_doiact_grav.h @@ -224,6 +224,9 @@ void runner_dopair_grav_pp_full(struct runner *r, struct cell *ci, const float h_inv_i = 1.f / h_i; const float h_inv3_i = h_inv_i * h_inv_i * h_inv_i; + ///* Can we use the multipole in cj ? */ + // if + /* Local accumulators for the acceleration */ float a_x = 0.f, a_y = 0.f, a_z = 0.f; @@ -1117,7 +1120,8 @@ void runner_dopair_grav(struct runner *r, struct cell *ci, struct cell *cj, * option... */ /* Can we use M-M interactions ? */ - if (gravity_multipole_accept(multi_i, multi_j, theta_crit2, r2)) { + if (gravity_multipole_accept(multi_i->r_max, multi_j->r_max, theta_crit2, + r2)) { /* MATTHIEU: make a symmetric M-M interaction function ! */ runner_dopair_grav_mm(r, ci, cj); @@ -1341,7 +1345,8 @@ void runner_do_grav_long_range(struct runner *r, struct cell *ci, int timer) { } /* Check the multipole acceptance criterion */ - if (gravity_multipole_accept(multi_i, multi_j, theta_crit2, r2)) { + if (gravity_multipole_accept(multi_i->r_max, multi_j->r_max, theta_crit2, + r2)) { /* Go for a (non-symmetric) M-M calculation */ runner_dopair_grav_mm(r, ci, cj); @@ -1364,8 +1369,9 @@ void runner_do_grav_long_range(struct runner *r, struct cell *ci, int timer) { const double r2_rebuild = dx * dx + dy * dy + dz * dz; /* Is the criterion violated now but was OK at the last rebuild ? */ - if (gravity_multipole_accept_rebuild(multi_i, multi_j, theta_crit2, - r2_rebuild)) { + if (gravity_multipole_accept(multi_i->r_max_rebuild, + multi_j->r_max_rebuild, theta_crit2, + r2_rebuild)) { /* Alright, we have to take charge of that pair in a different way. */ // MATTHIEU: We should actually open the tree-node here and recurse.