diff --git a/examples/main.c b/examples/main.c index 38847543bd25d6ce877b713b7e6e9ff614eb16bb..a15b62b06c13cb0f2aa0ca05345eb9a9ff670ede 100644 --- a/examples/main.c +++ b/examples/main.c @@ -777,19 +777,6 @@ int main(int argc, char *argv[]) { /* Verify that the fields to dump actually exist */ if (myrank == 0) io_check_output_fields(params, N_total); - /* Initialise the long-range gravity mesh */ - if (with_self_gravity && periodic) { -#ifdef HAVE_FFTW - pm_mesh_init(&mesh, &gravity_properties, dim); -#else - /* Need the FFTW library if periodic and self gravity. */ - error( - "No FFTW library found. Cannot compute periodic long-range forces."); -#endif - } else { - pm_mesh_init_no_mesh(&mesh, dim); - } - /* Initialize the space with these data. */ if (myrank == 0) clocks_gettime(&tic); space_init(&s, params, &cosmo, dim, parts, gparts, sparts, Ngas, Ngpart, @@ -803,6 +790,19 @@ int main(int argc, char *argv[]) { fflush(stdout); } + /* Initialise the long-range gravity mesh */ + if (with_self_gravity && periodic) { +#ifdef HAVE_FFTW + pm_mesh_init(&mesh, &gravity_properties, s.dim); +#else + /* Need the FFTW library if periodic and self gravity. */ + error( + "No FFTW library found. Cannot compute periodic long-range forces."); +#endif + } else { + pm_mesh_init_no_mesh(&mesh, s.dim); + } + /* Check that the matter content matches the cosmology given in the * parameter file. */ if (with_cosmology && with_self_gravity && !dry_run) diff --git a/src/mesh_gravity.c b/src/mesh_gravity.c index d49b850f98c3d011b4bc1f4a7462be31954e1e7f..4d78788a697b74b5e28061bc48cb03485caa89da 100644 --- a/src/mesh_gravity.c +++ b/src/mesh_gravity.c @@ -291,6 +291,8 @@ void pm_mesh_compute_potential(struct pm_mesh* mesh, const struct space* s, const double dim[3] = {s->dim[0], s->dim[1], s->dim[2]}; if (r_s <= 0.) error("Invalid value of a_smooth"); + if (mesh->dim[0] != dim[0] || mesh->dim[1] != dim[1] || mesh->dim[2] != dim[2]) + error("Domain size does not match the value stored in the space."); /* Some useful constants */ const int N = mesh->N;