diff --git a/doc/src/group.rst b/doc/src/group.rst index 1720ecfe1a..9d93e949f9 100644 --- a/doc/src/group.rst +++ b/doc/src/group.rst @@ -258,11 +258,17 @@ assignment is made at the beginning of the minimization, but not during the iterations of the minimizer. The point in the timestep at which atoms are assigned to a dynamic -group is after the initial stage of velocity Verlet time integration -has been performed, and before neighbor lists or forces are computed. -This is the point in the timestep where atom positions have just -changed due to the time integration, so the region criterion should be -accurate, if applied. +group is after interatomic forces have been computed, but before any +fixes which alter forces or otherwise update the system have been +invoked. This means that atom positions have been updated, neighbor +lists and ghost atoms are current, and both intermolecular and +intramolecular forces have been calculated based on the new +coordinates. Thus the region criterion, if applied, should be +accurate. Also, any computes invoked by an atom-style variable should +use updated information for that timestep, e.g. potential energy/atom +or coordination number/atom. Similarly, fixes or computes which are +invoked after that point in the timestep, should operate on the new +group of atoms. .. note:: diff --git a/src/compute_cluster_atom.cpp b/src/compute_cluster_atom.cpp index 07a0c700ee..bb4a7122e7 100644 --- a/src/compute_cluster_atom.cpp +++ b/src/compute_cluster_atom.cpp @@ -30,8 +30,6 @@ using namespace LAMMPS_NS; -enum { CLUSTER, MASK, COORDS }; - /* ---------------------------------------------------------------------- */ ComputeClusterAtom::ComputeClusterAtom(LAMMPS *lmp, int narg, char **arg) : @@ -44,7 +42,7 @@ ComputeClusterAtom::ComputeClusterAtom(LAMMPS *lmp, int narg, char **arg) : peratom_flag = 1; size_peratom_cols = 0; - comm_forward = 3; + comm_forward = 1; nmax = 0; } @@ -117,22 +115,6 @@ void ComputeClusterAtom::compute_peratom() numneigh = list->numneigh; firstneigh = list->firstneigh; - // if update->post_integrate set: - // a dynamic group in FixGroup is invoking a variable with this compute - // thus ghost atom coords need to be up-to-date after initial_integrate() - - if (update->post_integrate) { - commflag = COORDS; - comm->forward_comm(this); - } - - // if group is dynamic, insure ghost atom masks are current - - if (group->dynamic[igroup]) { - commflag = MASK; - comm->forward_comm(this); - } - // every atom starts in its own cluster, with clusterID = atomID tagint *tag = atom->tag; @@ -153,7 +135,6 @@ void ComputeClusterAtom::compute_peratom() // iterate until no changes in my atoms // then check if any proc made changes - commflag = CLUSTER; double **x = atom->x; int change, done, anychange; @@ -203,31 +184,15 @@ void ComputeClusterAtom::compute_peratom() /* ---------------------------------------------------------------------- */ -int ComputeClusterAtom::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, - int * /*pbc*/) +int ComputeClusterAtom::pack_forward_comm(int n, int *list, double *buf, + int /*pbc_flag*/, int * /*pbc*/) { int i, j, m; m = 0; - if (commflag == CLUSTER) { - for (i = 0; i < n; i++) { - j = list[i]; - buf[m++] = clusterID[j]; - } - } else if (commflag == MASK) { - int *mask = atom->mask; - for (i = 0; i < n; i++) { - j = list[i]; - buf[m++] = ubuf(mask[j]).d; - } - } else if (commflag == COORDS) { - double **x = atom->x; - for (i = 0; i < n; i++) { - j = list[i]; - buf[m++] = x[j][0]; - buf[m++] = x[j][1]; - buf[m++] = x[j][2]; - } + for (i = 0; i < n; i++) { + j = list[i]; + buf[m++] = clusterID[j]; } return m; @@ -241,19 +206,7 @@ void ComputeClusterAtom::unpack_forward_comm(int n, int first, double *buf) m = 0; last = first + n; - if (commflag == CLUSTER) { - for (i = first; i < last; i++) clusterID[i] = buf[m++]; - } else if (commflag == MASK) { - int *mask = atom->mask; - for (i = first; i < last; i++) mask[i] = (int) ubuf(buf[m++]).i; - } else if (commflag == COORDS) { - double **x = atom->x; - for (i = first; i < last; i++) { - x[i][0] = buf[m++]; - x[i][1] = buf[m++]; - x[i][2] = buf[m++]; - } - } + for (i = first; i < last; i++) clusterID[i] = buf[m++]; } /* ---------------------------------------------------------------------- diff --git a/src/compute_coord_atom.cpp b/src/compute_coord_atom.cpp index e7cd73df19..48646b3f38 100644 --- a/src/compute_coord_atom.cpp +++ b/src/compute_coord_atom.cpp @@ -259,14 +259,16 @@ void ComputeCoordAtom::compute_peratom() j = jlist[jj]; j &= NEIGHMASK; - jtype = type[j]; - delx = xtmp - x[j][0]; - dely = ytmp - x[j][1]; - delz = ztmp - x[j][2]; - rsq = delx * delx + dely * dely + delz * delz; - if (rsq < cutsq) { - for (m = 0; m < ncol; m++) - if (jtype >= typelo[m] && jtype <= typehi[m]) count[m] += 1.0; + if (mask[j] & jgroupbit) { + jtype = type[j]; + delx = xtmp - x[j][0]; + dely = ytmp - x[j][1]; + delz = ztmp - x[j][2]; + rsq = delx * delx + dely * dely + delz * delz; + if (rsq < cutsq) { + for (m = 0; m < ncol; m++) + if (jtype >= typelo[m] && jtype <= typehi[m]) count[m] += 1.0; + } } } } @@ -309,8 +311,8 @@ void ComputeCoordAtom::compute_peratom() /* ---------------------------------------------------------------------- */ -int ComputeCoordAtom::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, - int * /*pbc*/) +int ComputeCoordAtom::pack_forward_comm(int n, int *list, double *buf, + int /*pbc_flag*/, int * /*pbc*/) { int i, m = 0, j; for (i = 0; i < n; ++i) { diff --git a/src/fix_group.cpp b/src/fix_group.cpp index 1291721f0a..73143b2f4f 100644 --- a/src/fix_group.cpp +++ b/src/fix_group.cpp @@ -44,6 +44,8 @@ idregion(nullptr), idvar(nullptr), idprop(nullptr) gbit = group->bitmask[group->find(dgroupid)]; gbitinverse = group->inversemask[group->find(dgroupid)]; + comm_forward = 1; + // process optional args regionflag = 0; @@ -106,8 +108,6 @@ FixGroup::~FixGroup() int FixGroup::setmask() { int mask = 0; - mask |= POST_INTEGRATE; - mask |= POST_INTEGRATE_RESPA; return mask; } @@ -147,29 +147,6 @@ void FixGroup::init() if (iprop < 0 || cols) error->all(FLERR,"Group dynamic command custom property vector does not exist"); } - - // warn if any FixGroup is not at tail end of all post_integrate fixes - - Fix **fix = modify->fix; - int *fmask = modify->fmask; - int nfix = modify->nfix; - - int n = 0; - for (int i = 0; i < nfix; i++) if (POST_INTEGRATE & fmask[i]) n++; - int warn = 0; - for (int i = 0; i < nfix; i++) { - if (POST_INTEGRATE & fmask[i]) { - for (int j = i+1; j < nfix; j++) { - if (POST_INTEGRATE & fmask[j]) { - if (strstr(fix[j]->id,"GROUP_") != fix[j]->id) warn = 1; - } - } - } - } - - if (warn && comm->me == 0) - error->warning(FLERR,"One or more dynamic groups may not be " - "updated at correct point in timestep"); } /* ---------------------------------------------------------------------- @@ -183,7 +160,7 @@ void FixGroup::setup(int /*vflag*/) /* ---------------------------------------------------------------------- */ -void FixGroup::post_integrate() +void FixGroup::post_force(int /*vflag*/) { // only assign atoms to group on steps that are multiples of nevery @@ -192,9 +169,9 @@ void FixGroup::post_integrate() /* ---------------------------------------------------------------------- */ -void FixGroup::post_integrate_respa(int ilevel, int /*iloop*/) +void FixGroup::post_force_respa(int vflag, int ilevel, int /*iloop*/) { - if (ilevel == nlevels_respa-1) post_integrate(); + if (ilevel == nlevels_respa-1) post_force(vflag); } /* ---------------------------------------------------------------------- */ @@ -204,7 +181,6 @@ void FixGroup::set_group() int nlocal = atom->nlocal; // invoke atom-style variable if defined - // set post_integrate flag to 1, then unset after // this is for any compute to check if it needs to // operate differently due to invocation this early in timestep // e.g. perform ghost comm update due to atoms having just moved @@ -214,12 +190,10 @@ void FixGroup::set_group() double *dvector = nullptr; if (varflag) { - update->post_integrate = 1; modify->clearstep_compute(); - memory->create(var,nlocal,"fix/group:varvalue"); + memory->create(var,nlocal,"fix/group:var"); input->variable->compute_atom(ivar,igroup,var,1,0); modify->addstep_compute(update->ntimestep + nevery); - update->post_integrate = 0; } // set ptr to custom atom vector @@ -233,8 +207,6 @@ void FixGroup::set_group() // set mask for each atom // only in group if in parent group, in region, variable is non-zero - // if compute, fix, etc needs updated masks of ghost atoms, - // it must do forward_comm() to update them double **x = atom->x; int *mask = atom->mask; @@ -256,6 +228,42 @@ void FixGroup::set_group() } if (varflag) memory->destroy(var); + + // insure ghost atom masks are also updated + + comm->forward_comm(this); +} + +/* ---------------------------------------------------------------------- */ + +int FixGroup::pack_forward_comm(int n, int *list, double *buf, + int /*pbc_flag*/, int * /*pbc*/) +{ + int i, j, m; + + int *mask = atom->mask; + + m = 0; + for (i = 0; i < n; i++) { + j = list[i]; + buf[m++] = ubuf(mask[j]).d; + } + + return m; +} + +/* ---------------------------------------------------------------------- */ + +void FixGroup::unpack_forward_comm(int n, int first, double *buf) +{ + int i, m, last; + + m = 0; + last = first + n; + + int *mask = atom->mask; + + for (i = first; i < last; i++) mask[i] = (int) ubuf(buf[m++]).i; } /* ---------------------------------------------------------------------- */ diff --git a/src/fix_group.h b/src/fix_group.h index ffa2afcfc3..e4c4c9bd32 100644 --- a/src/fix_group.h +++ b/src/fix_group.h @@ -31,8 +31,10 @@ class FixGroup : public Fix { int setmask() override; void init() override; void setup(int) override; - void post_integrate() override; - void post_integrate_respa(int, int) override; + void post_force(int) override; + void post_force_respa(int, int, int) override; + int pack_forward_comm(int, int *, double *, int, int *) override; + void unpack_forward_comm(int, int, double *) override; void *extract(const char *, int &) override; private: diff --git a/src/modify.cpp b/src/modify.cpp index 7554079e2a..ec74b03554 100644 --- a/src/modify.cpp +++ b/src/modify.cpp @@ -58,6 +58,8 @@ Modify::Modify(LAMMPS *lmp) : Pointers(lmp) n_pre_force_respa = n_post_force_respa = n_final_integrate_respa = 0; n_min_pre_exchange = n_min_pre_force = n_min_pre_reverse = 0; n_min_post_force = n_min_energy = 0; + + n_post_force_group = 0; n_timeflag = -1; fix = nullptr; @@ -76,6 +78,8 @@ Modify::Modify(LAMMPS *lmp) : Pointers(lmp) end_of_step_every = nullptr; + list_post_force_group = nullptr; + list_timeflag = nullptr; nfix_restart_global = 0; @@ -156,6 +160,7 @@ Modify::~Modify() delete[] list_min_pre_reverse; delete[] list_min_post_force; delete[] list_min_energy; + delete[] list_post_force_group; delete[] end_of_step_every; delete[] list_timeflag; @@ -226,6 +231,7 @@ void Modify::init() list_init_energy_couple(n_energy_couple, list_energy_couple); list_init_energy_global(n_energy_global, list_energy_global); list_init_energy_atom(n_energy_atom, list_energy_atom); + list_init_post_force_group(n_post_force_group, list_post_force_group); list_init(INITIAL_INTEGRATE_RESPA, n_initial_integrate_respa, list_initial_integrate_respa); list_init(POST_INTEGRATE_RESPA, n_post_integrate_respa, list_post_integrate_respa); @@ -441,10 +447,16 @@ void Modify::pre_reverse(int eflag, int vflag) /* ---------------------------------------------------------------------- post_force call, only for relevant fixes + first call any instances of fix GROUP if they exist ------------------------------------------------------------------------- */ void Modify::post_force(int vflag) { + if (n_post_force_group) { + for (int i = 0; i < n_post_force_group; i++) + fix[list_post_force_group[i]]->post_force(vflag); + } + for (int i = 0; i < n_post_force; i++) fix[list_post_force[i]]->post_force(vflag); } @@ -585,10 +597,16 @@ void Modify::pre_force_respa(int vflag, int ilevel, int iloop) /* ---------------------------------------------------------------------- rRESPA post_force call, only for relevant fixes + first call any instances of fix GROUP if they exist ------------------------------------------------------------------------- */ void Modify::post_force_respa(int vflag, int ilevel, int iloop) { + if (n_post_force_group) { + for (int i = 0; i < n_post_force_group; i++) + fix[list_post_force_group[i]]->post_force_respa(vflag, ilevel, iloop); + } + for (int i = 0; i < n_post_force_respa; i++) fix[list_post_force_respa[i]]->post_force_respa(vflag, ilevel, iloop); } @@ -1716,6 +1734,26 @@ void Modify::list_init_energy_atom(int &n, int *&list) if (fix[i]->energy_peratom_flag && fix[i]->thermo_energy) list[n++] = i; } +/* ---------------------------------------------------------------------- + create list of fix indices for fix GROUP + are invoked first in post_force() or post_force_respa() +------------------------------------------------------------------------- */ + +void Modify::list_init_post_force_group(int &n, int *&list) +{ + delete[] list; + + n = 0; + for (int i = 0; i < nfix; i++) + if (strcmp(fix[i]->style,"GROUP") == 0) n++; + list = new int[n]; + + n = 0; + for (int i = 0; i < nfix; i++) + if (strcmp(fix[i]->style,"GROUP") == 0) + list[n++] = i; +} + /* ---------------------------------------------------------------------- create list of compute indices for computes which store invocation times ------------------------------------------------------------------------- */ diff --git a/src/modify.h b/src/modify.h index df8e52752f..da5d00a799 100644 --- a/src/modify.h +++ b/src/modify.h @@ -163,6 +163,9 @@ class Modify : protected Pointers { int *end_of_step_every; + int n_post_force_group; // list of fix GROUPs for post_force invocation + int *list_post_force_group; + int n_timeflag; // list of computes that store time invocation int *list_timeflag; @@ -187,6 +190,8 @@ class Modify : protected Pointers { void list_init_energy_couple(int &, int *&); void list_init_energy_global(int &, int *&); void list_init_energy_atom(int &, int *&); + void list_init_post_force_group(int &, int *&); + void list_init_post_force_respa_group(int &, int *&); void list_init_dofflag(int &, int *&); void list_init_compute(); diff --git a/src/update.cpp b/src/update.cpp index 3870f6d38e..297053e39d 100644 --- a/src/update.cpp +++ b/src/update.cpp @@ -60,7 +60,6 @@ Update::Update(LAMMPS *lmp) : Pointers(lmp) beginstep = endstep = 0; restrict_output = 0; setupflag = 0; - post_integrate = 0; multireplica = 0; eflag_global = vflag_global = -1; diff --git a/src/update.h b/src/update.h index afcafb87cb..c3e79b72a2 100644 --- a/src/update.h +++ b/src/update.h @@ -35,7 +35,6 @@ class Update : protected Pointers { int max_eval; // max force evaluations for minimizer int restrict_output; // 1 if output should not write dump/restart int setupflag; // set when setup() is computing forces - int post_integrate; // 1 if now at post_integrate() in timestep int multireplica; // 1 if min across replicas, else 0 int dt_default; // 1 if dt is at default value, else 0