diff --git a/src/change_box.cpp b/src/change_box.cpp
index add1fe1d0a49b57c128b945729c25cd8e03d8272..c6ec078a01f647d4059b692b03250bc5bd52012d 100644
--- a/src/change_box.cpp
+++ b/src/change_box.cpp
@@ -316,6 +316,9 @@ void ChangeBox::command(int narg, char **arg)
} else if (ops[m].style == REMAP) {
+ if (modify->check_rigid_group_overlap(groupbit))
+ error->warning(FLERR,"Attempting to remap atoms in rigid bodies");
+
// convert atoms to lamda coords, using last box state
// convert atoms back to box coords, using current box state
// save current box state
diff --git a/src/delete_atoms.cpp b/src/delete_atoms.cpp
index 116c2a2f1e057431644c82b2747d7d95eac6095c..85bd6f88ddf43b0d5a5a96612b726fb8fc2664ea 100644
--- a/src/delete_atoms.cpp
+++ b/src/delete_atoms.cpp
@@ -69,6 +69,15 @@ void DeleteAtoms::command(int narg, char **arg)
else if (strcmp(arg[0],"porosity") == 0) delete_porosity(narg,arg);
else error->all(FLERR,"Illegal delete_atoms command");
+ if (allflag) {
+ int igroup = group->find("all");
+ if ((igroup >= 0) && modify->check_rigid_group_overlap(group->bitmask[igroup]))
+ error->warning(FLERR,"Attempting to delete atoms in rigid bodies");
+ } else {
+ if (modify->check_rigid_list_overlap(dlist))
+ error->warning(FLERR,"Attempting to delete atoms in rigid bodies");
+ }
+
// if allflag = 1, just reset atom->nlocal
// else delete atoms one by one
@@ -89,16 +98,16 @@ void DeleteAtoms::command(int narg, char **arg)
int i = 0;
while (i < nlocal) {
if (dlist[i]) {
- avec->copy(nlocal-1,i,1);
- dlist[i] = dlist[nlocal-1];
- nlocal--;
+ avec->copy(nlocal-1,i,1);
+ dlist[i] = dlist[nlocal-1];
+ nlocal--;
} else i++;
}
-
+
atom->nlocal = nlocal;
memory->destroy(dlist);
}
-
+
// if non-molecular system and compress flag set,
// reset atom tags to be contiguous
// set all atom IDs to 0, call tag_extend()
@@ -201,7 +210,7 @@ void DeleteAtoms::delete_group(int narg, char **arg)
allflag = 1;
return;
}
-
+
// allocate and initialize deletion list
int nlocal = atom->nlocal;
diff --git a/src/displace_atoms.cpp b/src/displace_atoms.cpp
index 7db7e218391ff2c28f0a902b6a59e595c9877da8..a9aa5cf8558aa2c620d9e5ce99ab882d7eb7af4b 100644
--- a/src/displace_atoms.cpp
+++ b/src/displace_atoms.cpp
@@ -75,6 +75,9 @@ void DisplaceAtoms::command(int narg, char **arg)
if (igroup == -1) error->all(FLERR,"Could not find displace_atoms group ID");
groupbit = group->bitmask[igroup];
+ if (modify->check_rigid_group_overlap(groupbit))
+ error->warning(FLERR,"Attempting to displace atoms in rigid bodies");
+
int style = -1;
if (strcmp(arg[1],"move") == 0) style = MOVE;
else if (strcmp(arg[1],"ramp") == 0) style = RAMP;
diff --git a/src/fix_heat.cpp b/src/fix_heat.cpp
index d41aa4abea90d36adf7ee59f00441fa6d4ad47fd..97e0ed6a7fd76662087d4db5415d9808ec535c39 100644
--- a/src/fix_heat.cpp
+++ b/src/fix_heat.cpp
@@ -64,7 +64,7 @@ idregion(NULL), hstr(NULL), vheat(NULL), vscale(NULL)
// optional args
iregion = -1;
-
+
int iarg = 5;
while (iarg < narg) {
if (strcmp(arg[iarg],"region") == 0) {
@@ -126,6 +126,10 @@ void FixHeat::init()
else error->all(FLERR,"Variable for fix heat is invalid style");
}
+ // check for rigid bodies in region (done here for performance reasons)
+ if (modify->check_rigid_region_overlap(groupbit,domain->regions[iregion]))
+ error->warning(FLERR,"Cannot apply fix heat to atoms in rigid bodies");
+
// cannot have 0 atoms in group
if (group->count(igroup) == 0)
diff --git a/src/fix_temp_berendsen.cpp b/src/fix_temp_berendsen.cpp
index aff9a44977606a1b46716bed2f3e1b47f3b729ff..7b312cfb5f35dbc1777aeb9e63c3a530ed1e9d9c 100644
--- a/src/fix_temp_berendsen.cpp
+++ b/src/fix_temp_berendsen.cpp
@@ -128,6 +128,9 @@ void FixTempBerendsen::init()
error->all(FLERR,"Temperature ID for fix temp/berendsen does not exist");
temperature = modify->compute[icompute];
+ if (modify->check_rigid_group_overlap(groupbit))
+ error->warning(FLERR,"Cannot thermostat atoms in rigid bodies");
+
if (temperature->tempbias) which = BIAS;
else which = NOBIAS;
}
diff --git a/src/fix_temp_csld.cpp b/src/fix_temp_csld.cpp
index f24314ac80ef07d8a520816530deebb28d661d0d..63f27cdecb17968ecbc6f265ef78a45018a70535 100644
--- a/src/fix_temp_csld.cpp
+++ b/src/fix_temp_csld.cpp
@@ -155,6 +155,9 @@ void FixTempCSLD::init()
error->all(FLERR,"Temperature ID for fix temp/csld does not exist");
temperature = modify->compute[icompute];
+ if (modify->check_rigid_group_overlap(groupbit))
+ error->warning(FLERR,"Cannot thermostat atoms in rigid bodies");
+
if (temperature->tempbias) which = BIAS;
else which = NOBIAS;
}
diff --git a/src/modify.cpp b/src/modify.cpp
index 7af4576038c33d855325fb224630658fa0e7764c..01de6b59284676069de3980362accd8aa35badd6 100644
--- a/src/modify.cpp
+++ b/src/modify.cpp
@@ -23,6 +23,7 @@
#include "group.h"
#include "update.h"
#include "domain.h"
+#include "region.h"
#include "input.h"
#include "variable.h"
#include "memory.h"
@@ -995,6 +996,99 @@ int Modify::check_package(const char *package_fix_name)
return 1;
}
+
+/* ----------------------------------------------------------------------
+ check if the group indicated by groupbit overlaps with any
+ currently existing rigid fixes. return 1 in this case otherwise 0
+------------------------------------------------------------------------- */
+
+int Modify::check_rigid_group_overlap(int groupbit)
+{
+ const int * const mask = atom->mask;
+ const int nlocal = atom->nlocal;
+ int dim;
+
+ int n = 0;
+ for (int ifix = 0; ifix < nfix; ifix++) {
+ if (strncmp("rigid",fix[ifix]->style,5) == 0) {
+ const int * const body = (const int *)fix[ifix]->extract("body",dim);
+ if ((body == NULL) || (dim != 1)) break;
+
+ for (int i=0; (i < nlocal) && (n == 0); ++i)
+ if ((mask[i] & groupbit) && (body[i] >= 0)) ++n;
+ }
+ }
+
+ int n_all = 0;
+ MPI_Allreduce(&n,&n_all,1,MPI_INT,MPI_SUM,world);
+
+ if (n_all > 0) return 1;
+ return 0;
+}
+
+/* ----------------------------------------------------------------------
+ check if the atoms in the group indicated by groupbit _and_ region
+ indicated by regionid overlap with any currently existing rigid fixes.
+ return 1 in this case, otherwise 0
+------------------------------------------------------------------------- */
+
+int Modify::check_rigid_region_overlap(int groupbit, Region *reg)
+{
+ const int * const mask = atom->mask;
+ const double * const * const x = atom->x;
+ const int nlocal = atom->nlocal;
+ int dim;
+
+ int n = 0;
+ reg->prematch();
+ for (int ifix = 0; ifix < nfix; ifix++) {
+ if (strncmp("rigid",fix[ifix]->style,5) == 0) {
+ const int * const body = (const int *)fix[ifix]->extract("body",dim);
+ if ((body == NULL) || (dim != 1)) break;
+
+ for (int i=0; (i < nlocal) && (n == 0); ++i)
+ if ((mask[i] & groupbit) && (body[i] >= 0)
+ && reg->match(x[i][0],x[i][1],x[i][2])) ++n;
+ }
+ }
+
+ int n_all = 0;
+ MPI_Allreduce(&n,&n_all,1,MPI_INT,MPI_SUM,world);
+
+ if (n_all > 0) return 1;
+ return 0;
+}
+
+/* ----------------------------------------------------------------------
+ check if the atoms in the selection list (length atom->nlocal,
+ content: 1 if atom is contained, 0 if not) overlap with currently
+ existing rigid fixes. return 1 in this case otherwise 0
+------------------------------------------------------------------------- */
+
+int Modify::check_rigid_list_overlap(int *select)
+{
+ const int * const mask = atom->mask;
+ const int nlocal = atom->nlocal;
+ int dim;
+
+ int n = 0;
+ for (int ifix = 0; ifix < nfix; ifix++) {
+ if (strncmp("rigid",fix[ifix]->style,5) == 0) {
+ const int * const body = (const int *)fix[ifix]->extract("body",dim);
+ if ((body == NULL) || (dim != 1)) break;
+
+ for (int i=0; (i < nlocal) && (n == 0); ++i)
+ if ((body[i] >= 0) && select[i]) ++n;
+ }
+ }
+
+ int n_all = 0;
+ MPI_Allreduce(&n,&n_all,1,MPI_INT,MPI_SUM,world);
+
+ if (n_all > 0) return 1;
+ return 0;
+}
+
/* ----------------------------------------------------------------------
add a new compute
------------------------------------------------------------------------- */
diff --git a/src/modify.h b/src/modify.h
index 3ded3cbab6ea0d4a586358e6c850a03b8d9b4055..d825d5c4efbf543a709e766253a0a7e600cf8dcf 100644
--- a/src/modify.h
+++ b/src/modify.h
@@ -98,6 +98,9 @@ class Modify : protected Pointers {
int find_fix(const char *);
int find_fix_by_style(const char *);
int check_package(const char *);
+ int check_rigid_group_overlap(int);
+ int check_rigid_region_overlap(int, class Region *);
+ int check_rigid_list_overlap(int *);
void add_compute(int, char **, int trysuffix=1);
void modify_compute(int, char **);
diff --git a/src/set.cpp b/src/set.cpp
index 4ed07d423b112fd59c51ee8e3f24c73dbe8ab1b3..629f325d3539394db4085b146322e39a6a7acacb 100644
--- a/src/set.cpp
+++ b/src/set.cpp
@@ -585,6 +585,28 @@ void Set::set(int keyword)
}
}
+ // check if properties of atoms in rigid bodies are updated
+ // that are cached as per-body data.
+ switch (keyword) {
+ case X:
+ case Y:
+ case Z:
+ case MOLECULE:
+ case MASS:
+ case ANGMOM:
+ case SHAPE:
+ case DIAMETER:
+ case DENSITY:
+ case QUAT:
+ case IMAGE:
+ if (modify->check_rigid_list_overlap(select))
+ error->warning(FLERR,"Changing a property of atoms in rigid bodies "
+ "that has no effect unless rigid bodies are rebuild");
+ break;
+ default: // assume no conflict for all other properties
+ break;
+ }
+
// loop over selected atoms
AtomVecEllipsoid *avec_ellipsoid =
diff --git a/src/velocity.cpp b/src/velocity.cpp
index 82b6efbe1b8cadc2cb798feb301d7b2a3d9a146d..260a11bb4ee90dc2116907688e9b5bda4503e4cb 100644
--- a/src/velocity.cpp
+++ b/src/velocity.cpp
@@ -68,6 +68,12 @@ void Velocity::command(int narg, char **arg)
if (igroup == -1) error->all(FLERR,"Could not find velocity group ID");
groupbit = group->bitmask[igroup];
+ // check if velocities of atoms in rigid bodies are updated
+
+ if (modify->check_rigid_group_overlap(groupbit))
+ error->warning(FLERR,"Changing velocities of atoms in rigid bodies. "
+ "This has no effect unless rigid bodies are rebuild");
+
// identify style
if (strcmp(arg[1],"create") == 0) style = CREATE;