diff --git a/doc/src/compute_coord_atom.txt b/doc/src/compute_coord_atom.txt index ddc4cc82d3..893efc6cf6 100644 --- a/doc/src/compute_coord_atom.txt +++ b/doc/src/compute_coord_atom.txt @@ -15,8 +15,9 @@ compute ID group-ID coord/atom cstyle args ... :pre ID, group-ID are documented in "compute"_compute.html command :ulb,l coord/atom = style name of this compute command :l cstyle = {cutoff} or {orientorder} :l - {cutoff} args = cutoff typeN + {cutoff} args = cutoff \[group group2-ID\] typeN cutoff = distance within which to count coordination neighbors (distance units) + group {group2-ID} = select group-ID to restrict which atoms to cound for coordination number (optional) typeN = atom type for Nth coordination count (see asterisk form below) {orientorder} args = orientorderID threshold orientorderID = ID of an orientorder/atom compute @@ -28,6 +29,7 @@ cstyle = {cutoff} or {orientorder} :l compute 1 all coord/atom cutoff 2.0 compute 1 all coord/atom cutoff 6.0 1 2 compute 1 all coord/atom cutoff 6.0 2*4 5*8 * +compute 1 solute coord/atom cutoff 2.0 group solvent compute 1 all coord/atom orientorder 2 0.5 :pre [Description:] @@ -38,9 +40,14 @@ meaning of the resulting value depend on the {cstyle} keyword used. The {cutoff} cstyle calculates one or more traditional coordination numbers for each atom. A coordination number is defined as the number -of neighbor atoms with specified atom type(s) that are within the -specified cutoff distance from the central atom. Atoms not in the -specified group are included in the coordination number tally. +of neighbor atoms with specified atom type(s), and optionally within +the specified group, that are within the specified cutoff distance from +the central atom. The compute group selects only the central atoms; all +neighboring atoms, unless selected by type, type range, or group option, +are included in the coordinations number tally. + +The optional {group} keyword allows to specify from which group atoms +contribute to the coordination number. Default setting is group 'all'. The {typeN} keywords allow specification of which atom types contribute to each coordination number. One coordination number is @@ -122,7 +129,9 @@ explained above. "compute cluster/atom"_compute_cluster_atom.html "compute orientorder/atom"_compute_orientorder_atom.html -[Default:] none +[Default:] + +group = all :line diff --git a/src/compute_coord_atom.cpp b/src/compute_coord_atom.cpp index 33b318ea17..021d092bc4 100644 --- a/src/compute_coord_atom.cpp +++ b/src/compute_coord_atom.cpp @@ -25,6 +25,7 @@ #include "force.h" #include "pair.h" #include "comm.h" +#include "group.h" #include "memory.h" #include "error.h" @@ -37,10 +38,12 @@ using namespace LAMMPS_NS; ComputeCoordAtom::ComputeCoordAtom(LAMMPS *lmp, int narg, char **arg) : Compute(lmp, narg, arg), typelo(NULL), typehi(NULL), cvec(NULL), carray(NULL), - id_orientorder(NULL), normv(NULL) + group2(NULL), id_orientorder(NULL), normv(NULL) { if (narg < 5) error->all(FLERR,"Illegal compute coord/atom command"); + jgroup = group->find("all"); + jgroupbit = group->bitmask[jgroup]; cstyle = NONE; if (strcmp(arg[3],"cutoff") == 0) { @@ -48,18 +51,29 @@ ComputeCoordAtom::ComputeCoordAtom(LAMMPS *lmp, int narg, char **arg) : double cutoff = force->numeric(FLERR,arg[4]); cutsq = cutoff*cutoff; - ncol = narg-5 + 1; + int iarg = 5; + if ((narg > 6) && (strcmp(arg[5],"group") == 0)) { + int len = strlen(arg[6])+1; + group2 = new char[len]; + strcpy(group2,arg[6]); + iarg += 2; + jgroup = group->find(group2); + if (jgroup == -1) + error->all(FLERR,"Compute coord/atom group2 ID does not exist"); + jgroupbit = group->bitmask[jgroup]; + } + + ncol = narg-iarg + 1; int ntypes = atom->ntypes; typelo = new int[ncol]; typehi = new int[ncol]; - if (narg == 5) { + if (narg == iarg) { ncol = 1; typelo[0] = 1; typehi[0] = ntypes; } else { ncol = 0; - int iarg = 5; while (iarg < narg) { force->bounds(FLERR,arg[iarg],ntypes,typelo[ncol],typehi[ncol]); if (typelo[ncol] > typehi[ncol]) @@ -106,6 +120,7 @@ ComputeCoordAtom::ComputeCoordAtom(LAMMPS *lmp, int narg, char **arg) : ComputeCoordAtom::~ComputeCoordAtom() { + delete [] group2; delete [] typelo; delete [] typehi; memory->destroy(cvec); @@ -229,13 +244,15 @@ 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 && jtype >= typelo[0] && jtype <= typehi[0]) - n++; + 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 && jtype >= typelo[0] && jtype <= typehi[0]) + n++; + } } cvec[i] = n; diff --git a/src/compute_coord_atom.h b/src/compute_coord_atom.h index 2bbc1b6720..0c1cc90f84 100644 --- a/src/compute_coord_atom.h +++ b/src/compute_coord_atom.h @@ -45,6 +45,9 @@ class ComputeCoordAtom : public Compute { double *cvec; double **carray; + char *group2; + int jgroup,jgroupbit; + class ComputeOrientOrderAtom *c_orientorder; char *id_orientorder; double threshold;