From 69deec12f5036660f8b7f1ea0f476b8970c47d73 Mon Sep 17 00:00:00 2001
From: Matthieu Schaller <matthieu.schaller@durham.ac.uk>
Date: Fri, 22 Jun 2018 13:43:49 +0100
Subject: [PATCH] Pre-merge fixes
---
examples/EAGLE_12/eagle_12.yml | 8 ++-
examples/EAGLE_25/eagle_25.yml | 6 +-
.../ZeldovichPancake_3D/zeldovichPancake.yml | 2 +-
src/gravity_cache.h | 6 +-
src/kernel_gravity.h | 2 -
src/mesh_gravity.c | 72 ++-----------------
src/runner_doiact_grav.h | 8 +--
7 files changed, 22 insertions(+), 82 deletions(-)
diff --git a/examples/EAGLE_12/eagle_12.yml b/examples/EAGLE_12/eagle_12.yml
index 1403550036..24c4876f90 100644
--- a/examples/EAGLE_12/eagle_12.yml
+++ b/examples/EAGLE_12/eagle_12.yml
@@ -23,8 +23,9 @@ TimeIntegration:
dt_max: 1e-4 # The maximal time-step size of the simulation (in internal units).
Scheduler:
- max_top_level_cells: 15
-
+ max_top_level_cells: 8
+ cell_split_size: 64
+
# Parameters governing the snapshots
Snapshots:
basename: eagle # Common part of the name of output files
@@ -42,9 +43,10 @@ Statistics:
# Parameters for the self-gravity scheme
Gravity:
eta: 0.025 # Constant dimensionless multiplier for time integration.
- theta: 0.85 # Opening angle (Multipole acceptance criterion)
+ theta: 0.7 # Opening angle (Multipole acceptance criterion)
comoving_softening: 0.0026994 # Comoving softening length (in internal units).
max_physical_softening: 0.0007 # Physical softening length (in internal units).
+ mesh_side_length: 32
# Parameters for the hydrodynamics scheme
SPH:
diff --git a/examples/EAGLE_25/eagle_25.yml b/examples/EAGLE_25/eagle_25.yml
index 462ce37912..95cb53f665 100644
--- a/examples/EAGLE_25/eagle_25.yml
+++ b/examples/EAGLE_25/eagle_25.yml
@@ -23,7 +23,8 @@ TimeIntegration:
dt_max: 1e-4 # The maximal time-step size of the simulation (in internal units).
Scheduler:
- max_top_level_cells: 20
+ max_top_level_cells: 16
+ cell_split_size: 64
# Parameters governing the snapshots
Snapshots:
@@ -41,9 +42,10 @@ Statistics:
# Parameters for the self-gravity scheme
Gravity:
eta: 0.025 # Constant dimensionless multiplier for time integration.
- theta: 0.85 # Opening angle (Multipole acceptance criterion)
+ theta: 0.7 # Opening angle (Multipole acceptance criterion)
comoving_softening: 0.0026994 # Comoving softening length (in internal units).
max_physical_softening: 0.0007 # Physical softening length (in internal units).
+ mesh_side_length: 32
# Parameters for the hydrodynamics scheme
SPH:
diff --git a/examples/ZeldovichPancake_3D/zeldovichPancake.yml b/examples/ZeldovichPancake_3D/zeldovichPancake.yml
index 2319cf995c..d0bbe55222 100644
--- a/examples/ZeldovichPancake_3D/zeldovichPancake.yml
+++ b/examples/ZeldovichPancake_3D/zeldovichPancake.yml
@@ -46,7 +46,7 @@ Cosmology:
Gravity:
mesh_side_length: 16
eta: 0.025
- theta: 0.7
+ theta: 0.85
r_cut_max: 5.
comoving_softening: 0.0001
max_physical_softening: 0.0001
diff --git a/src/gravity_cache.h b/src/gravity_cache.h
index 09ddb3ef88..9101fc763d 100644
--- a/src/gravity_cache.h
+++ b/src/gravity_cache.h
@@ -136,7 +136,7 @@ static INLINE void gravity_cache_init(struct gravity_cache *c,
}
/**
- * @param Zero all the output fields (acceleration and potential) of a
+ * @brief Zero all the output fields (acceleration and potential) of a
* #gravity_cache.
*
* @param c The #gravity_cache to zero.
@@ -328,12 +328,16 @@ gravity_cache_populate_no_mpole(const timebin_t max_active_bin,
* the multi-pole.
*
* @param max_active_bin The largest active bin in the current time-step.
+ * @param periodic Are we using periodic BCs ?
+ * @param dim The size of the simulation volume along each dimension.
* @param c The #gravity_cache to fill.
* @param gparts The #gpart array to read from.
* @param gcount The number of particles to read.
* @param gcount_padded The number of particle to read padded to the next
* multiple of the vector length.
* @param cell The cell we play with (to get reasonable padding positions).
+ * @param CoM The position of the multipole.
+ * @param r_max2 The square of the multipole radius.
* @param grav_props The global gravity properties.
*/
__attribute__((always_inline)) INLINE static void
diff --git a/src/kernel_gravity.h b/src/kernel_gravity.h
index 208d397014..e23e43083d 100644
--- a/src/kernel_gravity.h
+++ b/src/kernel_gravity.h
@@ -26,8 +26,6 @@
#include "inline.h"
#include "minmax.h"
-//#define GADGET2_SOFTENING_CORRECTION
-
#ifdef GADGET2_SOFTENING_CORRECTION
/*! Conversion factor between Plummer softening and internal softening */
#define kernel_gravity_softening_plummer_equivalent 2.8
diff --git a/src/mesh_gravity.c b/src/mesh_gravity.c
index 04a1b51351..5c0016bc51 100644
--- a/src/mesh_gravity.c
+++ b/src/mesh_gravity.c
@@ -118,52 +118,6 @@ __attribute__((always_inline)) INLINE static void CIC_set(
mesh[row_major_id_periodic(i + 1, j + 1, k + 1, N)] += value * dx * dy * dz;
}
-/**
- * @brief Assigns a given multipole to a density mesh using the CIC method.
- *
- * @param m The #multipole.
- * @param rho The density mesh.
- * @param N the size of the mesh along one axis.
- * @param fac The width of a mesh cell.
- * @param dim The dimensions of the simulation box.
- */
-INLINE static void multipole_to_mesh_CIC(const struct gravity_tensors* m,
- double* rho, int N, double fac,
- const double dim[3]) {
- error("aa");
- /* Box wrap the multipole's position */
- const double CoM_x = box_wrap(m->CoM[0], 0., dim[0]);
- const double CoM_y = box_wrap(m->CoM[1], 0., dim[1]);
- const double CoM_z = box_wrap(m->CoM[2], 0., dim[2]);
-
- /* Workout the CIC coefficients */
- int i = (int)(fac * CoM_x);
- if (i >= N) i = N - 1;
- const double dx = fac * CoM_x - i;
- const double tx = 1. - dx;
-
- int j = (int)(fac * CoM_y);
- if (j >= N) j = N - 1;
- const double dy = fac * CoM_y - j;
- const double ty = 1. - dy;
-
- int k = (int)(fac * CoM_z);
- if (k >= N) k = N - 1;
- const double dz = fac * CoM_z - k;
- const double tz = 1. - dz;
-
-#ifdef SWIFT_DEBUG_CHECKS
- if (i < 0 || i >= N) error("Invalid multipole position in x");
- if (j < 0 || j >= N) error("Invalid multipole position in y");
- if (k < 0 || k >= N) error("Invalid multipole position in z");
-#endif
-
- const double mass = m->m_pole.M_000;
-
- /* CIC ! */
- CIC_set(rho, N, i, j, k, tx, ty, tz, dx, dy, dz, mass);
-}
-
/**
* @brief Assigns a given #gpart to a density mesh using the CIC method.
*
@@ -245,9 +199,9 @@ void mesh_to_gparts_CIC(struct gpart* gp, const double* pot, int N, double fac,
const double tz = 1. - dz;
#ifdef SWIFT_DEBUG_CHECKS
- if (i < 0 || i >= N) error("Invalid multipole position in x");
- if (j < 0 || j >= N) error("Invalid multipole position in y");
- if (k < 0 || k >= N) error("Invalid multipole position in z");
+ if (i < 0 || i >= N) error("Invalid gpart position in x");
+ if (j < 0 || j >= N) error("Invalid gpart position in y");
+ if (k < 0 || k >= N) error("Invalid gpart position in z");
#endif
#ifdef SWIFT_GRAVITY_FORCE_CHECKS
@@ -329,10 +283,8 @@ void pm_mesh_compute_potential(struct pm_mesh* mesh, const struct engine* e) {
const struct space* s = e->s;
const double r_s = mesh->r_s;
const double box_size = s->dim[0];
- const int cdim[3] = {s->cdim[0], s->cdim[1], s->cdim[2]};
const double dim[3] = {s->dim[0], s->dim[1], s->dim[2]};
- if (cdim[0] != cdim[1] || cdim[0] != cdim[2]) error("Non-square mesh");
if (r_s <= 0.) error("Invalid value of a_smooth");
/* Some useful constants */
@@ -340,20 +292,6 @@ void pm_mesh_compute_potential(struct pm_mesh* mesh, const struct engine* e) {
const int N_half = N / 2;
const double cell_fac = N / box_size;
-/* Recover the list of top-level multipoles */
-/* const int nr_cells = s->nr_cells; */
-/* struct gravity_tensors* restrict multipoles = s->multipoles_top; */
-
-#ifdef SWIFT_DEBUG_CHECKS
-/* const struct cell* cells = s->cells_top; */
-/* const integertime_t ti_current = e->ti_current; */
-
-/* /\* Make sure everything has been drifted to the current point *\/ */
-/* for (int i = 0; i < nr_cells; ++i) */
-/* if (cells[i].ti_old_multipole != ti_current) */
-/* error("Top-level multipole %d not drifted", i); */
-#endif
-
/* Use the memory allocated for the potential to temporarily store rho */
double* restrict rho = mesh->potential;
if (rho == NULL) error("Error allocating memory for density mesh");
@@ -371,11 +309,9 @@ void pm_mesh_compute_potential(struct pm_mesh* mesh, const struct engine* e) {
fftw_plan inverse_plan = fftw_plan_dft_c2r_3d(
N, N, N, frho, rho, FFTW_ESTIMATE | FFTW_DESTROY_INPUT);
- /* Do a CIC mesh assignment of the multipoles */
- /* for (int i = 0; i < nr_cells; ++i) */
- /* multipole_to_mesh_CIC(&multipoles[i], rho, N, cell_fac, dim); */
bzero(rho, N * N * N * sizeof(double));
+ /* Do a CIC mesh assignment of the gparts */
for (size_t i = 0; i < e->s->nr_gparts; ++i)
gpart_to_mesh_CIC(&e->s->gparts[i], rho, N, cell_fac, dim);
diff --git a/src/runner_doiact_grav.h b/src/runner_doiact_grav.h
index a36983eb43..b8a668bf98 100644
--- a/src/runner_doiact_grav.h
+++ b/src/runner_doiact_grav.h
@@ -1376,7 +1376,7 @@ static INLINE void runner_dopair_recursive_grav(struct runner *r,
} else if (!ci->split && !cj->split) {
/* We have two leaves. Go P-P. */
- runner_dopair_grav_pp(r, ci, cj, /*symmetric*/ 1, /*allow_mpoles*/ 0);
+ runner_dopair_grav_pp(r, ci, cj, /*symmetric*/ 1, /*allow_mpoles*/ 1);
} else {
@@ -1584,10 +1584,8 @@ static INLINE void runner_do_grav_long_range(struct runner *r, struct cell *ci,
/* Call the PM interaction fucntion on the active sub-cells of ci */
- runner_dopair_recursive_grav_pm(r, ci, cj);
- // runner_dopair_grav_mm(r, ci, cj);
- //runner_dopair_grav_pp(r, ci, cj, 0, 1);
-
+ //runner_dopair_recursive_grav_pm(r, ci, cj);
+ runner_dopair_grav_mm(r, ci, cj);
} /* We are in charge of this pair */
} /* Loop over top-level cells */
--
GitLab