Merge pull request #728 from danicholson/cluster-fragment-aggregate-fixes

Cluster/fragment/aggregate bugfixes

src/compute_aggregate_atom.cpp

@@ -52,6 +52,7 @@ ComputeAggregateAtom::ComputeAggregateAtom(LAMMPS *lmp, int narg, char **arg) :
   peratom_flag = 1;
   size_peratom_cols = 0;
   comm_forward = 1;
+  comm_reverse = 1;
 
   nmax = 0;
 }
@@ -120,8 +121,10 @@ void ComputeAggregateAtom::compute_peratom()
   }
 
   // invoke full neighbor list (will copy or build if necessary)
+  // on the first step of a run, set preflag to one in neighbor->build_one(...)
 
-  neighbor->build_one(list);
+  if (update->firststep == update->ntimestep) neighbor->build_one(list,1);
+  else neighbor->build_one(list);
 
   // if group is dynamic, insure ghost atom masks are current
 
@@ -164,6 +167,11 @@ void ComputeAggregateAtom::compute_peratom()
   while (1) {
     comm->forward_comm_compute(this);
 
+    // reverse communication when bonds are not stored on every processor
+
+    if (force->newton_bond)
+      comm->reverse_comm_compute(this);
+
     change = 0;
     while (1) {
       done = 1;
@@ -253,13 +261,50 @@ void ComputeAggregateAtom::unpack_forward_comm(int n, int first, double *buf)
   m = 0;
   last = first + n;
   if (commflag)
-    for (i = first; i < last; i++) aggregateID[i] = buf[m++];
+    for (i = first; i < last; i++) {
+      double x = buf[m++];
+
+      // only overwrite ghost IDs with values lower than current ones
+
+      aggregateID[i] = MIN(x,aggregateID[i]);
+    }
   else {
     int *mask = atom->mask;
     for (i = first; i < last; i++) mask[i] = (int) ubuf(buf[m++]).i;
   }
 }
 
+/* ---------------------------------------------------------------------- */
+
+int ComputeAggregateAtom::pack_reverse_comm(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    buf[m++] = aggregateID[i];
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
+void ComputeAggregateAtom::unpack_reverse_comm(int n, int *list, double *buf)
+{
+  int i,j,m;
+
+  m = 0;
+  for (i = 0; i < n; i++) {
+    j = list[i];
+    double x = buf[m++];
+
+    // only overwrite local IDs with values lower than current ones
+
+    aggregateID[j] = MIN(x,aggregateID[j]);
+  }
+}
+
 /* ----------------------------------------------------------------------
    memory usage of local atom-based array
 ------------------------------------------------------------------------- */

src/compute_aggregate_atom.h

@@ -33,6 +33,8 @@ class ComputeAggregateAtom : public Compute {
   void compute_peratom();
   int pack_forward_comm(int, int *, double *, int, int *);
   void unpack_forward_comm(int, int, double *);
+  int pack_reverse_comm(int, int, double *);
+  void unpack_reverse_comm(int, int *, double *);
   double memory_usage();
 
  private:

src/compute_cluster_atom.cpp

@@ -112,8 +112,10 @@ void ComputeClusterAtom::compute_peratom()
   }
 
   // invoke full neighbor list (will copy or build if necessary)
+  // on the first step of a run, set preflag to one in neighbor->build_one(...)
 
-  neighbor->build_one(list);
+  if (update->firststep == update->ntimestep) neighbor->build_one(list,1);
+  else neighbor->build_one(list);
 
   inum = list->inum;
   ilist = list->ilist;

src/compute_fragment_atom.cpp

@@ -44,6 +44,7 @@ ComputeFragmentAtom::ComputeFragmentAtom(LAMMPS *lmp, int narg, char **arg) :
   peratom_flag = 1;
   size_peratom_cols = 0;
   comm_forward = 1;
+  comm_reverse = 1;
 
   nmax = 0;
 }
@@ -122,6 +123,11 @@ void ComputeFragmentAtom::compute_peratom()
   while (1) {
     comm->forward_comm_compute(this);
 
+    // reverse communication when bonds are not stored on every processor
+
+    if (force->newton_bond)
+      comm->reverse_comm_compute(this);
+
     change = 0;
     while (1) {
       done = 1;
@@ -183,13 +189,50 @@ void ComputeFragmentAtom::unpack_forward_comm(int n, int first, double *buf)
   m = 0;
   last = first + n;
   if (commflag)
-    for (i = first; i < last; i++) fragmentID[i] = buf[m++];
+    for (i = first; i < last; i++) {
+      double x = buf[m++];
+
+      // only overwrite ghost IDs with values lower than current ones
+
+      fragmentID[i] = MIN(x,fragmentID[i]);
+    }
   else {
     int *mask = atom->mask;
     for (i = first; i < last; i++) mask[i] = (int) ubuf(buf[m++]).i;
   }
 }
 
+/* ---------------------------------------------------------------------- */
+
+int ComputeFragmentAtom::pack_reverse_comm(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    buf[m++] = fragmentID[i];
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
+void ComputeFragmentAtom::unpack_reverse_comm(int n, int *list, double *buf)
+{
+  int i,j,m;
+
+  m = 0;
+  for (i = 0; i < n; i++) {
+    j = list[i];
+    double x = buf[m++];
+
+    // only overwrite local IDs with values lower than current ones
+
+    fragmentID[j] = MIN(x,fragmentID[j]);
+  }
+}
+
 /* ----------------------------------------------------------------------
    memory usage of local atom-based array
 ------------------------------------------------------------------------- */

src/compute_fragment_atom.h

@@ -32,6 +32,8 @@ class ComputeFragmentAtom : public Compute {
   void compute_peratom();
   int pack_forward_comm(int, int *, double *, int, int *);
   void unpack_forward_comm(int, int, double *);
+  int pack_reverse_comm(int, int, double *);
+  void unpack_reverse_comm(int, int *, double *);
   double memory_usage();
 
  private:

src/compute_hexorder_atom.cpp

@@ -154,8 +154,10 @@ void ComputeHexOrderAtom::compute_peratom()
   }
 
   // invoke full neighbor list (will copy or build if necessary)
+  // on the first step of a run, set preflag to one in neighbor->build_one(...)
 
-  neighbor->build_one(list);
+  if (update->firststep == update->ntimestep) neighbor->build_one(list,1);
+  else neighbor->build_one(list);
 
   inum = list->inum;
   ilist = list->ilist;