From f89945e8e5c663ee4f4483eb3110f5059f298acb Mon Sep 17 00:00:00 2001
From: sjplimp <sjplimp@f3b2605a-c512-4ea7-a41b-209d697bcdaa>
Date: Thu, 14 Feb 2008 23:05:04 +0000
Subject: [PATCH] git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@1502
f3b2605a-c512-4ea7-a41b-209d697bcdaa
---
src/COLLOID/pair_colloid.cpp | 13 ++++++
src/COLLOID/pair_colloid.h | 3 +-
src/COLLOID/pair_lubricate.cpp | 2 -
src/comm.cpp | 67 +++++++++++++++++++++++-------
src/comm.h | 1 +
src/neigh_derive.cpp | 15 +++----
src/neigh_full.cpp | 12 +++---
src/neigh_gran.cpp | 22 +++++-----
src/neigh_half_bin.cpp | 12 +++---
src/neigh_half_multi.cpp | 12 +++---
src/neigh_half_nsq.cpp | 8 ++--
src/neigh_respa.cpp | 20 ++++-----
src/neighbor.cpp | 74 +++++++++++++++++++++-------------
src/neighbor.h | 3 ++
src/pair_lj_cut.cpp | 2 +-
15 files changed, 166 insertions(+), 100 deletions(-)
diff --git a/src/COLLOID/pair_colloid.cpp b/src/COLLOID/pair_colloid.cpp
index 59c6cc69d5..1e47215ae4 100644
--- a/src/COLLOID/pair_colloid.cpp
+++ b/src/COLLOID/pair_colloid.cpp
@@ -53,6 +53,7 @@ PairColloid::~PairColloid()
memory->destroy_2d_double_array(d2);
memory->destroy_2d_double_array(a1);
memory->destroy_2d_double_array(a2);
+ memory->destroy_2d_double_array(diameter);
memory->destroy_2d_double_array(cut);
memory->destroy_2d_double_array(offset);
memory->destroy_2d_double_array(sigma3);
@@ -229,6 +230,7 @@ void PairColloid::allocate()
d2 = memory->create_2d_double_array(n+1,n+1,"pair:d2");
a1 = memory->create_2d_double_array(n+1,n+1,"pair:a1");
a2 = memory->create_2d_double_array(n+1,n+1,"pair:a2");
+ diameter = memory->create_2d_double_array(n+1,n+1,"pair:diameter");
cut = memory->create_2d_double_array(n+1,n+1,"pair:cut");
offset = memory->create_2d_double_array(n+1,n+1,"pair:offset");
sigma3 = memory->create_2d_double_array(n+1,n+1,"pair:sigma3");
@@ -292,6 +294,7 @@ void PairColloid::coeff(int narg, char **arg)
error->all("Invalid d1 or d2 value for pair colloid coeff");
d1[i][j] = d1_one;
d2[i][j] = d2_one;
+ diameter[i][j] = 0.5*(d1_one+d2_one);
cut[i][j] = cut_one;
setflag[i][j] = 1;
count++;
@@ -312,6 +315,7 @@ double PairColloid::init_one(int i, int j)
sigma[i][j] = mix_distance(sigma[i][i],sigma[j][j]);
d1[i][j] = mix_distance(d1[i][i],d1[j][j]);
d2[i][j] = mix_distance(d2[i][i],d2[j][j]);
+ diameter[i][j] = 0.5 * (d1[i][j] + d2[i][j]);
cut[i][j] = mix_distance(cut[i][i],cut[j][j]);
}
@@ -340,6 +344,7 @@ double PairColloid::init_one(int i, int j)
sigma[j][i] = sigma[i][j];
sigma3[j][i] = sigma3[i][j];
sigma6[j][i] = sigma6[i][j];
+ diameter[j][i] = diameter[i][j];
cut[j][i] = cut[i][j];
cutsq[j][i] = cutsq[i][j] = cut[i][j] * cut[i][j];
@@ -516,3 +521,11 @@ double PairColloid::single(int i, int j, int itype, int jtype, double rsq,
return factor_lj*phi;
}
+
+/* ---------------------------------------------------------------------- */
+
+void *PairColloid::extract(char *str)
+{
+ if (strcmp(str,"diameter") == 0) return (void *) diameter;
+ return NULL;
+}
diff --git a/src/COLLOID/pair_colloid.h b/src/COLLOID/pair_colloid.h
index 87e10ef861..00f41b4b36 100644
--- a/src/COLLOID/pair_colloid.h
+++ b/src/COLLOID/pair_colloid.h
@@ -31,11 +31,12 @@ class PairColloid : public Pair {
void write_restart_settings(FILE *);
void read_restart_settings(FILE *);
double single(int, int, int, int, double, double, double, double &);
+ void *extract(char *);
private:
double cut_global;
double **cut;
- double **a12,**d1,**d2,**a1,**a2,**offset;
+ double **a12,**d1,**d2,**diameter,**a1,**a2,**offset;
double **sigma,**sigma3,**sigma6;
double **lj1,**lj2,**lj3,**lj4;
int **form;
diff --git a/src/COLLOID/pair_lubricate.cpp b/src/COLLOID/pair_lubricate.cpp
index b2323d2c93..78c7088fa8 100644
--- a/src/COLLOID/pair_lubricate.cpp
+++ b/src/COLLOID/pair_lubricate.cpp
@@ -352,8 +352,6 @@ void PairLubricate::init_style()
{
if (!atom->torque_flag || atom->shape == NULL)
error->all("Pair lubricate requires atom attributes torque, shape");
- if (domain->dimension != 3)
- error->all("Pair lubricate only available for 3d");
// insure all particle shapes are spherical
diff --git a/src/comm.cpp b/src/comm.cpp
index 414f70c152..51f74e9b0f 100644
--- a/src/comm.cpp
+++ b/src/comm.cpp
@@ -29,6 +29,7 @@
#include "neighbor.h"
#include "modify.h"
#include "fix.h"
+#include "group.h"
#include "compute.h"
#include "error.h"
#include "memory.h"
@@ -57,6 +58,7 @@ Comm::Comm(LAMMPS *lmp) : Pointers(lmp)
user_procgrid[0] = user_procgrid[1] = user_procgrid[2] = 0;
grid2proc = NULL;
+ igroup = 0;
style = SINGLE;
multilo = multihi = NULL;
cutghostmulti = NULL;
@@ -170,6 +172,7 @@ void Comm::init()
{
triclinic = domain->triclinic;
map_style = atom->map_style;
+ groupbit = group->bitmask[igroup];
// comm_only = 1 if only x,f are exchanged in forward/reverse comm
@@ -605,7 +608,7 @@ void Comm::borders()
{
int i,n,itype,iswap,dim,ineed,maxneed,nsend,nrecv,nfirst,nlast,smax,rmax;
double lo,hi;
- int *type;
+ int *type,*mask;
double **x;
double *buf,*mlo,*mhi;
MPI_Request request;
@@ -648,22 +651,47 @@ void Comm::borders()
}
nsend = 0;
- if (style == SINGLE) {
- for (i = nfirst; i < nlast; i++)
- if (x[i][dim] >= lo && x[i][dim] <= hi) {
- if (nsend == maxsendlist[iswap]) grow_list(iswap,nsend);
- sendlist[iswap][nsend++] = i;
- }
+
+ // SINGLE vs MULTI, all atoms versus only atoms in group
+
+ if (igroup) {
+ mask = atom->mask;
+ if (style == SINGLE) {
+ for (i = nfirst; i < nlast; i++)
+ if (mask[i] & groupbit) {
+ if (x[i][dim] >= lo && x[i][dim] <= hi) {
+ if (nsend == maxsendlist[iswap]) grow_list(iswap,nsend);
+ sendlist[iswap][nsend++] = i;
+ }
+ }
+ } else {
+ for (i = nfirst; i < nlast; i++)
+ if (mask[i] & groupbit) {
+ itype = type[i];
+ if (x[i][dim] >= mlo[itype] && x[i][dim] <= mhi[itype]) {
+ if (nsend == maxsendlist[iswap]) grow_list(iswap,nsend);
+ sendlist[iswap][nsend++] = i;
+ }
+ }
+ }
} else {
- for (i = nfirst; i < nlast; i++) {
- itype = type[i];
- if (x[i][dim] >= mlo[itype] && x[i][dim] <= mhi[itype]) {
- if (nsend == maxsendlist[iswap]) grow_list(iswap,nsend);
- sendlist[iswap][nsend++] = i;
+ if (style == SINGLE) {
+ for (i = nfirst; i < nlast; i++)
+ if (x[i][dim] >= lo && x[i][dim] <= hi) {
+ if (nsend == maxsendlist[iswap]) grow_list(iswap,nsend);
+ sendlist[iswap][nsend++] = i;
+ }
+ } else {
+ for (i = nfirst; i < nlast; i++) {
+ itype = type[i];
+ if (x[i][dim] >= mlo[itype] && x[i][dim] <= mhi[itype]) {
+ if (nsend == maxsendlist[iswap]) grow_list(iswap,nsend);
+ sendlist[iswap][nsend++] = i;
+ }
}
}
}
-
+
// pack up list of border atoms
if (nsend*size_border > maxsend)
@@ -1466,11 +1494,22 @@ void Comm::free_multi()
void Comm::set(int narg, char **arg)
{
- if (narg != 1) error->all("Illegal communicate command");
+ if (narg < 1) error->all("Illegal communicate command");
if (strcmp(arg[0],"single") == 0) style = SINGLE;
else if (strcmp(arg[0],"multi") == 0) style = MULTI;
else error->all("Illegal communicate command");
+
+ int iarg = 1;
+ while (iarg < narg) {
+ if (strcmp(arg[iarg],"group") == 0) {
+ if (iarg+2 > narg) error->all("Illegal communicate command");
+ igroup = group->find(arg[iarg+1]);
+ if (igroup < 0)
+ error->all("Invalid group ID in communicate command");
+ iarg += 2;
+ } else error->all("Illegal communicate command");
+ }
}
/* ----------------------------------------------------------------------
diff --git a/src/comm.h b/src/comm.h
index 09cf8046ef..b6b8ef0576 100644
--- a/src/comm.h
+++ b/src/comm.h
@@ -71,6 +71,7 @@ class Comm : protected Pointers {
int comm_x_only,comm_f_only; // 1 if only exchange x,f in for/rev comm
int map_style; // non-0 if global->local mapping is done
int ***grid2proc; // which proc owns i,j,k loc in 3d grid
+ int igroup,groupbit; // only communicate this group
int *firstrecv; // where to put 1st recv atom in each swap
int **sendlist; // list of atoms to send in each swap
diff --git a/src/neigh_derive.cpp b/src/neigh_derive.cpp
index d2c11f1292..cb1a2a6709 100644
--- a/src/neigh_derive.cpp
+++ b/src/neigh_derive.cpp
@@ -67,10 +67,9 @@ void Neighbor::half_from_full_no_newton(NeighList *list)
if (j > i) neighptr[n++] = j;
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -145,10 +144,9 @@ void Neighbor::half_from_full_newton(NeighList *list)
neighptr[n++] = joriginal;
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -217,10 +215,9 @@ void Neighbor::skip_from(NeighList *list)
neighptr[n++] = joriginal;
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -311,12 +308,11 @@ void Neighbor::skip_from_granular(NeighList *list)
shearptr[nn++] = shearptr_skip[3*jj+2];
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
firsttouch[i] = touchptr;
firstshear[i] = shearptr;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -454,10 +450,9 @@ void Neighbor::skip_from_respa(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
diff --git a/src/neigh_full.cpp b/src/neigh_full.cpp
index 878ab0819c..81c6c9586a 100644
--- a/src/neigh_full.cpp
+++ b/src/neigh_full.cpp
@@ -47,6 +47,7 @@ void Neighbor::full_nsq(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -83,10 +84,9 @@ void Neighbor::full_nsq(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -132,6 +132,7 @@ void Neighbor::full_bin(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -173,10 +174,9 @@ void Neighbor::full_bin(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -225,6 +225,7 @@ void Neighbor::full_multi(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -271,10 +272,9 @@ void Neighbor::full_multi(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
diff --git a/src/neigh_gran.cpp b/src/neigh_gran.cpp
index 915caa8e49..8adbc4501c 100644
--- a/src/neigh_gran.cpp
+++ b/src/neigh_gran.cpp
@@ -74,6 +74,7 @@ void Neighbor::granular_nsq_no_newton(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -142,15 +143,13 @@ void Neighbor::granular_nsq_no_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
if (fix_history) {
firsttouch[i] = touchptr;
firstshear[i] = shearptr;
}
-
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -194,6 +193,7 @@ void Neighbor::granular_nsq_newton(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -239,10 +239,9 @@ void Neighbor::granular_nsq_newton(NeighList *list)
if (rsq <= cutsq) neighptr[n++] = j;
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -314,6 +313,7 @@ void Neighbor::granular_bin_no_newton(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -389,15 +389,13 @@ void Neighbor::granular_bin_no_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
if (fix_history) {
firsttouch[i] = touchptr;
firstshear[i] = shearptr;
}
-
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -446,6 +444,7 @@ void Neighbor::granular_bin_newton(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -501,10 +500,9 @@ void Neighbor::granular_bin_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -553,6 +551,7 @@ void Neighbor::granular_bin_newton_tri(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -592,10 +591,9 @@ void Neighbor::granular_bin_newton_tri(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
diff --git a/src/neigh_half_bin.cpp b/src/neigh_half_bin.cpp
index c0e5ac7214..3f19e4cd6e 100644
--- a/src/neigh_half_bin.cpp
+++ b/src/neigh_half_bin.cpp
@@ -57,6 +57,7 @@ void Neighbor::half_bin_no_newton(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -100,10 +101,9 @@ void Neighbor::half_bin_no_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -150,6 +150,7 @@ void Neighbor::half_bin_newton(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -214,10 +215,9 @@ void Neighbor::half_bin_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -264,6 +264,7 @@ void Neighbor::half_bin_newton_tri(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -308,10 +309,9 @@ void Neighbor::half_bin_newton_tri(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
diff --git a/src/neigh_half_multi.cpp b/src/neigh_half_multi.cpp
index 8ecad243a2..cda4965ebc 100644
--- a/src/neigh_half_multi.cpp
+++ b/src/neigh_half_multi.cpp
@@ -60,6 +60,7 @@ void Neighbor::half_multi_no_newton(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -108,10 +109,9 @@ void Neighbor::half_multi_no_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -161,6 +161,7 @@ void Neighbor::half_multi_newton(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -232,10 +233,9 @@ void Neighbor::half_multi_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -285,6 +285,7 @@ void Neighbor::half_multi_newton_tri(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -334,10 +335,9 @@ void Neighbor::half_multi_newton_tri(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
diff --git a/src/neigh_half_nsq.cpp b/src/neigh_half_nsq.cpp
index ea7ae68e02..aaee116dcf 100644
--- a/src/neigh_half_nsq.cpp
+++ b/src/neigh_half_nsq.cpp
@@ -48,6 +48,7 @@ void Neighbor::half_nsq_no_newton(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -82,10 +83,9 @@ void Neighbor::half_nsq_no_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -125,6 +125,7 @@ void Neighbor::half_nsq_newton(NeighList *list)
int npnt = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -175,10 +176,9 @@ void Neighbor::half_nsq_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
diff --git a/src/neigh_respa.cpp b/src/neigh_respa.cpp
index 38bb3232e7..125bfecc99 100644
--- a/src/neigh_respa.cpp
+++ b/src/neigh_respa.cpp
@@ -68,6 +68,7 @@ void Neighbor::respa_nsq_no_newton(NeighList *list)
int npnt_middle = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -131,10 +132,9 @@ void Neighbor::respa_nsq_no_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -209,6 +209,7 @@ void Neighbor::respa_nsq_newton(NeighList *list)
int npnt_middle = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -288,10 +289,9 @@ void Neighbor::respa_nsq_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -373,6 +373,7 @@ void Neighbor::respa_bin_no_newton(NeighList *list)
int npnt_middle = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -445,10 +446,9 @@ void Neighbor::respa_bin_no_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -529,6 +529,7 @@ void Neighbor::respa_bin_newton(NeighList *list)
int npnt_middle = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -634,10 +635,9 @@ void Neighbor::respa_bin_newton(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
@@ -718,6 +718,7 @@ void Neighbor::respa_bin_newton_tri(NeighList *list)
int npnt_middle = 0;
for (i = 0; i < nlocal; i++) {
+ if (include_group && !(mask[i] & include_groupbit)) continue;
if (pgsize - npnt < oneatom) {
npnt = 0;
@@ -792,10 +793,9 @@ void Neighbor::respa_bin_newton_tri(NeighList *list)
}
}
- ilist[inum] = i;
+ ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
- inum++;
npnt += n;
if (npnt >= pgsize)
error->one("Neighbor list overflow, boost neigh_modify one or page");
diff --git a/src/neighbor.cpp b/src/neighbor.cpp
index b9af0647b3..0c0b6b2b84 100644
--- a/src/neighbor.cpp
+++ b/src/neighbor.cpp
@@ -90,6 +90,8 @@ Neighbor::Neighbor(LAMMPS *lmp) : Pointers(lmp)
// pair exclusion list info
+ include_group = 0;
+
nex_type = maxex_type = 0;
ex1_type = ex2_type = NULL;
ex_type = NULL;
@@ -338,6 +340,8 @@ void Neighbor::init()
n = atom->ntypes;
+ include_groupbit = group->bitmask[include_group];
+
if (nex_type == 0 && nex_group == 0 && nex_mol == 0) exclude = 0;
else exclude = 1;
@@ -920,15 +924,28 @@ int Neighbor::check_distance()
{
double delx,dely,delz,rsq;
- int nlocal = atom->nlocal;
double **x = atom->x;
+ int *mask = atom->mask;
+ int nlocal = atom->nlocal;
int flag = 0;
- for (int i = 0; i < nlocal; i++) {
- delx = x[i][0] - xhold[i][0];
- dely = x[i][1] - xhold[i][1];
- delz = x[i][2] - xhold[i][2];
- rsq = delx*delx + dely*dely + delz*delz;
- if (rsq > triggersq) flag = 1;
+
+ if (include_group) {
+ for (int i = 0; i < nlocal; i++)
+ if (mask[i] & include_groupbit) {
+ delx = x[i][0] - xhold[i][0];
+ dely = x[i][1] - xhold[i][1];
+ delz = x[i][2] - xhold[i][2];
+ rsq = delx*delx + dely*dely + delz*delz;
+ if (rsq > triggersq) flag = 1;
+ }
+ } else {
+ for (int i = 0; i < nlocal; i++) {
+ delx = x[i][0] - xhold[i][0];
+ dely = x[i][1] - xhold[i][1];
+ delz = x[i][2] - xhold[i][2];
+ rsq = delx*delx + dely*dely + delz*delz;
+ if (rsq > triggersq) flag = 1;
+ }
}
int flagall;
@@ -1288,6 +1305,12 @@ void Neighbor::modify_params(int narg, char **arg)
if (binsize_user <= 0.0) binsizeflag = 0;
else binsizeflag = 1;
iarg += 2;
+ } else if (strcmp(arg[iarg],"include") == 0) {
+ if (iarg+2 > narg) error->all("Illegal neigh_modify command");
+ include_group = group->find(arg[iarg+1]);
+ if (include_group == -1)
+ error->all("Invalid group ID in neigh_modify command");
+ iarg += 2;
} else if (strcmp(arg[iarg],"exclude") == 0) {
if (iarg+2 > narg) error->all("Illegal neigh_modify command");
@@ -1395,35 +1418,30 @@ int Neighbor::find_special(int i, int j)
void Neighbor::bin_atoms()
{
int i,ibin,nall;
- double **x;
-
- nall = atom->nlocal + atom->nghost;
- x = atom->x;
for (i = 0; i < mbins; i++) binhead[i] = -1;
// bin in reverse order so linked list will be in forward order
// also puts ghost atoms at end of list, which is necessary
- for (i = nall-1; i >= 0; i--) {
- ibin = coord2bin(x[i]);
- bins[i] = binhead[ibin];
- binhead[ibin] = i;
- }
-
- /*
- for (i = nlocal; i < nall; i++) {
- ibin = coord2bin(x[i]);
- bins[i] = binhead[ibin];
- binhead[ibin] = i;
- }
+ double **x = atom->x;
+ int *mask = atom->mask;
+ nall = atom->nlocal + atom->nghost;
- for (i = 0; i < nlocal; i++) {
- ibin = coord2bin(x[i]);
- bins[i] = binhead[ibin];
- binhead[ibin] = i;
+ if (include_group) {
+ for (i = nall-1; i >= 0; i--)
+ if (mask[i] & include_groupbit) {
+ ibin = coord2bin(x[i]);
+ bins[i] = binhead[ibin];
+ binhead[ibin] = i;
+ }
+ } else {
+ for (i = nall-1; i >= 0; i--) {
+ ibin = coord2bin(x[i]);
+ bins[i] = binhead[ibin];
+ binhead[ibin] = i;
+ }
}
- */
}
/* ----------------------------------------------------------------------
diff --git a/src/neighbor.h b/src/neighbor.h
index 39f5788986..3fc62c4421 100644
--- a/src/neighbor.h
+++ b/src/neighbor.h
@@ -123,6 +123,9 @@ class Neighbor : protected Pointers {
int special_flag[4]; // flags for 1-2, 1-3, 1-4 neighbors
+ int include_group; // only form neighbor lists on this group
+ int include_groupbit;
+
int exclude; // 0 if no type/group exclusions, 1 if yes
int nex_type; // # of entries in type exclusion list
diff --git a/src/pair_lj_cut.cpp b/src/pair_lj_cut.cpp
index 6ab2b0bbe6..24985e302b 100644
--- a/src/pair_lj_cut.cpp
+++ b/src/pair_lj_cut.cpp
@@ -717,6 +717,6 @@ double PairLJCut::single(int i, int j, int itype, int jtype, double rsq,
void *PairLJCut::extract(char *str)
{
- if (strcmp(str,"sigma") == 0) return (void *) sigma;
+ if (strcmp(str,"diameter") == 0) return (void *) sigma;
return NULL;
}
--
GitLab