From 6c0d5eb89902e70ce21b14a111ca55baccad08ce Mon Sep 17 00:00:00 2001
From: Matthieu Schaller <schaller@strw.leidenuniv.nl>
Date: Sat, 21 Mar 2020 23:16:57 +0100
Subject: [PATCH] Make the gravity tasks only update the particles atomically

---
 src/gravity/MultiSoftening/gravity.h |  3 ++-
 src/gravity_cache.h                  |  6 +++---
 src/mesh_gravity.c                   |  6 +++---
 src/multipole.h                      |  6 +++---
 src/runner_doiact_grav.c             | 20 +++++++++++---------
 5 files changed, 22 insertions(+), 19 deletions(-)

diff --git a/src/gravity/MultiSoftening/gravity.h b/src/gravity/MultiSoftening/gravity.h
index 07ed0e78b1..7cb2c50434 100644
--- a/src/gravity/MultiSoftening/gravity.h
+++ b/src/gravity/MultiSoftening/gravity.h
@@ -23,6 +23,7 @@
 #include <float.h>
 
 /* Local includes. */
+#include "atomic.h"
 #include "cosmology.h"
 #include "error.h"
 #include "gravity_properties.h"
@@ -77,7 +78,7 @@ __attribute__((always_inline)) INLINE static float gravity_get_softening(
 __attribute__((always_inline)) INLINE static void
 gravity_add_comoving_potential(struct gpart* restrict gp, float pot) {
 
-  gp->potential += pot;
+  atomic_add_f(&gp->potential, pot);
 }
 
 /**
diff --git a/src/gravity_cache.h b/src/gravity_cache.h
index 912961d2e0..6e55e3f443 100644
--- a/src/gravity_cache.h
+++ b/src/gravity_cache.h
@@ -491,9 +491,9 @@ __attribute__((always_inline)) INLINE static void gravity_cache_write_back(
   /* Write stuff back to the particles */
   for (int i = 0; i < gcount; ++i) {
     if (active[i]) {
-      gparts[i].a_grav[0] += a_x[i];
-      gparts[i].a_grav[1] += a_y[i];
-      gparts[i].a_grav[2] += a_z[i];
+      atomic_add_f(&gparts[i].a_grav[0], a_x[i]);
+      atomic_add_f(&gparts[i].a_grav[1], a_y[i]);
+      atomic_add_f(&gparts[i].a_grav[2], a_z[i]);
       gravity_add_comoving_potential(&gparts[i], pot[i]);
     }
   }
diff --git a/src/mesh_gravity.c b/src/mesh_gravity.c
index 0efec0bc80..7cabd7bbae 100644
--- a/src/mesh_gravity.c
+++ b/src/mesh_gravity.c
@@ -335,9 +335,9 @@ void mesh_to_gparts_CIC(struct gpart* gp, const double* pot, const int N,
 
   /* Store things back */
   gravity_add_comoving_potential(gp, p);
-  gp->a_grav[0] += fac * a[0];
-  gp->a_grav[1] += fac * a[1];
-  gp->a_grav[2] += fac * a[2];
+  atomic_add_f(&gp->a_grav[0], fac * a[0]);
+  atomic_add_f(&gp->a_grav[1], fac * a[1]);
+  atomic_add_f(&gp->a_grav[2], fac * a[2]);
 #ifdef SWIFT_GRAVITY_FORCE_CHECKS
   gp->potential_PM = p;
   gp->a_grav_PM[0] = fac * a[0];
diff --git a/src/multipole.h b/src/multipole.h
index e03e302c41..178c8088a6 100644
--- a/src/multipole.h
+++ b/src/multipole.h
@@ -2589,9 +2589,9 @@ INLINE static void gravity_L2P(const struct grav_tensor *lb,
 #endif
 
   /* Update the particle */
-  gp->a_grav[0] += a_grav[0];
-  gp->a_grav[1] += a_grav[1];
-  gp->a_grav[2] += a_grav[2];
+  atomic_add_f(&gp->a_grav[0], a_grav[0]);
+  atomic_add_f(&gp->a_grav[1], a_grav[1]);
+  atomic_add_f(&gp->a_grav[2], a_grav[2]);
   gravity_add_comoving_potential(gp, pot);
 
 #ifdef SWIFT_GRAVITY_FORCE_CHECKS
diff --git a/src/runner_doiact_grav.c b/src/runner_doiact_grav.c
index 043b662b10..0d9b2c9d67 100644
--- a/src/runner_doiact_grav.c
+++ b/src/runner_doiact_grav.c
@@ -264,7 +264,7 @@ static INLINE void runner_dopair_grav_pp_full(
 #ifdef SWIFT_DEBUG_CHECKS
       /* Update the interaction counter if it's not a padded gpart */
       if (pjd < gcount_j && !gpart_is_inhibited(&gparts_j[pjd], e))
-        gparts_i[pid].num_interacted++;
+        atomic_inc(&gparts_i[pid].num_interacted);
 #endif
 
 #ifdef SWIFT_GRAVITY_FORCE_CHECKS
@@ -420,7 +420,7 @@ static INLINE void runner_dopair_grav_pp_truncated(
 #ifdef SWIFT_DEBUG_CHECKS
       /* Update the interaction counter if it's not a padded gpart */
       if (pjd < gcount_j && !gpart_is_inhibited(&gparts_j[pjd], e))
-        gparts_i[pid].num_interacted++;
+        atomic_inc(&gparts_i[pid].num_interacted);
 #endif
 
 #ifdef SWIFT_GRAVITY_FORCE_CHECKS
@@ -565,7 +565,8 @@ static INLINE void runner_dopair_grav_pm_full(
 #ifdef SWIFT_DEBUG_CHECKS
     /* Update the interaction counter */
     if (pid < gcount_i)
-      gparts_i[pid].num_interacted += cj->grav.multipole->m_pole.num_gpart;
+      atomic_add(&gparts_i[pid].num_interacted,
+                 cj->grav.multipole->m_pole.num_gpart);
 #endif
 
 #ifdef SWIFT_GRAVITY_FORCE_CHECKS
@@ -706,7 +707,8 @@ static INLINE void runner_dopair_grav_pm_truncated(
 #ifdef SWIFT_DEBUG_CHECKS
     /* Update the interaction counter */
     if (pid < gcount_i)
-      gparts_i[pid].num_interacted += cj->grav.multipole->m_pole.num_gpart;
+      atomic_add(&gparts_i[pid].num_interacted,
+                 cj->grav.multipole->m_pole.num_gpart);
 #endif
 
 #ifdef SWIFT_GRAVITY_FORCE_CHECKS
@@ -1041,7 +1043,7 @@ static INLINE void runner_doself_grav_pp_full(
 #ifdef SWIFT_DEBUG_CHECKS
       /* Update the interaction counter if it's not a padded gpart */
       if (pjd < gcount && !gpart_is_inhibited(&gparts[pjd], e))
-        gparts[pid].num_interacted++;
+        atomic_inc(&gparts[pid].num_interacted);
 #endif
 
 #ifdef SWIFT_GRAVITY_FORCE_CHECKS
@@ -1180,7 +1182,7 @@ static INLINE void runner_doself_grav_pp_truncated(
 #ifdef SWIFT_DEBUG_CHECKS
       /* Update the interaction counter if it's not a padded gpart */
       if (pjd < gcount && !gpart_is_inhibited(&gparts[pjd], e))
-        gparts[pid].num_interacted++;
+        atomic_inc(&gparts[pid].num_interacted);
 #endif
 
 #ifdef SWIFT_GRAVITY_FORCE_CHECKS
@@ -1638,9 +1640,9 @@ void runner_dopair_recursive_grav(struct runner *r, struct cell *ci,
 
 #ifdef SWIFT_DEBUG_CHECKS
     if (cell_is_active_gravity(ci, e))
-      multi_i->pot.num_interacted += multi_j->m_pole.num_gpart;
+      atomic_add(&multi_i->pot.num_interacted, multi_j->m_pole.num_gpart);
     if (cell_is_active_gravity(cj, e))
-      multi_j->pot.num_interacted += multi_i->m_pole.num_gpart;
+      atomic_add(&multi_j->pot.num_interacted, multi_i->m_pole.num_gpart);
 #endif
 
 #ifdef SWIFT_GRAVITY_FORCE_CHECKS
@@ -1854,7 +1856,7 @@ void runner_do_grav_long_range(struct runner *r, struct cell *ci, int timer) {
 
 #ifdef SWIFT_DEBUG_CHECKS
         /* Need to account for the interactions we missed */
-        multi_i->pot.num_interacted += multi_j->m_pole.num_gpart;
+        atomic_add(&multi_i->pot.num_interacted, multi_j->m_pole.num_gpart);
 #endif
 
 #ifdef SWIFT_GRAVITY_FORCE_CHECKS
-- 
GitLab