Merge #3288 #3290

3288: core: Simplify particle removal r=RudolfWeeber a=fweik

Fixes #3286.

Description of changes:
 - Simplify particle removal by doing one linear pass over all particles
   instead using index


3290: fix writevtk() bug where types=all was ignored r=jngrad a=christophlohrmann

Description of changes:
 - fixed bug where the argument ``types='all'`` was ignored by ``system.part.writevtk()`` 


PR Checklist
------------
 - [ ] Tests?
   - [ ] Interface
   - [ ] Core 
 - [ ] Docs?


Co-authored-by: Florian Weik <fweik@icp.uni-stuttgart.de>
Co-authored-by: Christoph Lohrmann <clohrmann@icp.uni-stuttgart.de>
Co-authored-by: christophlohrmann <47663906+christophlohrmann@users.noreply.github.com>

src/core/communication.cpp

@@ -280,10 +280,11 @@ void mpi_remove_particle_slave(int pnode, int part) {
   if (part != -1) {
     n_part--;
 
-    if (pnode == this_node)
+    if (pnode == this_node) {
       local_remove_particle(part);
-
-    remove_all_bonds_to(part);
+    } else {
+      remove_all_bonds_to(part);
+    }
   } else
     local_remove_all_particles();
 

src/core/particle_data.cpp

@@ -1036,50 +1036,54 @@ int remove_particle(int p_id) {
   return ES_OK;
 }
 
-namespace {
-std::pair<Cell *, size_t> find_particle(Particle *p, Cell *c) {
-  for (int i = 0; i < c->n; ++i) {
-    if ((c->part + i) == p) {
-      return {c, i};
-    }
+/**
+ * @brief Remove all bonds on particle involing other particle.
+ *
+ * @param p Particle whose bond list is modified.
+ * @param id Bonds involving this id are removed.
+ */
+static void remove_all_bonds_to(Particle &p, int id) {
+  IntList *bl = &p.bl;
+  int i, j, partners;
+
+  for (i = 0; i < bl->n;) {
+    partners = bonded_ia_params[bl->e[i]].num;
+    for (j = 1; j <= partners; j++)
+      if (bl->e[i + j] == id)
+        break;
+    if (j <= partners) {
+      bl->erase(bl->begin() + i, bl->begin() + i + 1 + partners);
+    } else
+      i += 1 + partners;
   }
-  return {nullptr, 0};
+  assert(i == bl->n);
 }
 
-std::pair<Cell *, size_t> find_particle(Particle *p, CellPList cells) {
-  for (auto &c : cells) {
-    auto res = find_particle(p, c);
-    if (res.first) {
-      return res;
-    }
+void remove_all_bonds_to(int identity) {
+  for (auto &p : local_cells.particles()) {
+    remove_all_bonds_to(p, identity);
   }
-
-  return {nullptr, 0};
 }
-} // namespace
 
 void local_remove_particle(int part) {
-  Particle *p = local_particles[part];
-  assert(p);
-  assert(not p->l.ghost);
-
-  /* If the particles are sorted we can use the
-   * cell system to find the cell containing the
-   * particle. Otherwise we do a brute force search
-   * of the cells. */
   Cell *cell = nullptr;
-  size_t n = 0;
-  if (Cells::RESORT_NONE == get_resort_particles()) {
-    std::tie(cell, n) = find_particle(p, find_current_cell(*p));
-  }
-
-  if (not cell) {
-    std::tie(cell, n) = find_particle(p, local_cells);
+  int position = -1;
+  for (auto c : local_cells) {
+    for (int i = 0; i < c->n; i++) {
+      auto &p = c->part[i];
+
+      if (p.identity() == part) {
+        cell = c;
+        position = i;
+      } else {
+        remove_all_bonds_to(p, i);
+      }
+    }
   }
 
-  assert(cell && cell->part && (n < cell->n) && ((cell->part + n) == p));
+  assert(cell && (position >= 0));
 
-  Particle p_destroy = extract_indexed_particle(cell, n);
+  extract_indexed_particle(cell, position);
 }
 
 Particle *local_place_particle(int id, const Utils::Vector3d &pos, int _new) {
@@ -1193,31 +1197,6 @@ int try_delete_bond(Particle *part, const int *bond) {
   return ES_ERROR;
 }
 
-void remove_all_bonds_to(int identity) {
-  for (auto &p : local_cells.particles()) {
-    IntList *bl = &p.bl;
-    int i, j, partners;
-
-    for (i = 0; i < bl->n;) {
-      partners = bonded_ia_params[bl->e[i]].num;
-      for (j = 1; j <= partners; j++)
-        if (bl->e[i + j] == identity)
-          break;
-      if (j <= partners) {
-        bl->erase(bl->begin() + i, bl->begin() + i + 1 + partners);
-      } else
-        i += 1 + partners;
-    }
-    if (i != bl->n) {
-      fprintf(stderr,
-              "%d: INTERNAL ERROR: bond information corrupt for "
-              "particle %d, exiting...\n",
-              this_node, p.p.identity);
-      errexit();
-    }
-  }
-}
-
 #ifdef EXCLUSIONS
 void local_change_exclusion(int part1, int part2, int _delete) {
   if (part1 == -1 && part2 == -1) {

src/python/espressomd/particle_data.pyx

@@ -1931,9 +1931,8 @@ Set quat and scalar dipole moment (dipm) instead.")
 
         n = 0
         for p in self:
-            for t in types:
-                if p.type == t or t == "all":
-                    n += 1
+            if types == 'all' or p.type in types:
+                n += 1
 
         with open(fname, "w") as vtk:
             vtk.write("# vtk DataFile Version 2.0\n")
@@ -1942,17 +1941,15 @@ Set quat and scalar dipole moment (dipm) instead.")
             vtk.write("DATASET UNSTRUCTURED_GRID\n")
             vtk.write("POINTS {} floats\n".format(n))
             for p in self:
-                for t in types:
-                    if p.type == t or t == "all":
-                        vtk.write("{} {} {}\n".format(*(p.pos_folded)))
+                if types == 'all' or p.type in types:
+                    vtk.write("{} {} {}\n".format(*(p.pos_folded)))
 
             vtk.write("POINT_DATA {}\n".format(n))
             vtk.write("SCALARS velocity float 3\n")
             vtk.write("LOOKUP_TABLE default\n")
             for p in self:
-                for t in types:
-                    if p.type == t or t == "all":
-                        vtk.write("{} {} {}\n".format(*p.v))
+                if types == 'all' or p.type in types:
+                    vtk.write("{} {} {}\n".format(*p.v))
 
     property highest_particle_id:
         """