WIP: Remove ghost shifts

[fweik]

Removes the shifts on the ghosts. So far normal MD works, LB coupling still needs fixing, and the layered cell system does not work.

[fweik]

@RudolfWeeber do you have an idea why the virtual sites relative is broken?

[bindgens1]

Thanks a lot!

[espresso-ci]

Your pull request does not meet our code formatting rules. Specifically, I suggest you make the following changes:

diff --git a/src/core/cells.cpp b/src/core/cells.cpp
index e9050a5..f70dcce 100644
--- a/src/core/cells.cpp
+++ b/src/core/cells.cpp
@@ -167,8 +167,9 @@ static void topology_release(int cs) {
     layered_topology_release();
     break;
   default:
-    fprintf(stderr, "INTERNAL ERROR: attempting to sort the particles in an "
-                    "unknown way (%d)\n",
+    fprintf(stderr,
+            "INTERNAL ERROR: attempting to sort the particles in an "
+            "unknown way (%d)\n",
             cs);
     errexit();
   }
@@ -195,8 +196,9 @@ void topology_init(int cs, CellPList *local) {
     layered_topology_init(local, node_grid);
     break;
   default:
-    fprintf(stderr, "INTERNAL ERROR: attempting to sort the particles in an "
-                    "unknown way (%d)\n",
+    fprintf(stderr,
+            "INTERNAL ERROR: attempting to sort the particles in an "
+            "unknown way (%d)\n",
             cs);
     errexit();
   }
diff --git a/src/core/domain_decomposition.cpp b/src/core/domain_decomposition.cpp
index f45beec..2d8c743 100644
--- a/src/core/domain_decomposition.cpp
+++ b/src/core/domain_decomposition.cpp
@@ -460,57 +460,59 @@ void dd_init_cell_interactions() {
   }
 
   /* loop all local cells */
-    for (int o = 1; o < dd.cell_grid[2] + 1; o++)
-      for (int n = 1; n < dd.cell_grid[1] + 1; n++)
-        for (int m = 1; m < dd.cell_grid[0] + 1; m++) {
-          using boost::algorithm::clamp;
+  for (int o = 1; o < dd.cell_grid[2] + 1; o++)
+    for (int n = 1; n < dd.cell_grid[1] + 1; n++)
+      for (int m = 1; m < dd.cell_grid[0] + 1; m++) {
+        using boost::algorithm::clamp;
 
-    auto const gind1 = global_index({m, n, o});
-    auto const lind1 = get_linear_index(m, n, o, dd.ghost_cell_grid);
+        auto const gind1 = global_index({m, n, o});
+        auto const lind1 = get_linear_index(m, n, o, dd.ghost_cell_grid);
 
         std::vector<Cell *> red_neighbors;
         std::vector<Cell *> black_neighbors;
 
-    /* loop all neighbor cells */
-    int lower_index[3] = {m - 1, n - 1, o - 1};
-    int upper_index[3] = {m + 1, n + 1, o + 1};
+        /* loop all neighbor cells */
+        int lower_index[3] = {m - 1, n - 1, o - 1};
+        int upper_index[3] = {m + 1, n + 1, o + 1};
 
-    for (int i = 0; i < 3; i++) {
-      if (dd.fully_connected[i]) {
-        lower_index[i] = 0;
-        upper_index[i] = dd.ghost_cell_grid[i] - 1;
-      }
+        for (int i = 0; i < 3; i++) {
+          if (dd.fully_connected[i]) {
+            lower_index[i] = 0;
+            upper_index[i] = dd.ghost_cell_grid[i] - 1;
+          }
 
-      /* Cut of ghost layer at boundary if not periodic. */
-      if ((!PERIODIC(i)) && boundary[2 * i]) {
-        lower_index[i] = clamp(lower_index[i], 0, dd.ghost_cell_grid[i] - 1);
-      }
-      if ((!PERIODIC(i)) && boundary[2 * i + 1]) {
-        upper_index[i] = clamp(upper_index[i], 0, dd.ghost_cell_grid[i] - 1);
-      }
-    }
+          /* Cut of ghost layer at boundary if not periodic. */
+          if ((!PERIODIC(i)) && boundary[2 * i]) {
+            lower_index[i] =
+                clamp(lower_index[i], 0, dd.ghost_cell_grid[i] - 1);
+          }
+          if ((!PERIODIC(i)) && boundary[2 * i + 1]) {
+            upper_index[i] =
+                clamp(upper_index[i], 0, dd.ghost_cell_grid[i] - 1);
+          }
+        }
 
-    for (int p = lower_index[2]; p <= upper_index[2]; p++)
-      for (int q = lower_index[1]; q <= upper_index[1]; q++)
-        for (int r = lower_index[0]; r <= upper_index[0]; r++) {
-          auto const gind2 = global_index({r, q, p});
+        for (int p = lower_index[2]; p <= upper_index[2]; p++)
+          for (int q = lower_index[1]; q <= upper_index[1]; q++)
+            for (int r = lower_index[0]; r <= upper_index[0]; r++) {
+              auto const gind2 = global_index({r, q, p});
 
-          if (gind1 == gind2) {
-            continue;
-          }
+              if (gind1 == gind2) {
+                continue;
+              }
 
-          auto const lind2 = get_linear_index(r, q, p, dd.ghost_cell_grid);
+              auto const lind2 = get_linear_index(r, q, p, dd.ghost_cell_grid);
 
-          if (gind2 > gind1) {
-            red_neighbors.push_back(&cells[lind2]);
-          } else {
-            black_neighbors.push_back(&cells[lind2]);
-          }
-        }
+              if (gind2 > gind1) {
+                red_neighbors.push_back(&cells[lind2]);
+              } else {
+                black_neighbors.push_back(&cells[lind2]);
+              }
+            }
 
-    cells[lind1].m_neighbors =
-        Neighbors<Cell *>(red_neighbors, black_neighbors);
-  }
+        cells[lind1].m_neighbors =
+            Neighbors<Cell *>(red_neighbors, black_neighbors);
+      }
 }
 
 /*************************************************/
diff --git a/testsuite/python/virtual_sites_relative.py b/testsuite/python/virtual_sites_relative.py
index bcd7b9e..7fcd696 100644
--- a/testsuite/python/virtual_sites_relative.py
+++ b/testsuite/python/virtual_sites_relative.py
@@ -302,13 +302,13 @@ class VirtualSites(ut.TestCase):
         system.cell_system.skin = 0.4
         system.cell_system.set_n_square(use_verlet_lists=True)
         self.run_test_lj()
-	print("n2")
+        print("n2")
         system.cell_system.set_domain_decomposition(use_verlet_lists=True)
         self.run_test_lj()
-	print("dd-vl")
+        print("dd-vl")
         system.cell_system.set_domain_decomposition(use_verlet_lists=False)
         self.run_test_lj()
-	print("dd")
+        print("dd")
 
     @utx.skipIfMissingFeatures("EXTERNAL_FORCES")
     def test_zz_stress_tensor(self):

To apply these changes, please do one of the following:

• You can download a patch with my suggested changes here, inspect it and make changes manually.
• You can directly apply it to your repository by running curl https://gitlab.icp.uni-stuttgart.de/espressomd/espresso/-/jobs/123987/artifacts/raw/style.patch | git apply -.
• You can run maintainer/CI/fix_style.sh to automatically fix your coding style. This is the same command that I have executed to generate the patch above, but it requires certain tools to be installed on your computer.

You can run gitlab-runner exec docker style afterwards to check if your changes worked out properly.

Please note that there are often multiple ways to correctly format code. As I am just a robot, I sometimes fail to identify the most aesthetically pleasing way. So please look over my suggested changes and adapt them where the style does not make sense.

[pdebuyl]

Hi @fweik

I have trouble making sense of this full PR. Would it be possible at all to have a PR against the current python branch for which only the ghost shift changes are visible? Or is there a subset of the commits here that I could read? Thanks!

[fweik]

I'm not sure that I follow, everything in this PR is related to the ghost shifts (arguably except for the removal of the centor of mass virtual sites, but you can just ignore those two filed).

[pdebuyl]

Hi @fweik I did look at the PR at some point in time where it was not in sync so the changes were huge.

The current version is ok, thanks for the reply :-)

[pdebuyl]

One of the failures in the tests is testsuite/python/cellsystem.py. With two MPI processes, I get

x: array([[0, 1],
[0, 1],
[0, 1],
[0, 1]])

for get_pairs_. From my reading of the code, this means that link_cell goes over all cell pairs (hence one finds i,j and j,i) and the results are collected across nodes (so that any node having both i and j (one as a normal particle the other as a ghost) will return the pairs).

I don't know enough about the cell iterators internals to go further unfortunately.