git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@3355 f3b2605a-c512-4ea7-a41b-209d697bcdaa

src/ASPHERE/atom_vec_ellipsoid.cpp

@@ -35,12 +35,15 @@ using namespace LAMMPS_NS;
 AtomVecEllipsoid::AtomVecEllipsoid(LAMMPS *lmp, int narg, char **arg) :
   AtomVec(lmp, narg, arg)
 {
+  molecular = 0;
   mass_type = 1;
   shape_type = 1;
+
   comm_x_only = comm_f_only = 0;
-  size_comm = 7;
+  size_forward = 7;
   size_reverse = 6;
   size_border = 10;
+  size_velocity = 6;
   size_data_atom = 9;
   size_data_vel = 7;
   xcol_data = 3;
@@ -158,6 +161,62 @@ int AtomVecEllipsoid::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecEllipsoid::pack_comm_vel(int n, int *list, double *buf,
+				    int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = quat[j][0];
+      buf[m++] = quat[j][1];
+      buf[m++] = quat[j][2];
+      buf[m++] = quat[j][3];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      buf[m++] = angmom[j][0];
+      buf[m++] = angmom[j][1];
+      buf[m++] = angmom[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = quat[j][0];
+      buf[m++] = quat[j][1];
+      buf[m++] = quat[j][2];
+      buf[m++] = quat[j][3];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      buf[m++] = angmom[j][0];
+      buf[m++] = angmom[j][1];
+      buf[m++] = angmom[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecEllipsoid::pack_comm_one(int i, double *buf)
 {
   buf[0] = quat[i][0];
@@ -188,6 +247,31 @@ void AtomVecEllipsoid::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecEllipsoid::unpack_comm_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    quat[i][0] = buf[m++];
+    quat[i][1] = buf[m++];
+    quat[i][2] = buf[m++];
+    quat[i][3] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+    angmom[i][0] = buf[m++];
+    angmom[i][1] = buf[m++];
+    angmom[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecEllipsoid::unpack_comm_one(int i, double *buf)
 {
   quat[i][0] = buf[0];
@@ -306,6 +390,68 @@ int AtomVecEllipsoid::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecEllipsoid::pack_border_vel(int n, int *list, double *buf,
+				      int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = quat[j][0];
+      buf[m++] = quat[j][1];
+      buf[m++] = quat[j][2];
+      buf[m++] = quat[j][3];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      buf[m++] = angmom[j][0];
+      buf[m++] = angmom[j][1];
+      buf[m++] = angmom[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = quat[j][0];
+      buf[m++] = quat[j][1];
+      buf[m++] = quat[j][2];
+      buf[m++] = quat[j][3];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      buf[m++] = angmom[j][0];
+      buf[m++] = angmom[j][1];
+      buf[m++] = angmom[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecEllipsoid::pack_border_one(int i, double *buf)
 {
   buf[0] = quat[i][0];
@@ -340,6 +486,35 @@ void AtomVecEllipsoid::unpack_border(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecEllipsoid::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    quat[i][0] = buf[m++];
+    quat[i][1] = buf[m++];
+    quat[i][2] = buf[m++];
+    quat[i][3] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+    angmom[i][0] = buf[m++];
+    angmom[i][1] = buf[m++];
+    angmom[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecEllipsoid::unpack_border_one(int i, double *buf)
 {
   quat[i][0] = buf[0];

src/ASPHERE/atom_vec_ellipsoid.h

@@ -25,16 +25,20 @@ class AtomVecEllipsoid : public AtomVec {
   void grow(int);
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   int pack_comm_one(int, double *);
   void unpack_comm(int, int, double *);
+  void unpack_comm_vel(int, int, double *);
   int unpack_comm_one(int, double *);
   int pack_reverse(int, int, double *);
   int pack_reverse_one(int, double *);
   void unpack_reverse(int, int *, double *);
   int unpack_reverse_one(int, double *);
   int pack_border(int, int *, double *, int, int *);
+  int pack_border_vel(int, int *, double *, int, int *);
   int pack_border_one(int, double *);
   void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int unpack_border_one(int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);

src/COLLOID/atom_vec_colloid.cpp

@@ -30,13 +30,16 @@ using namespace LAMMPS_NS;
 AtomVecColloid::AtomVecColloid(LAMMPS *lmp, int narg, char **arg) :
   AtomVec(lmp, narg, arg)
 {
+  molecular = 0;
   mass_type = 1;
   shape_type = 1;
-  comm_x_only = comm_f_only = 0;
+
-  ghost_velocity = 1;
+  comm_x_only = 1;
-  size_comm = 9;
+  comm_f_only = 0;
+  size_forward = 3;
   size_reverse = 6;
-  size_border = 12;
+  size_border = 6;
+  size_velocity = 6;
   size_data_atom = 5;
   size_data_vel = 7;
   xcol_data = 3;
@@ -110,6 +113,42 @@ int AtomVecColloid::pack_comm(int n, int *list, double *buf,
   int i,j,m;
   double dx,dy,dz;
 
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
+int AtomVecColloid::pack_comm_vel(int n, int *list, double *buf,
+				  int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
   m = 0;
   if (pbc_flag == 0) {
     for (i = 0; i < n; i++) {
@@ -152,20 +191,23 @@ int AtomVecColloid::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
-int AtomVecColloid::pack_comm_one(int i, double *buf)
+void AtomVecColloid::unpack_comm(int n, int first, double *buf)
 {
-  buf[0] = v[i][0];
+  int i,m,last;
-  buf[1] = v[i][1];
+
-  buf[2] = v[i][2];
+  m = 0;
-  buf[3] = omega[i][0];
+  last = first + n;
-  buf[4] = omega[i][1];
+  for (i = first; i < last; i++) {
-  buf[5] = omega[i][2];
+    x[i][0] = buf[m++];
-  return 6;
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+  }
 }
 
+
 /* ---------------------------------------------------------------------- */
 
-void AtomVecColloid::unpack_comm(int n, int first, double *buf)
+void AtomVecColloid::unpack_comm_vel(int n, int first, double *buf)
 {
   int i,m,last;
 
@@ -186,19 +228,6 @@ void AtomVecColloid::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
-int AtomVecColloid::unpack_comm_one(int i, double *buf)
-{
-  v[i][0] = buf[0];
-  v[i][1] = buf[1];
-  v[i][2] = buf[2];
-  omega[i][0] = buf[3];
-  omega[i][1] = buf[4];
-  omega[i][2] = buf[5];
-  return 6;
-}
-
-/* ---------------------------------------------------------------------- */
-
 int AtomVecColloid::pack_reverse(int n, int first, double *buf)
 {
   int i,m,last;
@@ -262,6 +291,48 @@ int AtomVecColloid::pack_border(int n, int *list, double *buf,
   int i,j,m;
   double dx,dy,dz;
 
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
+int AtomVecColloid::pack_border_vel(int n, int *list, double *buf,
+				    int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
   m = 0;
   if (pbc_flag == 0) {
     for (i = 0; i < n; i++) {
@@ -310,20 +381,26 @@ int AtomVecColloid::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
-int AtomVecColloid::pack_border_one(int i, double *buf)
+void AtomVecColloid::unpack_border(int n, int first, double *buf)
 {
-  buf[0] = v[i][0];
+  int i,m,last;
-  buf[1] = v[i][1];
+
-  buf[2] = v[i][2];
+  m = 0;
-  buf[3] = omega[i][0];
+  last = first + n;
-  buf[4] = omega[i][1];
+  for (i = first; i < last; i++) {
-  buf[5] = omega[i][2];
+    if (i == nmax) grow(0);
-  return 6;
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+  }
 }
 
 /* ---------------------------------------------------------------------- */
 
-void AtomVecColloid::unpack_border(int n, int first, double *buf)
+void AtomVecColloid::unpack_border_vel(int n, int first, double *buf)
 {
   int i,m,last;
 
@@ -346,19 +423,6 @@ void AtomVecColloid::unpack_border(int n, int first, double *buf)
   }
 }
 
-/* ---------------------------------------------------------------------- */
-
-int AtomVecColloid::unpack_border_one(int i, double *buf)
-{
-  v[i][0] = buf[0];
-  v[i][1] = buf[1];
-  v[i][2] = buf[2];
-  omega[i][0] = buf[3];
-  omega[i][1] = buf[4];
-  omega[i][2] = buf[5];
-  return 6;
-}
-
 /* ----------------------------------------------------------------------
    pack data for atom I for sending to another proc
    xyz must be 1st 3 values, so comm::exchange() can test on them

src/COLLOID/atom_vec_colloid.h

@@ -25,17 +25,17 @@ class AtomVecColloid : public AtomVec {
   void grow(int);
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
-  int pack_comm_one(int, double *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   void unpack_comm(int, int, double *);
-  int unpack_comm_one(int, double *);
+  void unpack_comm_vel(int, int, double *);
   int pack_reverse(int, int, double *);
   int pack_reverse_one(int, double *);
   void unpack_reverse(int, int *, double *);
   int unpack_reverse_one(int, double *);
   int pack_border(int, int *, double *, int, int *);
-  int pack_border_one(int, double *);
+  int pack_border_vel(int, int *, double *, int, int *);
   void unpack_border(int, int, double *);
-  int unpack_border_one(int, double *);
+  void unpack_border_vel(int, int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);
   int size_restart();

src/COLLOID/pair_lubricate.cpp

@@ -390,6 +390,8 @@ void PairLubricate::coeff(int narg, char **arg)
 
 void PairLubricate::init_style()
 {
+  if (comm->ghost_velocity == 0)
+    error->all("Pair lubricate requires ghost atoms store velocity");
   if (!atom->torque_flag || !atom->avec->shape_type)
     error->all("Pair lubricate requires atom attributes torque and shape");
   if (atom->radius_flag || atom->rmass_flag)

src/DIPOLE/atom_vec_dipole.cpp

@@ -29,13 +29,16 @@ using namespace LAMMPS_NS;
 AtomVecDipole::AtomVecDipole(LAMMPS *lmp, int narg, char **arg) :
   AtomVec(lmp, narg, arg)
 {
+  molecular = 0;
   mass_type = 1;
   shape_type = 1;
   dipole_type = 1;
+
   comm_x_only = comm_f_only = 0;
-  size_comm = 6;
+  size_forward = 6;
   size_reverse = 6;
   size_border = 10;
+  size_velocity = 6;
   size_data_atom = 9;
   size_data_vel = 7;
   xcol_data = 4;
@@ -153,6 +156,60 @@ int AtomVecDipole::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecDipole::pack_comm_vel(int n, int *list, double *buf,
+				 int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = mu[j][0];
+      buf[m++] = mu[j][1];
+      buf[m++] = mu[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      buf[m++] = omega[j][0];
+      buf[m++] = omega[j][1];
+      buf[m++] = omega[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = mu[j][0];
+      buf[m++] = mu[j][1];
+      buf[m++] = mu[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      buf[m++] = omega[j][0];
+      buf[m++] = omega[j][1];
+      buf[m++] = omega[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecDipole::pack_comm_one(int i, double *buf)
 {
   buf[0] = mu[i][0];
@@ -181,6 +238,30 @@ void AtomVecDipole::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecDipole::unpack_comm_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    mu[i][0] = buf[m++];
+    mu[i][1] = buf[m++];
+    mu[i][2] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+    omega[i][0] = buf[m++];
+    omega[i][1] = buf[m++];
+    omega[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecDipole::unpack_comm_one(int i, double *buf)
 {
   mu[i][0] = buf[0];
@@ -298,6 +379,68 @@ int AtomVecDipole::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecDipole::pack_border_vel(int n, int *list, double *buf,
+				   int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = q[j];
+      buf[m++] = mu[j][0];
+      buf[m++] = mu[j][1];
+      buf[m++] = mu[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      buf[m++] = omega[j][0];
+      buf[m++] = omega[j][1];
+      buf[m++] = omega[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = q[j];
+      buf[m++] = mu[j][0];
+      buf[m++] = mu[j][1];
+      buf[m++] = mu[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      buf[m++] = omega[j][0];
+      buf[m++] = omega[j][1];
+      buf[m++] = omega[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecDipole::pack_border_one(int i, double *buf)
 {
   buf[0] = q[i];
@@ -332,6 +475,35 @@ void AtomVecDipole::unpack_border(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecDipole::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    q[i] = buf[m++];
+    mu[i][0] = buf[m++];
+    mu[i][1] = buf[m++];
+    mu[i][2] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+    omega[i][0] = buf[m++];
+    omega[i][1] = buf[m++];
+    omega[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecDipole::unpack_border_one(int i, double *buf)
 {
   q[i] = buf[0];

src/DIPOLE/atom_vec_dipole.h

@@ -24,16 +24,20 @@ class AtomVecDipole : public AtomVec {
   void grow(int);
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   int pack_comm_one(int, double *);
   void unpack_comm(int, int, double *);
+  void unpack_comm_vel(int, int, double *);
   int unpack_comm_one(int, double *);
   int pack_reverse(int, int, double *);
   int pack_reverse_one(int, double *);
   void unpack_reverse(int, int *, double *);
   int unpack_reverse_one(int, double *);
   int pack_border(int, int *, double *, int, int *);
+  int pack_border_vel(int, int *, double *, int, int *);
   int pack_border_one(int, double *);
   void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int unpack_border_one(int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);

src/DPD/Install.csh

@@ -1,24 +0,0 @@
-# Install/unInstall package classes in LAMMPS
-
-if ($1 == 1) then
-
-  cp style_dpd.h ..
-
-  cp atom_vec_dpd.cpp ..
-  cp pair_dpd.cpp ..
-
-  cp atom_vec_dpd.h ..
-  cp pair_dpd.h ..
-
-else if ($1 == 0) then
-
-  rm ../style_dpd.h
-  touch ../style_dpd.h
-
-  rm ../atom_vec_dpd.cpp
-  rm ../pair_dpd.cpp
-
-  rm ../atom_vec_dpd.h
-  rm ../pair_dpd.h
-
-endif

src/DPD/atom_vec_dpd.cpp

@@ -1,211 +0,0 @@
-/* ----------------------------------------------------------------------
-   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
-   http://lammps.sandia.gov, Sandia National Laboratories
-   Steve Plimpton, sjplimp@sandia.gov
-
-   Copyright (2003) Sandia Corporation.  Under the terms of Contract
-   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
-   certain rights in this software.  This software is distributed under 
-   the GNU General Public License.
-
-   See the README file in the top-level LAMMPS directory.
-------------------------------------------------------------------------- */
-
-#include "stdlib.h"
-#include "atom_vec_dpd.h"
-#include "atom.h"
-#include "domain.h"
-#include "modify.h"
-#include "fix.h"
-#include "memory.h"
-#include "error.h"
-
-using namespace LAMMPS_NS;
-
-#define DELTA 10000
-
-/* ---------------------------------------------------------------------- */
-
-AtomVecDPD::AtomVecDPD(LAMMPS *lmp, int narg, char **arg) :
-  AtomVecAtomic(lmp, narg, arg)
-{
-  mass_type = 1;
-  comm_x_only = 0;
-  ghost_velocity = 1;
-  size_comm = 6;
-  size_reverse = 3;
-  size_border = 9;
-  size_data_atom = 5;
-  size_data_vel = 4;
-  xcol_data = 3;
-}
-
-/* ---------------------------------------------------------------------- */
-
-int AtomVecDPD::pack_comm(int n, int *list, double *buf,
-			  int pbc_flag, int *pbc)
-{
-  int i,j,m;
-  double dx,dy,dz;
-
-  m = 0;
-  if (pbc_flag == 0) {
-    for (i = 0; i < n; i++) {
-      j = list[i];
-      buf[m++] = x[j][0];
-      buf[m++] = x[j][1];
-      buf[m++] = x[j][2];
-      buf[m++] = v[j][0];
-      buf[m++] = v[j][1];
-      buf[m++] = v[j][2];
-    }
-  } else {
-    if (domain->triclinic == 0) {
-      dx = pbc[0]*domain->xprd;
-      dy = pbc[1]*domain->yprd;
-      dz = pbc[2]*domain->zprd;
-    } else {
-      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
-      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
-      dz = pbc[2]*domain->zprd;
-    }
-    for (i = 0; i < n; i++) {
-      j = list[i];
-      buf[m++] = x[j][0] + dx;
-      buf[m++] = x[j][1] + dy;
-      buf[m++] = x[j][2] + dz;
-      buf[m++] = v[j][0];
-      buf[m++] = v[j][1];
-      buf[m++] = v[j][2];
-    }
-  }
-  return m;
-}
-
-/* ---------------------------------------------------------------------- */
-
-int AtomVecDPD::pack_comm_one(int i, double *buf)
-{
-  buf[0] = v[i][0];
-  buf[1] = v[i][1];
-  buf[2] = v[i][2];
-  return 3;
-}
-
-/* ---------------------------------------------------------------------- */
-
-void AtomVecDPD::unpack_comm(int n, int first, double *buf)
-{
-  int i,m,last;
-
-  m = 0;
-  last = first + n;
-  for (i = first; i < last; i++) {
-    x[i][0] = buf[m++];
-    x[i][1] = buf[m++];
-    x[i][2] = buf[m++];
-    v[i][0] = buf[m++];
-    v[i][1] = buf[m++];
-    v[i][2] = buf[m++];
-  }
-}
-
-/* ---------------------------------------------------------------------- */
-
-int AtomVecDPD::unpack_comm_one(int i, double *buf)
-{
-  v[i][0] = buf[0];
-  v[i][1] = buf[1];
-  v[i][2] = buf[2];
-  return 3;
-}
-
-/* ---------------------------------------------------------------------- */
-
-int AtomVecDPD::pack_border(int n, int *list, double *buf,
-			    int pbc_flag, int *pbc)
-{
-  int i,j,m;
-  double dx,dy,dz;
-
-  m = 0;
-  if (pbc_flag == 0) {
-    for (i = 0; i < n; i++) {
-      j = list[i];
-      buf[m++] = x[j][0];
-      buf[m++] = x[j][1];
-      buf[m++] = x[j][2];
-      buf[m++] = tag[j];
-      buf[m++] = type[j];
-      buf[m++] = mask[j];
-      buf[m++] = v[j][0];
-      buf[m++] = v[j][1];
-      buf[m++] = v[j][2];
-   }
-  } else {
-    if (domain->triclinic == 0) {
-      dx = pbc[0]*domain->xprd;
-      dy = pbc[1]*domain->yprd;
-      dz = pbc[2]*domain->zprd;
-    } else {
-      dx = pbc[0];
-      dy = pbc[1];
-      dz = pbc[2];
-    }
-    for (i = 0; i < n; i++) {
-      j = list[i];
-      buf[m++] = x[j][0] + dx;
-      buf[m++] = x[j][1] + dy;
-      buf[m++] = x[j][2] + dz;
-      buf[m++] = tag[j];
-      buf[m++] = type[j];
-      buf[m++] = mask[j];
-      buf[m++] = v[j][0];
-      buf[m++] = v[j][1];
-      buf[m++] = v[j][2];
-    }
-  }
-  return m;
-}
-
-/* ---------------------------------------------------------------------- */
-
-int AtomVecDPD::pack_border_one(int i, double *buf)
-{
-  buf[0] = v[i][0];
-  buf[1] = v[i][1];
-  buf[2] = v[i][2];
-  return 3;
-}
-
-/* ---------------------------------------------------------------------- */
-
-void AtomVecDPD::unpack_border(int n, int first, double *buf)
-{
-  int i,m,last;
-
-  m = 0;
-  last = first + n;
-  for (i = first; i < last; i++) {
-    if (i == nmax) grow(0);
-    x[i][0] = buf[m++];
-    x[i][1] = buf[m++];
-    x[i][2] = buf[m++];
-    tag[i] = static_cast<int> (buf[m++]);
-    type[i] = static_cast<int> (buf[m++]);
-    mask[i] = static_cast<int> (buf[m++]);
-    v[i][0] = buf[m++];
-    v[i][1] = buf[m++];
-    v[i][2] = buf[m++];
-  }
-}
-
-/* ---------------------------------------------------------------------- */
-
-int AtomVecDPD::unpack_border_one(int i, double *buf)
-{
-  v[i][0] = buf[0];
-  v[i][1] = buf[1];
-  v[i][2] = buf[2];
-  return 3;
-}

src/DPD/atom_vec_dpd.h

@@ -1,37 +0,0 @@
-/* ----------------------------------------------------------------------
-   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
-   http://lammps.sandia.gov, Sandia National Laboratories
-   Steve Plimpton, sjplimp@sandia.gov
-
-   Copyright (2003) Sandia Corporation.  Under the terms of Contract
-   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
-   certain rights in this software.  This software is distributed under 
-   the GNU General Public License.
-
-   See the README file in the top-level LAMMPS directory.
-------------------------------------------------------------------------- */
-
-#ifndef ATOM_VEC_DPD_H
-#define ATOM_VEC_DPD_H
-
-#include "atom_vec_atomic.h"
-
-namespace LAMMPS_NS {
-
-class AtomVecDPD : public AtomVecAtomic {
- public:
-  AtomVecDPD(class LAMMPS *, int, char **);
-  void zero_ghost(int, int);
-  int pack_comm(int, int *, double *, int, int *);
-  int pack_comm_one(int, double *);
-  void unpack_comm(int, int, double *);
-  int unpack_comm_one(int, double *);
-  int pack_border(int, int *, double *, int, int *);
-  int pack_border_one(int, double *);
-  void unpack_border(int, int, double *);
-  int unpack_border_one(int, double *);
-};
-
-}
-
-#endif

src/DPD/style_dpd.h

@@ -1,28 +0,0 @@
-/* ----------------------------------------------------------------------
-   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
-   http://lammps.sandia.gov, Sandia National Laboratories
-   Steve Plimpton, sjplimp@sandia.gov
-
-   Copyright (2003) Sandia Corporation.  Under the terms of Contract
-   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
-   certain rights in this software.  This software is distributed under 
-   the GNU General Public License.
-
-   See the README file in the top-level LAMMPS directory.
-------------------------------------------------------------------------- */
-
-#ifdef AtomInclude
-#include "atom_vec_dpd.h"
-#endif
-
-#ifdef AtomClass
-AtomStyle(dpd,AtomVecDPD)
-#endif
-
-#ifdef PairInclude
-#include "pair_dpd.h"
-#endif
-
-#ifdef PairClass
-PairStyle(dpd,PairDPD)
-#endif

src/GRANULAR/atom_vec_granular.cpp

@@ -31,11 +31,14 @@ using namespace LAMMPS_NS;
 AtomVecGranular::AtomVecGranular(LAMMPS *lmp, int narg, char **arg) :
   AtomVec(lmp, narg, arg)
 {
-  comm_x_only = comm_f_only = 0;
+  molecular = 0;
-  ghost_velocity = 1;
+
-  size_comm = 9;
+  comm_x_only = 1;
+  comm_f_only = 0;
+  size_forward = 3;
   size_reverse = 6;
-  size_border = 14;
+  size_border = 8;
+  size_velocity = 6;
   size_data_atom = 7;
   size_data_vel = 7;
   xcol_data = 5;
@@ -122,6 +125,42 @@ int AtomVecGranular::pack_comm(int n, int *list, double *buf,
   int i,j,m;
   double dx,dy,dz;
 
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
+int AtomVecGranular::pack_comm_vel(int n, int *list, double *buf,
+				   int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
   m = 0;
   if (pbc_flag == 0) {
     for (i = 0; i < n; i++) {
@@ -164,20 +203,22 @@ int AtomVecGranular::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
-int AtomVecGranular::pack_comm_one(int i, double *buf)
+void AtomVecGranular::unpack_comm(int n, int first, double *buf)
 {
-  buf[0] = v[i][0];
+  int i,m,last;
-  buf[1] = v[i][1];
+
-  buf[2] = v[i][2];
+  m = 0;
-  buf[3] = omega[i][0];
+  last = first + n;
-  buf[4] = omega[i][1];
+  for (i = first; i < last; i++) {
-  buf[5] = omega[i][2];
+    x[i][0] = buf[m++];
-  return 6;
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+  }
 }
 
 /* ---------------------------------------------------------------------- */
 
-void AtomVecGranular::unpack_comm(int n, int first, double *buf)
+void AtomVecGranular::unpack_comm_vel(int n, int first, double *buf)
 {
   int i,m,last;
 
@@ -198,19 +239,6 @@ void AtomVecGranular::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
-int AtomVecGranular::unpack_comm_one(int i, double *buf)
-{
-  v[i][0] = buf[0];
-  v[i][1] = buf[1];
-  v[i][2] = buf[2];
-  omega[i][0] = buf[3];
-  omega[i][1] = buf[4];
-  omega[i][2] = buf[5];
-  return 6;
-}
-
-/* ---------------------------------------------------------------------- */
-
 int AtomVecGranular::pack_reverse(int n, int first, double *buf)
 {
   int i,m,last;
@@ -284,11 +312,57 @@ int AtomVecGranular::pack_border(int n, int *list, double *buf,
       buf[m++] = tag[j];
       buf[m++] = type[j];
       buf[m++] = mask[j];
+      buf[m++] = radius[j];
+      buf[m++] = rmass[j];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = radius[j];
+      buf[m++] = rmass[j];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
+int AtomVecGranular::pack_border_vel(int n, int *list, double *buf,
+				     int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = radius[j];
+      buf[m++] = rmass[j];
       buf[m++] = v[j][0];
       buf[m++] = v[j][1];
       buf[m++] = v[j][2];
-      buf[m++] = radius[j];
-      buf[m++] = rmass[j];
       buf[m++] = omega[j][0];
       buf[m++] = omega[j][1];
       buf[m++] = omega[j][2];
@@ -311,11 +385,11 @@ int AtomVecGranular::pack_border(int n, int *list, double *buf,
       buf[m++] = tag[j];
       buf[m++] = type[j];
       buf[m++] = mask[j];
+      buf[m++] = radius[j];
+      buf[m++] = rmass[j];
       buf[m++] = v[j][0];
       buf[m++] = v[j][1];
       buf[m++] = v[j][2];
-      buf[m++] = radius[j];
-      buf[m++] = rmass[j];
       buf[m++] = omega[j][0];
       buf[m++] = omega[j][1];
       buf[m++] = omega[j][2];
@@ -328,15 +402,9 @@ int AtomVecGranular::pack_border(int n, int *list, double *buf,
 
 int AtomVecGranular::pack_border_one(int i, double *buf)
 {
-  buf[0] = v[i][0];
+  buf[0] = radius[i];
-  buf[1] = v[i][1];
+  buf[1] = rmass[i];
-  buf[2] = v[i][2];
+  return 2;
-  buf[3] = radius[i];
-  buf[4] = rmass[i];
-  buf[5] = omega[i][0];
-  buf[6] = omega[i][1];
-  buf[7] = omega[i][2];
-  return 8;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -355,11 +423,33 @@ void AtomVecGranular::unpack_border(int n, int first, double *buf)
     tag[i] = static_cast<int> (buf[m++]);
     type[i] = static_cast<int> (buf[m++]);
     mask[i] = static_cast<int> (buf[m++]);
+    radius[i] = buf[m++];
+    rmass[i] = buf[m++];
+  }
+}
+
+
+/* ---------------------------------------------------------------------- */
+
+void AtomVecGranular::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    radius[i] = buf[m++];
+    rmass[i] = buf[m++];
     v[i][0] = buf[m++];
     v[i][1] = buf[m++];
     v[i][2] = buf[m++];
-    radius[i] = buf[m++];
-    rmass[i] = buf[m++];
     omega[i][0] = buf[m++];
     omega[i][1] = buf[m++];
     omega[i][2] = buf[m++];
@@ -370,15 +460,9 @@ void AtomVecGranular::unpack_border(int n, int first, double *buf)
 
 int AtomVecGranular::unpack_border_one(int i, double *buf)
 {
-  v[i][0] = buf[0];
+  radius[i] = buf[0];
-  v[i][1] = buf[1];
+  rmass[i] = buf[1];
-  v[i][2] = buf[2];
+  return 2;
-  radius[i] = buf[3];
-  rmass[i] = buf[4];
-  omega[i][0] = buf[5];
-  omega[i][1] = buf[6];
-  omega[i][2] = buf[7];
-  return 8;
 }
 
 /* ----------------------------------------------------------------------

src/GRANULAR/atom_vec_granular.h

@@ -25,16 +25,18 @@ class AtomVecGranular : public AtomVec {
   void grow(int);
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
-  int pack_comm_one(int, double *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   void unpack_comm(int, int, double *);
-  int unpack_comm_one(int, double *);
+  void unpack_comm_vel(int, int, double *);
   int pack_reverse(int, int, double *);
   int pack_reverse_one(int, double *);
   void unpack_reverse(int, int *, double *);
   int unpack_reverse_one(int, double *);
   int pack_border(int, int *, double *, int, int *);
+  int pack_border_vel(int, int *, double *, int, int *);
   int pack_border_one(int, double *);
   void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int unpack_border_one(int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);

src/GRANULAR/pair_gran_hooke_history.cpp

@@ -343,7 +343,7 @@ void PairGranHookeHistory::init_style()
 
   if (!atom->radius_flag || !atom->omega_flag || !atom->torque_flag)
     error->all("Pair granular requires atom attributes radius, omega, torque");
-  if (atom->avec->ghost_velocity == 0)
+  if (comm->ghost_velocity == 0)
     error->all("Pair granular requires ghost atoms store velocity");
 
   // need a half neigh list and optionally a granular history neigh list

src/MOLECULE/atom_vec_angle.cpp

@@ -33,12 +33,14 @@ AtomVecAngle::AtomVecAngle(LAMMPS *lmp, int narg, char **arg) :
   AtomVec(lmp, narg, arg)
 {
   molecular = 1;
-  bonds_allow = 1;
+  bonds_allow = angles_allow = 1;
-  angles_allow = 1;
   mass_type = 1;
-  size_comm = 3;
+
+  comm_x_only = comm_f_only = 1;
+  size_forward = 3;
   size_reverse = 3;
   size_border = 7;
+  size_velocity = 3;
   size_data_atom = 6;
   size_data_vel = 4;
   xcol_data = 4;
@@ -166,6 +168,48 @@ int AtomVecAngle::pack_comm(int n, int *list, double *buf,
   int i,j,m;
   double dx,dy,dz;
 
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
+int AtomVecAngle::pack_comm_vel(int n, int *list, double *buf,
+				int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
   m = 0;
   if (pbc_flag == 0) {
     for (i = 0; i < n; i++) {
@@ -211,6 +255,24 @@ void AtomVecAngle::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecAngle::unpack_comm_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecAngle::pack_reverse(int n, int first, double *buf)
 {
   int i,m,last;
@@ -286,6 +348,56 @@ int AtomVecAngle::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecAngle::pack_border_vel(int n, int *list, double *buf,
+				  int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = molecule[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = molecule[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecAngle::pack_border_one(int i, double *buf)
 {
   buf[0] = molecule[i];
@@ -314,6 +426,29 @@ void AtomVecAngle::unpack_border(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecAngle::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    molecule[i] = static_cast<int> (buf[m++]);
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecAngle::unpack_border_one(int i, double *buf)
 {
   molecule[i] = static_cast<int> (buf[0]);

src/MOLECULE/atom_vec_angle.h

@@ -25,12 +25,16 @@ class AtomVecAngle : public AtomVec {
   void reset_special();
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   void unpack_comm(int, int, double *);
+  void unpack_comm_vel(int, int, double *);
   int pack_reverse(int, int, double *);
   void unpack_reverse(int, int *, double *);
   int pack_border(int, int *, double *, int, int *);
+  int pack_border_vel(int, int *, double *, int, int *);
   int pack_border_one(int, double *);
   void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int unpack_border_one(int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);

src/MOLECULE/atom_vec_bond.cpp

@@ -35,9 +35,12 @@ AtomVecBond::AtomVecBond(LAMMPS *lmp, int narg, char **arg) :
   molecular = 1;
   bonds_allow = 1;
   mass_type = 1;
-  size_comm = 3;
+
+  comm_x_only = comm_f_only = 1;
+  size_forward = 3;
   size_reverse = 3;
   size_border = 7;
+  size_velocity = 3;
   size_data_atom = 6;
   size_data_vel = 4;
   xcol_data = 4;
@@ -172,6 +175,48 @@ int AtomVecBond::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecBond::pack_comm_vel(int n, int *list, double *buf,
+			       int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 void AtomVecBond::unpack_comm(int n, int first, double *buf)
 {
   int i,m,last;
@@ -187,6 +232,24 @@ void AtomVecBond::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecBond::unpack_comm_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecBond::pack_reverse(int n, int first, double *buf)
 {
   int i,m,last;
@@ -262,6 +325,56 @@ int AtomVecBond::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecBond::pack_border_vel(int n, int *list, double *buf,
+				 int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = molecule[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = molecule[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecBond::pack_border_one(int i, double *buf)
 {
   buf[0] = molecule[i];
@@ -290,6 +403,29 @@ void AtomVecBond::unpack_border(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecBond::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    molecule[i] = static_cast<int> (buf[m++]);
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecBond::unpack_border_one(int i, double *buf)
 {
   molecule[i] = static_cast<int> (buf[0]);

src/MOLECULE/atom_vec_bond.h

@@ -25,12 +25,16 @@ class AtomVecBond : public AtomVec {
   void reset_special();
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   void unpack_comm(int, int, double *);
+  void unpack_comm_vel(int, int, double *);
   int pack_reverse(int, int, double *);
   void unpack_reverse(int, int *, double *);
   int pack_border(int, int *, double *, int, int *);
+  int pack_border_vel(int, int *, double *, int, int *);
   int pack_border_one(int, double *);
   void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int unpack_border_one(int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);

src/MOLECULE/atom_vec_full.cpp

@@ -33,14 +33,14 @@ AtomVecFull::AtomVecFull(LAMMPS *lmp, int narg, char **arg) :
   AtomVec(lmp, narg, arg)
 {
   molecular = 1;
-  bonds_allow = 1;
+  bonds_allow = angles_allow = dihedrals_allow = impropers_allow = 1;
-  angles_allow = 1;
-  dihedrals_allow = 1;
-  impropers_allow = 1;
   mass_type = 1;
-  size_comm = 3;
+
+  comm_x_only = comm_f_only = 1;
+  size_forward = 3;
   size_reverse = 3;
   size_border = 8;
+  size_velocity = 3;
   size_data_atom = 7;
   size_data_vel = 4;
   xcol_data = 5;
@@ -255,6 +255,48 @@ int AtomVecFull::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecFull::pack_comm_vel(int n, int *list, double *buf,
+			       int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 void AtomVecFull::unpack_comm(int n, int first, double *buf)
 {
   int i,m,last;
@@ -270,6 +312,24 @@ void AtomVecFull::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecFull::unpack_comm_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecFull::pack_reverse(int n, int first, double *buf)
 {
   int i,m,last;
@@ -347,6 +407,58 @@ int AtomVecFull::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecFull::pack_border_vel(int n, int *list, double *buf,
+				 int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = q[j];
+      buf[m++] = molecule[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = q[j];
+      buf[m++] = molecule[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecFull::pack_border_one(int i, double *buf)
 {
   buf[0] = q[i];
@@ -377,6 +489,30 @@ void AtomVecFull::unpack_border(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecFull::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    q[i] = buf[m++];
+    molecule[i] = static_cast<int> (buf[m++]);
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecFull::unpack_border_one(int i, double *buf)
 {
   q[i] = buf[0];

src/MOLECULE/atom_vec_full.h

@@ -25,12 +25,16 @@ class AtomVecFull : public AtomVec {
   void reset_special();
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   void unpack_comm(int, int, double *);
+  void unpack_comm_vel(int, int, double *);
   int pack_reverse(int, int, double *);
   void unpack_reverse(int, int *, double *);
   int pack_border(int, int *, double *, int, int *);
+  int pack_border_vel(int, int *, double *, int, int *);
   int pack_border_one(int, double *);
   void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int unpack_border_one(int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);

src/MOLECULE/atom_vec_molecular.cpp

@@ -33,14 +33,14 @@ AtomVecMolecular::AtomVecMolecular(LAMMPS *lmp, int narg, char **arg) :
   AtomVec(lmp, narg, arg)
 {
   molecular = 1;
-  bonds_allow = 1;
+  bonds_allow = angles_allow = dihedrals_allow = impropers_allow = 1;
-  angles_allow = 1;
-  dihedrals_allow = 1;
-  impropers_allow = 1;
   mass_type = 1;
-  size_comm = 3;
+
+  comm_x_only = comm_f_only = 1;
+  size_forward = 3;
   size_reverse = 3;
   size_border = 7;
+  size_velocity = 3;
   size_data_atom = 6;
   size_data_vel = 4;
   xcol_data = 4;
@@ -252,6 +252,48 @@ int AtomVecMolecular::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecMolecular::pack_comm_vel(int n, int *list, double *buf, 
+				    int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 void AtomVecMolecular::unpack_comm(int n, int first, double *buf)
 {
   int i,m,last;
@@ -267,6 +309,24 @@ void AtomVecMolecular::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecMolecular::unpack_comm_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecMolecular::pack_reverse(int n, int first, double *buf)
 {
   int i,m,last;
@@ -342,6 +402,56 @@ int AtomVecMolecular::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecMolecular::pack_border_vel(int n, int *list, double *buf,
+				      int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = molecule[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = molecule[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecMolecular::pack_border_one(int i, double *buf)
 {
   buf[0] = molecule[i];
@@ -370,6 +480,29 @@ void AtomVecMolecular::unpack_border(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecMolecular::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    molecule[i] = static_cast<int> (buf[m++]);
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecMolecular::unpack_border_one(int i, double *buf)
 {
   molecule[i] = static_cast<int> (buf[0]);

src/MOLECULE/atom_vec_molecular.h

@@ -25,12 +25,16 @@ class AtomVecMolecular : public AtomVec {
   void reset_special();
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   void unpack_comm(int, int, double *);
+  void unpack_comm_vel(int, int, double *);
   int pack_reverse(int, int, double *);
   void unpack_reverse(int, int *, double *);
   int pack_border(int, int *, double *, int, int *);
+  int pack_border_vel(int, int *, double *, int, int *);
   int pack_border_one(int, double *);
   void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int unpack_border_one(int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);

src/Makefile

@@ -13,7 +13,7 @@ OBJ = 	$(SRC:.cpp=.o)
 
 # Package variables
 
-PACKAGE = asphere class2 colloid dipole dpd dsmc gpu granular \
+PACKAGE = asphere class2 colloid dipole dsmc gpu granular \
 	  kspace manybody meam molecule opt peri poems prd reax xtc
 
 PACKUSER = user-ackland user-atc user-cd-eam user-cg-cmm user-ewaldn \

src/PERI/atom_vec_peri.cpp

@@ -32,12 +32,16 @@ using namespace LAMMPS_NS;
 /* ---------------------------------------------------------------------- */
 
 AtomVecPeri::AtomVecPeri(LAMMPS *lmp, int narg, char **arg) : 
-AtomVec(lmp, narg, arg) 
+  AtomVec(lmp, narg, arg) 
 {
+  molecular = 0;
+
   comm_x_only = 0;
-  size_comm = 4;
+  comm_f_only = 1;
+  size_forward = 4;
   size_reverse = 3;
   size_border = 11;
+  size_velocity = 3;
   size_data_atom = 7;
   size_data_vel = 4;
   xcol_data = 5;
@@ -153,6 +157,51 @@ int AtomVecPeri::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecPeri::pack_comm_vel(int n, int *list, double *buf,
+			       int pbc_flag, int *pbc)
+
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = s0[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = s0[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecPeri::pack_comm_one(int i, double *buf)
 {
   buf[0] = s0[i];
@@ -177,6 +226,25 @@ void AtomVecPeri::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecPeri::unpack_comm_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    s0[i] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecPeri::unpack_comm_one(int i, double *buf)
 {
   s0[i] = buf[0];
@@ -268,6 +336,64 @@ int AtomVecPeri::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecPeri::pack_border_vel(int n, int *list, double *buf,
+				 int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = vfrac[j];
+      buf[m++] = s0[j];
+      buf[m++] = x0[j][0];
+      buf[m++] = x0[j][1];
+      buf[m++] = x0[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = vfrac[j];
+      buf[m++] = s0[j];
+      buf[m++] = x0[j][0];
+      buf[m++] = x0[j][1];
+      buf[m++] = x0[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecPeri::pack_border_one(int i, double *buf)
 {
   buf[0] = vfrac[i];
@@ -304,6 +430,33 @@ void AtomVecPeri::unpack_border(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecPeri::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    vfrac[i] = buf[m++];
+    s0[i] = buf[m++];
+    x0[i][0] = buf[m++];
+    x0[i][1] = buf[m++];
+    x0[i][2] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+ 
 int AtomVecPeri::unpack_border_one(int i, double *buf)
 {
   vfrac[i] = buf[0];

src/PERI/atom_vec_peri.h

@@ -24,14 +24,18 @@ class AtomVecPeri : public AtomVec {
   void grow(int);
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   int pack_comm_one(int, double *);
   void unpack_comm(int, int, double *);
+  void unpack_comm_vel(int, int, double *);
   int unpack_comm_one(int, double *);
   int pack_reverse(int, int, double *);
   void unpack_reverse(int, int *, double *);
   int pack_border(int, int *, double *, int, int *);
+  int pack_border_vel(int, int *, double *, int, int *);
   int pack_border_one(int, double *);
   void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int unpack_border_one(int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);

src/atom_vec.cpp

@@ -22,12 +22,8 @@ using namespace LAMMPS_NS;
 AtomVec::AtomVec(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp)
 {
   nmax = 0;
-  molecular = 0;
   bonds_allow = angles_allow = dihedrals_allow = impropers_allow = 0;
   mass_type = shape_type = dipole_type = 0;
-  comm_x_only = comm_f_only = 1;
-  ghost_velocity = 0;
-  size_comm = size_reverse = size_border = 0;
 }
 
 /* ----------------------------------------------------------------------

src/atom_vec.h

@@ -26,13 +26,14 @@ class AtomVec : protected Pointers {
   int mass_type;                       // 1 if per-type masses
   int shape_type;                      // 1 if per-type shape array
   int dipole_type;                     // 1 if per-type dipole moments
+
   int comm_x_only;                     // 1 if only exchange x in forward comm
   int comm_f_only;                     // 1 if only exchange f in reverse comm
-  int ghost_velocity;                  // 1 if ghost atoms store velocity
 
-  int size_comm;                       // # of values per atom in comm       
+  int size_forward;                    // # of values per atom in comm       
   int size_reverse;                    // # in reverse comm
   int size_border;                     // # in border comm
+  int size_velocity;                   // # of velocity based quantities
   int size_data_atom;                  // number of values in Atom line
   int size_data_vel;                   // number of values in Velocity line
   int xcol_data;                       // column (1-N) where x is in Atom line
@@ -46,8 +47,10 @@ class AtomVec : protected Pointers {
   virtual void copy(int, int) = 0;
 
   virtual int pack_comm(int, int *, double *, int, int *) = 0;
+  virtual int pack_comm_vel(int, int *, double *, int, int *) = 0;
   virtual int pack_comm_one(int, double *) {return 0;}
   virtual void unpack_comm(int, int, double *) = 0;
+  virtual void unpack_comm_vel(int, int, double *) = 0;
   virtual int unpack_comm_one(int, double *) {return 0;}
 
   virtual int pack_reverse(int, int, double *) = 0;
@@ -56,8 +59,10 @@ class AtomVec : protected Pointers {
   virtual int unpack_reverse_one(int, double *) {return 0;}
 
   virtual int pack_border(int, int *, double *, int, int *) = 0;
+  virtual int pack_border_vel(int, int *, double *, int, int *) = 0;
   virtual int pack_border_one(int, double *) {return 0;}
   virtual void unpack_border(int, int, double *) = 0;
+  virtual void unpack_border_vel(int, int, double *) = 0;
   virtual int unpack_border_one(int, double *) {return 0;}
 
   virtual int pack_exchange(int, double *) = 0;

src/atom_vec_atomic.cpp

@@ -29,10 +29,14 @@ using namespace LAMMPS_NS;
 AtomVecAtomic::AtomVecAtomic(LAMMPS *lmp, int narg, char **arg) :
   AtomVec(lmp, narg, arg)
 {
+  molecular = 0;
   mass_type = 1;
-  size_comm = 3;
+
+  comm_x_only = comm_f_only = 1;
+  size_forward = 3;
   size_reverse = 3;
   size_border = 6;
+  size_velocity = 3;
   size_data_atom = 5;
   size_data_vel = 4;
   xcol_data = 3;
@@ -125,6 +129,48 @@ int AtomVecAtomic::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecAtomic::pack_comm_vel(int n, int *list, double *buf,
+				 int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 void AtomVecAtomic::unpack_comm(int n, int first, double *buf)
 {
   int i,m,last;
@@ -140,6 +186,24 @@ void AtomVecAtomic::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecAtomic::unpack_comm_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecAtomic::pack_reverse(int n, int first, double *buf)
 {
   int i,m,last;
@@ -213,6 +277,54 @@ int AtomVecAtomic::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecAtomic::pack_border_vel(int n, int *list, double *buf,
+				   int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 void AtomVecAtomic::unpack_border(int n, int first, double *buf)
 {
   int i,m,last;
@@ -230,6 +342,28 @@ void AtomVecAtomic::unpack_border(int n, int first, double *buf)
   }
 }
 
+/* ---------------------------------------------------------------------- */
+
+void AtomVecAtomic::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
 /* ----------------------------------------------------------------------
    pack data for atom I for sending to another proc
    xyz must be 1st 3 values, so comm::exchange() can test on them 

src/atom_vec_atomic.h

@@ -21,15 +21,19 @@ namespace LAMMPS_NS {
 class AtomVecAtomic : public AtomVec {
  public:
   AtomVecAtomic(class LAMMPS *, int, char **);
-  virtual ~AtomVecAtomic() {}
+  ~AtomVecAtomic() {}
   void grow(int);
   void copy(int, int);
-  virtual int pack_comm(int, int *, double *, int, int *);
+  int pack_comm(int, int *, double *, int, int *);
-  virtual void unpack_comm(int, int, double *);
+  int pack_comm_vel(int, int *, double *, int, int *);
+  void unpack_comm(int, int, double *);
+  void unpack_comm_vel(int, int, double *);
   int pack_reverse(int, int, double *);
   void unpack_reverse(int, int *, double *);
-  virtual int pack_border(int, int *, double *, int, int *);
+  int pack_border(int, int *, double *, int, int *);
-  virtual void unpack_border(int, int, double *);
+  int pack_border_vel(int, int *, double *, int, int *);
+  void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);
   int size_restart();
@@ -40,7 +44,7 @@ class AtomVecAtomic : public AtomVec {
   int data_atom_hybrid(int, char **);
   double memory_usage();
 
- protected:
+ private:
   int *tag,*type,*mask,*image;
   double **x,**v,**f;
 };

src/atom_vec_charge.cpp

@@ -29,10 +29,14 @@ using namespace LAMMPS_NS;
 AtomVecCharge::AtomVecCharge(LAMMPS *lmp, int narg, char **arg) : 
   AtomVec(lmp, narg, arg)
 {
+  molecular = 0;
   mass_type = 1;
-  size_comm = 3;
+
+  comm_x_only = comm_f_only = 1;
+  size_forward = 3;
   size_reverse = 3;
   size_border = 7;
+  size_velocity = 3;
   size_data_atom = 6;
   size_data_vel = 4;
   xcol_data = 4;
@@ -132,6 +136,48 @@ int AtomVecCharge::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecCharge::pack_comm_vel(int n, int *list, double *buf,
+				 int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 void AtomVecCharge::unpack_comm(int n, int first, double *buf)
 {
   int i,m,last;
@@ -147,6 +193,24 @@ void AtomVecCharge::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecCharge::unpack_comm_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecCharge::pack_reverse(int n, int first, double *buf)
 {
   int i,m,last;
@@ -222,6 +286,56 @@ int AtomVecCharge::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecCharge::pack_border_vel(int n, int *list, double *buf,
+				   int pbc_flag, int *pbc)
+{
+  int i,j,m;
+  double dx,dy,dz;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = q[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = q[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecCharge::pack_border_one(int i, double *buf)
 {
   buf[0] = q[i];
@@ -250,6 +364,29 @@ void AtomVecCharge::unpack_border(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecCharge::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,m,last;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    q[i] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecCharge::unpack_border_one(int i, double *buf)
 {
   q[i] = buf[0];

src/atom_vec_charge.h

@@ -24,12 +24,16 @@ class AtomVecCharge : public AtomVec {
   void grow(int);
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   void unpack_comm(int, int, double *);
+  void unpack_comm_vel(int, int, double *);
   int pack_reverse(int, int, double *);
   void unpack_reverse(int, int *, double *);
   int pack_border(int, int *, double *, int, int *);
+  int pack_border_vel(int, int *, double *, int, int *);
   int pack_border_one(int, double *);
   void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int unpack_border_one(int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);

src/atom_vec_hybrid.cpp

@@ -57,7 +57,7 @@ AtomVecHybrid::AtomVecHybrid(LAMMPS *lmp, int narg, char **arg) :
   // hybrid settings are MAX or MIN of sub-style settings
   // hybrid sizes are minimial values plus extra values for each sub-style
 
-  size_comm = 3;
+  size_forward = 3;
   size_reverse = 3;
   size_border = 6;
   size_data_atom = 5;
@@ -73,16 +73,19 @@ AtomVecHybrid::AtomVecHybrid(LAMMPS *lmp, int narg, char **arg) :
     mass_type = MAX(mass_type,styles[k]->mass_type);
     shape_type = MAX(shape_type,styles[k]->shape_type);
     dipole_type = MAX(dipole_type,styles[k]->dipole_type);
+
     comm_x_only = MIN(comm_x_only,styles[k]->comm_x_only);
     comm_f_only = MIN(comm_f_only,styles[k]->comm_f_only);
-    ghost_velocity = MAX(ghost_velocity,styles[k]->ghost_velocity);
+    size_forward += styles[k]->size_forward - 3;
-
-    size_comm += styles[k]->size_comm - 3;
     size_reverse += styles[k]->size_reverse - 3;
     size_border += styles[k]->size_border - 6;
     size_data_atom += styles[k]->size_data_atom - 5;
     size_data_vel += styles[k]->size_data_vel - 4;
   }
+
+  size_velocity = 3;
+  if (atom->omega_flag) size_velocity += 3;
+  if (atom->angmom_flag) size_velocity += 3;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -122,6 +125,9 @@ void AtomVecHybrid::grow(int n)
   v = atom->v;
   f = atom->f;
 
+  omega = atom->omega;
+  angmom = atom->angmom;
+
   if (atom->nextra_grow)
     for (int iextra = 0; iextra < atom->nextra_grow; iextra++) 
       modify->fix[atom->extra_grow[iextra]]->grow_arrays(nmax);
@@ -192,6 +198,74 @@ int AtomVecHybrid::pack_comm(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecHybrid::pack_comm_vel(int n, int *list, double *buf,
+				 int pbc_flag, int *pbc)
+{
+  int i,j,k,m;
+  double dx,dy,dz;
+  int omega_flag = atom->omega_flag;
+  int angmom_flag = atom->angmom_flag;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      if (omega_flag) {
+	buf[m++] = omega[j][0];
+	buf[m++] = omega[j][1];
+	buf[m++] = omega[j][2];
+      }
+      if (angmom_flag) {
+	buf[m++] = angmom[j][0];
+	buf[m++] = angmom[j][1];
+	buf[m++] = angmom[j][2];
+      }
+      for (k = 0; k < nstyles; k++)
+	m += styles[k]->pack_comm_one(j,&buf[m]);
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0]*domain->xprd + pbc[5]*domain->xy + pbc[4]*domain->xz;
+      dy = pbc[1]*domain->yprd + pbc[3]*domain->yz;
+      dz = pbc[2]*domain->zprd;
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      if (omega_flag) {
+	buf[m++] = omega[j][0];
+	buf[m++] = omega[j][1];
+	buf[m++] = omega[j][2];
+      }
+      if (angmom_flag) {
+	buf[m++] = angmom[j][0];
+	buf[m++] = angmom[j][1];
+	buf[m++] = angmom[j][2];
+      }
+      for (k = 0; k < nstyles; k++)
+	m += styles[k]->pack_comm_one(j,&buf[m]);
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 void AtomVecHybrid::unpack_comm(int n, int first, double *buf)
 {
   int i,k,last;
@@ -209,6 +283,38 @@ void AtomVecHybrid::unpack_comm(int n, int first, double *buf)
 
 /* ---------------------------------------------------------------------- */
 
+void AtomVecHybrid::unpack_comm_vel(int n, int first, double *buf)
+{
+  int i,k,last;
+  int omega_flag = atom->omega_flag;
+  int angmom_flag = atom->angmom_flag;
+
+  int m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+    if (omega_flag) {
+      omega[i][0] = buf[m++];
+      omega[i][1] = buf[m++];
+      omega[i][2] = buf[m++];
+    }
+    if (angmom_flag) {
+      angmom[i][0] = buf[m++];
+      angmom[i][1] = buf[m++];
+      angmom[i][2] = buf[m++];
+    }
+    for (k = 0; k < nstyles; k++)
+      m += styles[k]->unpack_comm_one(i,&buf[m]);
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
 int AtomVecHybrid::pack_reverse(int n, int first, double *buf)
 {
   int i,k,last;
@@ -290,6 +396,80 @@ int AtomVecHybrid::pack_border(int n, int *list, double *buf,
 
 /* ---------------------------------------------------------------------- */
 
+int AtomVecHybrid::pack_border_vel(int n, int *list, double *buf,
+				   int pbc_flag, int *pbc)
+{
+  int i,j,k,m;
+  double dx,dy,dz;
+  int omega_flag = atom->omega_flag;
+  int angmom_flag = atom->angmom_flag;
+
+  m = 0;
+  if (pbc_flag == 0) {
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0];
+      buf[m++] = x[j][1];
+      buf[m++] = x[j][2];
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      if (omega_flag) {
+	buf[m++] = omega[j][0];
+	buf[m++] = omega[j][1];
+	buf[m++] = omega[j][2];
+      }
+      if (angmom_flag) {
+	buf[m++] = angmom[j][0];
+	buf[m++] = angmom[j][1];
+	buf[m++] = angmom[j][2];
+      }
+      for (k = 0; k < nstyles; k++)
+	m += styles[k]->pack_border_one(j,&buf[m]);
+    }
+  } else {
+    if (domain->triclinic == 0) {
+      dx = pbc[0]*domain->xprd;
+      dy = pbc[1]*domain->yprd;
+      dz = pbc[2]*domain->zprd;
+    } else {
+      dx = pbc[0];
+      dy = pbc[1];
+      dz = pbc[2];
+    }
+    for (i = 0; i < n; i++) {
+      j = list[i];
+      buf[m++] = x[j][0] + dx;
+      buf[m++] = x[j][1] + dy;
+      buf[m++] = x[j][2] + dz;
+      buf[m++] = tag[j];
+      buf[m++] = type[j];
+      buf[m++] = mask[j];
+      buf[m++] = v[j][0];
+      buf[m++] = v[j][1];
+      buf[m++] = v[j][2];
+      if (omega_flag) {
+	buf[m++] = omega[j][0];
+	buf[m++] = omega[j][1];
+	buf[m++] = omega[j][2];
+      }
+      if (angmom_flag) {
+	buf[m++] = angmom[j][0];
+	buf[m++] = angmom[j][1];
+	buf[m++] = angmom[j][2];
+      }
+      for (k = 0; k < nstyles; k++)
+	m += styles[k]->pack_border_one(j,&buf[m]);
+    }
+  }
+  return m;
+}
+
+/* ---------------------------------------------------------------------- */
+
 void AtomVecHybrid::unpack_border(int n, int first, double *buf)
 {
   int i,k,m,last;
@@ -309,6 +489,42 @@ void AtomVecHybrid::unpack_border(int n, int first, double *buf)
   }
 }
 
+/* ---------------------------------------------------------------------- */
+
+void AtomVecHybrid::unpack_border_vel(int n, int first, double *buf)
+{
+  int i,k,m,last;
+  int omega_flag = atom->omega_flag;
+  int angmom_flag = atom->angmom_flag;
+
+  m = 0;
+  last = first + n;
+  for (i = first; i < last; i++) {
+    if (i == nmax) grow(0);
+    x[i][0] = buf[m++];
+    x[i][1] = buf[m++];
+    x[i][2] = buf[m++];
+    tag[i] = static_cast<int> (buf[m++]);
+    type[i] = static_cast<int> (buf[m++]);
+    mask[i] = static_cast<int> (buf[m++]);
+    v[i][0] = buf[m++];
+    v[i][1] = buf[m++];
+    v[i][2] = buf[m++];
+    if (omega_flag) {
+      omega[i][0] = buf[m++];
+      omega[i][1] = buf[m++];
+      omega[i][2] = buf[m++];
+    }
+    if (angmom_flag) {
+      angmom[i][0] = buf[m++];
+      angmom[i][1] = buf[m++];
+      angmom[i][2] = buf[m++];
+    }
+    for (k = 0; k < nstyles; k++)
+      m += styles[k]->unpack_border_one(i,&buf[m]);
+  }
+}
+
 /* ----------------------------------------------------------------------
    pack data for atom I for sending to another proc
    pack each sub-style one after the other

src/atom_vec_hybrid.h

@@ -30,11 +30,15 @@ class AtomVecHybrid : public AtomVec {
   void reset_special();
   void copy(int, int);
   int pack_comm(int, int *, double *, int, int *);
+  int pack_comm_vel(int, int *, double *, int, int *);
   void unpack_comm(int, int, double *);
+  void unpack_comm_vel(int, int, double *);
   int pack_reverse(int, int, double *);
   void unpack_reverse(int, int *, double *);
   int pack_border(int, int *, double *, int, int *);
+  int pack_border_vel(int, int *, double *, int, int *);
   void unpack_border(int, int, double *);
+  void unpack_border_vel(int, int, double *);
   int pack_exchange(int, double *);
   int unpack_exchange(double *);
   int size_restart();
@@ -49,6 +53,7 @@ class AtomVecHybrid : public AtomVec {
  private:
   int *tag,*type,*mask,*image;
   double **x,**v,**f;
+  double **omega,**angmom;
 };
 
 }

src/comm.cpp

@@ -63,6 +63,7 @@ Comm::Comm(LAMMPS *lmp) : Pointers(lmp)
   multilo = multihi = NULL;
   cutghostmulti = NULL;
   cutghostuser = 0.0;
+  ghost_velocity = 0;
 
   // initialize comm buffers & exchange memory
 
@@ -175,16 +176,28 @@ void Comm::init()
   map_style = atom->map_style;
 
   // comm_only = 1 if only x,f are exchanged in forward/reverse comm
+  // comm_x_only = 0 if ghost_velocity since velocities are added
 
   comm_x_only = atom->avec->comm_x_only;
   comm_f_only = atom->avec->comm_f_only;
+  if (ghost_velocity) comm_x_only = 0;
+
+  // set per-atom sizes for forward/reverse/border comm
+  // augment by velocity quantities if needed
+
+  size_forward = atom->avec->size_forward;
+  size_reverse = atom->avec->size_reverse;
+  size_border = atom->avec->size_border;
+
+  if (ghost_velocity) size_forward += atom->avec->size_velocity;
+  if (ghost_velocity) size_border += atom->avec->size_velocity;
 
   // maxforward = # of datums in largest forward communication
   // maxreverse = # of datums in largest reverse communication
   // query pair,fix,compute for their requirements
 
-  maxforward = MAX(atom->avec->size_comm,atom->avec->size_border);
+  maxforward = MAX(size_forward,size_border);
-  maxreverse = atom->avec->size_reverse;
+  maxreverse = size_reverse;
 
   if (force->pair) maxforward = MAX(maxforward,force->pair->comm_forward);
   if (force->pair) maxreverse = MAX(maxreverse,force->pair->comm_reverse);
@@ -391,7 +404,7 @@ void Comm::setup()
 
 /* ----------------------------------------------------------------------
    communication of atom coords every timestep
-   other stuff may also be sent via pack/unpack routines
+   other per-atom attributes may also be sent via pack/unpack routines
 ------------------------------------------------------------------------- */
 
 void Comm::communicate()
@@ -410,16 +423,24 @@ void Comm::communicate()
   for (int iswap = 0; iswap < nswap; iswap++) {
     if (sendproc[iswap] != me) {
       if (comm_x_only) {
-	if (size_comm_recv[iswap]) buf = x[firstrecv[iswap]];
+	if (size_forward_recv[iswap]) buf = x[firstrecv[iswap]];
 	else buf = NULL;
-	MPI_Irecv(buf,size_comm_recv[iswap],MPI_DOUBLE,
+	MPI_Irecv(buf,size_forward_recv[iswap],MPI_DOUBLE,
 		  recvproc[iswap],0,world,&request);
 	n = avec->pack_comm(sendnum[iswap],sendlist[iswap],
 			    buf_send,pbc_flag[iswap],pbc[iswap]);
 	MPI_Send(buf_send,n,MPI_DOUBLE,sendproc[iswap],0,world);
 	MPI_Wait(&request,&status);
+      } else if (ghost_velocity) {
+	MPI_Irecv(buf_recv,size_forward_recv[iswap],MPI_DOUBLE,
+		  recvproc[iswap],0,world,&request);
+	n = avec->pack_comm_vel(sendnum[iswap],sendlist[iswap],
+				buf_send,pbc_flag[iswap],pbc[iswap]);
+	MPI_Send(buf_send,n,MPI_DOUBLE,sendproc[iswap],0,world);
+	MPI_Wait(&request,&status);
+	avec->unpack_comm_vel(recvnum[iswap],firstrecv[iswap],buf_recv);
       } else {
-	MPI_Irecv(buf_recv,size_comm_recv[iswap],MPI_DOUBLE,
+	MPI_Irecv(buf_recv,size_forward_recv[iswap],MPI_DOUBLE,
 		  recvproc[iswap],0,world,&request);
 	n = avec->pack_comm(sendnum[iswap],sendlist[iswap],
 			    buf_send,pbc_flag[iswap],pbc[iswap]);
@@ -430,13 +451,19 @@ void Comm::communicate()
 
     } else {
       if (comm_x_only) {
-	if (sendnum[iswap])
+	if (sendnum[iswap]) {
 	  n = avec->pack_comm(sendnum[iswap],sendlist[iswap],
-			      x[firstrecv[iswap]],pbc_flag[iswap],pbc[iswap]);
+			      x[firstrecv[iswap]],pbc_flag[iswap],
-      } else {
+			      pbc[iswap]);
-	n = avec->pack_comm(sendnum[iswap],sendlist[iswap],
+	} else if (ghost_velocity) {
-			    buf_send,pbc_flag[iswap],pbc[iswap]);
+	  n = avec->pack_comm_vel(sendnum[iswap],sendlist[iswap],
-	avec->unpack_comm(recvnum[iswap],firstrecv[iswap],buf_send);
+				  buf_send,pbc_flag[iswap],pbc[iswap]);
+	  avec->unpack_comm_vel(recvnum[iswap],firstrecv[iswap],buf_send);
+	} else {
+	  n = avec->pack_comm(sendnum[iswap],sendlist[iswap],
+			      buf_send,pbc_flag[iswap],pbc[iswap]);
+	  avec->unpack_comm(recvnum[iswap],firstrecv[iswap],buf_send);
+	}
       }
     }
   }
@@ -444,7 +471,7 @@ void Comm::communicate()
 
 /* ----------------------------------------------------------------------
    reverse communication of forces on atoms every timestep 
-   other stuff can also be sent via pack/unpack routines
+   other per-atom attributes may also be sent via pack/unpack routines
 ------------------------------------------------------------------------- */
       
 void Comm::reverse_communicate()
@@ -623,7 +650,6 @@ void Comm::borders()
   MPI_Request request;
   MPI_Status status;
   AtomVec *avec = atom->avec;
-  int size_border = avec->size_border;
 
   // clear old ghosts
 
@@ -719,8 +745,12 @@ void Comm::borders()
 
       if (nsend*size_border > maxsend)
 	grow_send(nsend*size_border,0);
-      n = avec->pack_border(nsend,sendlist[iswap],buf_send,
+      if (ghost_velocity)
-			    pbc_flag[iswap],pbc[iswap]);
+	n = avec->pack_border_vel(nsend,sendlist[iswap],buf_send,
+				  pbc_flag[iswap],pbc[iswap]);
+      else
+	n = avec->pack_border(nsend,sendlist[iswap],buf_send,
+			      pbc_flag[iswap],pbc[iswap]);
       
       // swap atoms with other proc
       // put incoming ghosts at end of my atom arrays
@@ -743,7 +773,10 @@ void Comm::borders()
 
       // unpack buffer
 
-      avec->unpack_border(nrecv,atom->nlocal+atom->nghost,buf);
+      if (ghost_velocity)
+	avec->unpack_border_vel(nrecv,atom->nlocal+atom->nghost,buf);
+      else
+	avec->unpack_border(nrecv,atom->nlocal+atom->nghost,buf);
 
       // set all pointers & counters
 
@@ -751,9 +784,9 @@ void Comm::borders()
       rmax = MAX(rmax,nrecv);
       sendnum[iswap] = nsend;
       recvnum[iswap] = nrecv;
-      size_comm_recv[iswap] = nrecv * avec->size_comm;
+      size_forward_recv[iswap] = nrecv*size_forward;
-      size_reverse_send[iswap] = nrecv * avec->size_reverse;
+      size_reverse_send[iswap] = nrecv*size_reverse;
-      size_reverse_recv[iswap] = nsend * avec->size_reverse;
+      size_reverse_recv[iswap] = nsend*size_reverse;
       firstrecv[iswap] = atom->nlocal + atom->nghost;
       atom->nghost += nrecv;
       iswap++;
@@ -1469,7 +1502,7 @@ void Comm::allocate_swap(int n)
   recvnum = (int *) memory->smalloc(n*sizeof(int),"comm:recvnum");
   sendproc = (int *) memory->smalloc(n*sizeof(int),"comm:sendproc");
   recvproc = (int *) memory->smalloc(n*sizeof(int),"comm:recvproc");
-  size_comm_recv = (int *) memory->smalloc(n*sizeof(int),"comm:size");
+  size_forward_recv = (int *) memory->smalloc(n*sizeof(int),"comm:size");
   size_reverse_send = (int *) memory->smalloc(n*sizeof(int),"comm:size");
   size_reverse_recv = (int *) memory->smalloc(n*sizeof(int),"comm:size");
   slablo = (double *) memory->smalloc(n*sizeof(double),"comm:slablo");
@@ -1500,7 +1533,7 @@ void Comm::free_swap()
   memory->sfree(recvnum);
   memory->sfree(sendproc);
   memory->sfree(recvproc);
-  memory->sfree(size_comm_recv);
+  memory->sfree(size_forward_recv);
   memory->sfree(size_reverse_send);
   memory->sfree(size_reverse_recv);
   memory->sfree(slablo);
@@ -1549,6 +1582,12 @@ void Comm::set(int narg, char **arg)
       if (cutghostuser < 0.0) 
 	error->all("Invalid cutoff in communicate command");
       iarg += 2;
+    } else if (strcmp(arg[iarg],"vel") == 0) {
+      if (iarg+2 > narg) error->all("Illegal communicate command");
+      if (strcmp(arg[iarg+1],"yes") == 0) ghost_velocity = 1;
+      else if (strcmp(arg[iarg+1],"yes") == 0) ghost_velocity = 0;
+      else error->all("Illegal communicate command");
+      iarg += 2;
     } else error->all("Illegal communicate command");
   }
 }

src/comm.h

@@ -28,6 +28,7 @@ class Comm : protected Pointers {
   int procneigh[3][2];              // my 6 neighboring procs
   int nswap;                        // # of swaps to perform
   int need[3];                      // procs I need atoms from in each dim
+  int ghost_velocity;               // 1 if ghost atoms have velocity, 0 if not
   int maxforward_fix;               // comm sizes called from Fix,Pair
   int maxreverse_fix;
   int maxforward_pair;
@@ -35,6 +36,7 @@ class Comm : protected Pointers {
   double cutghost[3];               // cutoffs used for acquiring ghost atoms
   double cutghostuser;              // user-specified ghost cutoff
   
+
   Comm(class LAMMPS *);
   ~Comm();
 
@@ -61,9 +63,12 @@ class Comm : protected Pointers {
  private:
   int triclinic;                    // 0 if domain is orthog, 1 if triclinic
   int maxswap;                      // max # of swaps memory is allocated for
+  int size_forward;                 // # of per-atom datums in forward comm
+  int size_reverse;                 // # of datums in reverse comm
+  int size_border;                  // # of datums in forward border comm
   int *sendnum,*recvnum;            // # of atoms to send/recv in each swap
   int *sendproc,*recvproc;          // proc to send/recv to/from at each swap
-  int *size_comm_recv;              // # of values to recv in each forward comm
+  int *size_forward_recv;           // # of values to recv in each forward comm
   int *size_reverse_send;           // # to send in each reverse comm
   int *size_reverse_recv;           // # to recv in each reverse comm
   double *slablo,*slabhi;           // bounds of slab to send at each swap

src/compute_heat_flux.cpp

@@ -24,6 +24,7 @@
 #include "atom_vec.h"
 #include "update.h"
 #include "force.h"
+#include "comm.h"
 #include "pair.h"
 #include "modify.h"
 #include "group.h"
@@ -75,7 +76,7 @@ void ComputeHeatFlux::init()
 {
   // error checks
 
-  if (atom->avec->ghost_velocity == 0)
+  if (comm->ghost_velocity == 0)
     error->all("Compute heat/flux requires ghost atoms store velocity");
 
   if (force->pair == NULL || force->pair->single_enable == 0)

src/pair_dpd.cpp

@@ -251,7 +251,7 @@ void PairDPD::coeff(int narg, char **arg)
 
 void PairDPD::init_style()
 {
-  if (atom->avec->ghost_velocity == 0)
+  if (comm->ghost_velocity == 0)
     error->all("Pair dpd requires ghost atoms store velocity");
 
   // if newton off, forces between atoms ij will be double computed

src/style.h

@@ -312,6 +312,7 @@ MinimizeStyle(sd,MinSD)
 #include "pair_buck_coul_cut.h"
 #include "pair_coul_cut.h"
 #include "pair_coul_debye.h"
+#include "pair_dpd.h"
 #include "pair_hybrid.h"
 #include "pair_hybrid_overlay.h"
 #include "pair_lj_cut.h"
@@ -333,6 +334,7 @@ PairStyle(buck,PairBuck)
 PairStyle(buck/coul/cut,PairBuckCoulCut)
 PairStyle(coul/cut,PairCoulCut)
 PairStyle(coul/debye,PairCoulDebye)
+PairStyle(dpd,PairDPD)
 PairStyle(hybrid,PairHybrid)
 PairStyle(hybrid/overlay,PairHybridOverlay)
 PairStyle(lj/cut,PairLJCut)

src/style_packages.h

@@ -17,7 +17,6 @@
 #include "style_class2.h"
 #include "style_colloid.h"
 #include "style_dipole.h"
-#include "style_dpd.h"
 #include "style_dsmc.h"
 #include "style_gpu.h"
 #include "style_granular.h"