From a264c94c998d7da8b6b0de6d55b6ea4bf1df85b2 Mon Sep 17 00:00:00 2001
From: Matthieu Schaller <matthieu.schaller@durham.ac.uk>
Date: Wed, 16 May 2018 15:38:52 +0100
Subject: [PATCH] Correctly define the Ewald factor also when reading the table
from a file.
---
examples/plot_gravity_checks.py | 5 ++++-
src/gravity.c | 8 ++++----
2 files changed, 8 insertions(+), 5 deletions(-)
diff --git a/examples/plot_gravity_checks.py b/examples/plot_gravity_checks.py
index ac6d1699b7..ca4e17c7c8 100644
--- a/examples/plot_gravity_checks.py
+++ b/examples/plot_gravity_checks.py
@@ -81,7 +81,10 @@ plt.figure()
count = 0
# Get the Gadget-2 data if existing
-gadget2_file_list = glob.glob("forcetest_gadget2.txt")
+if periodic:
+ gadget2_file_list = glob.glob("forcetest_gadget2_periodic.txt")
+else:
+ gadget2_file_list = glob.glob("forcetest_gadget2.txt")
if len(gadget2_file_list) != 0:
gadget2_data = np.loadtxt(gadget2_file_list[0])
diff --git a/src/gravity.c b/src/gravity.c
index 42cf1dc4bb..b21c104d06 100644
--- a/src/gravity.c
+++ b/src/gravity.c
@@ -145,9 +145,6 @@ void gravity_exact_force_ewald_init(double boxSize) {
const float factor_cos = 2.f * M_PI;
const float factor_pot = M_PI / alpha2;
- /* Ewald factor to access the table */
- ewald_fac = (double)(2 * Newald) / boxSize;
-
/* Zero everything */
bzero(fewald_x, (Newald + 1) * (Newald + 1) * (Newald + 1) * sizeof(float));
bzero(fewald_y, (Newald + 1) * (Newald + 1) * (Newald + 1) * sizeof(float));
@@ -295,6 +292,9 @@ void gravity_exact_force_ewald_init(double boxSize) {
clocks_from_ticks(getticks() - tic), clocks_getunit());
}
+ /* Ewald factor to access the table */
+ ewald_fac = (double)(2 * Newald) / boxSize;
+
/* Apply the box-size correction */
for (int i = 0; i <= Newald; ++i) {
for (int j = 0; j <= Newald; ++j) {
@@ -508,7 +508,7 @@ void gravity_exact_force_compute_mapper(void *map_data, int nr_gparts,
/* Interact it with all other particles in the space.*/
for (int j = 0; j < (int)s->nr_gparts; ++j) {
- struct gpart *gpj = &s->gparts[j];
+ const struct gpart *gpj = &s->gparts[j];
/* No self interaction */
if (gpi == gpj) continue;
--
GitLab