From 15068b63547898c62806153c59605c914196207d Mon Sep 17 00:00:00 2001
From: sjplimp <sjplimp@f3b2605a-c512-4ea7-a41b-209d697bcdaa>
Date: Thu, 2 Jun 2016 14:08:20 +0000
Subject: [PATCH] git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@15100
 f3b2605a-c512-4ea7-a41b-209d697bcdaa

---
 src/MANYBODY/pair_airebo.cpp                  |   7 +-
 src/MOLECULE/improper_umbrella.cpp            |   5 +-
 src/Makefile                                  |   2 +-
 src/Purge.list                                | 134 ++++++++++++++++++
 src/USER-MANIFOLD/fix_nve_manifold_rattle.cpp |  18 +--
 src/USER-MANIFOLD/fix_nvt_manifold_rattle.cpp |  22 +--
 src/USER-OMP/pair_airebo_omp.cpp              |   7 +-
 src/USER-REAXC/pair_reax_c.cpp                |   2 +-
 src/VORONOI/compute_voronoi_atom.cpp          |  21 ++-
 src/VORONOI/compute_voronoi_atom.h            |   2 +-
 src/compute_com.cpp                           |   1 +
 src/compute_gyration.cpp                      |   1 +
 src/compute_msd.cpp                           |   1 +
 src/info.cpp                                  |   6 +-
 14 files changed, 193 insertions(+), 36 deletions(-)

diff --git a/src/MANYBODY/pair_airebo.cpp b/src/MANYBODY/pair_airebo.cpp
index d4e57c2c6c..2bb1c1da22 100644
--- a/src/MANYBODY/pair_airebo.cpp
+++ b/src/MANYBODY/pair_airebo.cpp
@@ -592,6 +592,8 @@ void PairAIREBO::FLJ(int eflag, int vflag)
       // if outside of 2-path cutoff but inside 4-path cutoff,
       //   best = 0.0, test 3-,4-paths
       // if inside 2-path cutoff, best = wij, only test 3-,4-paths if best < 1
+      npath = testpath = done = 0;
+      best = 0.0;
 
       if (rijsq >= cutljsq[itype][jtype]) continue;
       rij = sqrt(rijsq);
@@ -608,7 +610,6 @@ void PairAIREBO::FLJ(int eflag, int vflag)
         else testpath = 0;
       }
 
-      done = 0;
       if (testpath) {
 
         // test all 3-body paths = I-K-J
@@ -629,7 +630,7 @@ void PairAIREBO::FLJ(int eflag, int vflag)
           if (rsq < rcmaxsq[itype][ktype]) {
             rik = sqrt(rsq);
             wik = Sp(rik,rcmin[itype][ktype],rcmax[itype][ktype],dwik);
-          } else wik = 0.0;
+          } else { dwik = wik = 0.0; rikS = rik = 1.0; }
 
           if (wik > best) {
             deljk[0] = x[j][0] - x[k][0];
@@ -679,7 +680,7 @@ void PairAIREBO::FLJ(int eflag, int vflag)
               if (rsq < rcmaxsq[ktype][mtype]) {
                 rkm = sqrt(rsq);
                 wkm = Sp(rkm,rcmin[ktype][mtype],rcmax[ktype][mtype],dwkm);
-              } else wkm = 0.0;
+              } else { dwkm = wkm = 0.0; rkmS = rkm = 1.0; }
 
               if (wik*wkm > best) {
                 deljm[0] = x[j][0] - x[m][0];
diff --git a/src/MOLECULE/improper_umbrella.cpp b/src/MOLECULE/improper_umbrella.cpp
index 14551e6c5c..30a62ffb18 100644
--- a/src/MOLECULE/improper_umbrella.cpp
+++ b/src/MOLECULE/improper_umbrella.cpp
@@ -37,7 +37,10 @@ using namespace MathConst;
 
 /* ---------------------------------------------------------------------- */
 
-ImproperUmbrella::ImproperUmbrella(LAMMPS *lmp) : Improper(lmp) {}
+ImproperUmbrella::ImproperUmbrella(LAMMPS *lmp) : Improper(lmp)
+{
+  writedata = 1;
+}
 
 /* ---------------------------------------------------------------------- */
 
diff --git a/src/Makefile b/src/Makefile
index 85b212d539..7cf9d1730f 100755
--- a/src/Makefile
+++ b/src/Makefile
@@ -54,7 +54,7 @@ PACKUSER = user-atc user-awpmd user-cg-cmm user-colvars \
 	   user-vtk
 
 PACKLIB = compress gpu kim kokkos meam poems python reax voronoi \
-	  user-atc user-awpmd user-colvars user-cuda user-h5md user-molfile \
+	  user-atc user-awpmd user-colvars user-h5md user-molfile \
 	  user-qmmm user-quip user-vtk
 
 PACKALL = $(PACKAGE) $(PACKUSER)
diff --git a/src/Purge.list b/src/Purge.list
index 79ad3c0c42..37d137e838 100644
--- a/src/Purge.list
+++ b/src/Purge.list
@@ -13,6 +13,140 @@ style_kspace.h
 style_minimize.h
 style_pair.h
 style_region.h
+# deleted on 31 May 2016
+fix_ave_spatial_sphere.cpp
+fix_ave_spatial_sphere.h
+atom_vec_angle_cuda.cpp
+atom_vec_angle_cuda.h
+atom_vec_atomic_cuda.cpp
+atom_vec_atomic_cuda.h
+atom_vec_charge_cuda.cpp
+atom_vec_charge_cuda.h
+atom_vec_full_cuda.cpp
+atom_vec_full_cuda.h
+comm_cuda.cpp
+comm_cuda.h
+compute_pe_cuda.cpp
+compute_pe_cuda.h
+compute_pressure_cuda.cpp
+compute_pressure_cuda.h
+compute_temp_cuda.cpp
+compute_temp_cuda.h
+compute_temp_partial_cuda.cpp
+compute_temp_partial_cuda.h
+cuda.cpp
+cuda_data.h
+cuda_modify_flags.h
+cuda_neigh_list.cpp
+cuda_neigh_list.h
+domain_cuda.cpp
+domain_cuda.h
+fft3d_cuda.cpp
+fft3d_cuda.h
+fft3d_wrap_cuda.cpp
+fft3d_wrap_cuda.h
+fix_addforce_cuda.cpp
+fix_addforce_cuda.h
+fix_aveforce_cuda.cpp
+fix_aveforce_cuda.h
+fix_enforce2d_cuda.cpp
+fix_enforce2d_cuda.h
+fix_freeze_cuda.cpp
+fix_freeze_cuda.h
+fix_gravity_cuda.cpp
+fix_gravity_cuda.h
+fix_nh_cuda.cpp
+fix_nh_cuda.h
+fix_npt_cuda.cpp
+fix_npt_cuda.h
+fix_nve_cuda.cpp
+fix_nve_cuda.h
+fix_nvt_cuda.cpp
+fix_nvt_cuda.h
+fix_set_force_cuda.cpp
+fix_set_force_cuda.h
+fix_shake_cuda.cpp
+fix_shake_cuda.h
+fix_temp_berendsen_cuda.cpp
+fix_temp_berendsen_cuda.h
+fix_temp_rescale_cuda.cpp
+fix_temp_rescale_cuda.h
+fix_temp_rescale_limit_cuda.cpp
+fix_temp_rescale_limit_cuda.h
+fix_viscous_cuda.cpp
+fix_viscous_cuda.h
+modify_cuda.cpp
+modify_cuda.h
+neighbor_cuda.cpp
+neighbor_cuda.h
+neigh_full_cuda.cpp
+pair_born_coul_long_cuda.cpp
+pair_born_coul_long_cuda.h
+pair_buck_coul_cut_cuda.cpp
+pair_buck_coul_cut_cuda.h
+pair_buck_coul_long_cuda.cpp
+pair_buck_coul_long_cuda.h
+pair_buck_cuda.cpp
+pair_buck_cuda.h
+pair_eam_alloy_cuda.cpp
+pair_eam_alloy_cuda.h
+pair_eam_cuda.cpp
+pair_eam_cuda.h
+pair_eam_fs_cuda.cpp
+pair_eam_fs_cuda.h
+pair_gran_hooke_cuda.cpp
+pair_gran_hooke_cuda.h
+pair_lj96_cut_cuda.cpp
+pair_lj96_cut_cuda.h
+pair_lj_charmm_coul_charmm_cuda.cpp
+pair_lj_charmm_coul_charmm_cuda.h
+pair_lj_charmm_coul_charmm_implicit_cuda.cpp
+pair_lj_charmm_coul_charmm_implicit_cuda.h
+pair_lj_charmm_coul_long_cuda.cpp
+pair_lj_charmm_coul_long_cuda.h
+pair_lj_class2_coul_cut_cuda.cpp
+pair_lj_class2_coul_cut_cuda.h
+pair_lj_class2_coul_long_cuda.cpp
+pair_lj_class2_coul_long_cuda.h
+pair_lj_class2_cuda.cpp
+pair_lj_class2_cuda.h
+pair_lj_cut_coul_cut_cuda.cpp
+pair_lj_cut_coul_cut_cuda.h
+pair_lj_cut_coul_debye_cuda.cpp
+pair_lj_cut_coul_debye_cuda.h
+pair_lj_cut_coul_long_cuda.cpp
+pair_lj_cut_coul_long_cuda.h
+pair_lj_cut_cuda.cpp
+pair_lj_cut_cuda.h
+pair_lj_cut_experimental_cuda.cpp
+pair_lj_cut_experimental_cuda.h
+pair_lj_expand_cuda.cpp
+pair_lj_expand_cuda.h
+pair_lj_gromacs_coul_gromacs_cuda.cpp
+pair_lj_gromacs_coul_gromacs_cuda.h
+pair_lj_gromacs_cuda.cpp
+pair_lj_gromacs_cuda.h
+pair_lj_sdk_coul_long_cuda.cpp
+pair_lj_sdk_coul_long_cuda.h
+pair_lj_sdk_cuda.cpp
+pair_lj_sdk_cuda.h
+pair_lj_smooth_cuda.cpp
+pair_lj_smooth_cuda.h
+pair_morse_cuda.cpp
+pair_morse_cuda.h
+pair_sw_cuda.cpp
+pair_sw_cuda.h
+pair_tersoff_cuda.cpp
+pair_tersoff_cuda.h
+pair_tersoff_zbl_cuda.cpp
+pair_tersoff_zbl_cuda.h
+pppm_cuda.cpp
+pppm_cuda.h
+pppm_old.cpp
+pppm_old.h
+user_cuda.h
+verlet_cuda.cpp
+verlet_cuda.h
 # deleted on 11 May 2016
 pair_dpd_conservative.cpp
 pair_dpd_conservative.h
diff --git a/src/USER-MANIFOLD/fix_nve_manifold_rattle.cpp b/src/USER-MANIFOLD/fix_nve_manifold_rattle.cpp
index d2d73de6bd..f485c4b7b4 100644
--- a/src/USER-MANIFOLD/fix_nve_manifold_rattle.cpp
+++ b/src/USER-MANIFOLD/fix_nve_manifold_rattle.cpp
@@ -64,14 +64,16 @@ enum { CONST, EQUAL }; // For treating the variables.
 
 static const char* cite_fix_nve_manifold_rattle =
   "fix nve/manifold/rattle command:\n\n"
-  "@article{paquay-2016,\n"
-  " author = {Paquay, Stefan and Kusters, Remy}, \n"
-  " eprint = {arXiv:1411.3019}, \n"
-  " title = {A method for molecular dynamics on curved surfaces}, \n"
-  " url = {http://arxiv.org/abs/1411.3019}, \n"
-  " year = 2016, \n"
-  " journal = {Arxiv preprint, \\url{http://arxiv.org/abs/1411.3019}}, \n"
-  " note = {To be published in Biophys. J.}, \n"
+  "   author        = {Paquay, Stefan and Kusters, Remy},\n"
+  "   doi           = {10.1016/j.bpj.2016.02.017},\n"
+  "   issn          = {0006-3495},\n"
+  "   journal       = {Biophysical Journal},\n"
+  "   month         = apr,\n"
+  "   number        = {6},\n"
+  "   pages         = {1226--1233},\n"
+  "   title         = {{A Method for Molecular Dynamics on Curved Surfaces}},\n"
+  "   volume        = {110},\n"
+  "   year          = {2016}\n"
   "}\n\n";
 
 
diff --git a/src/USER-MANIFOLD/fix_nvt_manifold_rattle.cpp b/src/USER-MANIFOLD/fix_nvt_manifold_rattle.cpp
index 3801508336..c304a8718d 100644
--- a/src/USER-MANIFOLD/fix_nvt_manifold_rattle.cpp
+++ b/src/USER-MANIFOLD/fix_nvt_manifold_rattle.cpp
@@ -65,18 +65,18 @@ enum {NOBIAS,BIAS};
 
 
 static const char* cite_fix_nvt_manifold_rattle =
-  "fix nve/manifold/rattle command:\n\n"
-  "@article{paquay-2016,\n"
-  " author = {Paquay, Stefan and Kusters, Remy}, \n"
-  " eprint = {arXiv:1411.3019}, \n"
-  " title = {A method for molecular dynamics on curved surfaces}, \n"
-  " url = {http://arxiv.org/abs/1411.3019}, \n"
-  " year = 2016, \n"
-  " journal = {Arxiv preprint, \\url{http://arxiv.org/abs/1411.3019}}, \n"
-  " note = {To be published in Biophys. J.}, \n"
+  "fix nvt/manifold/rattle command:\n\n"
+  "   author        = {Paquay, Stefan and Kusters, Remy},\n"
+  "   doi           = {10.1016/j.bpj.2016.02.017},\n"
+  "   issn          = {0006-3495},\n"
+  "   journal       = {Biophysical Journal},\n"
+  "   month         = apr,\n"
+  "   number        = {6},\n"
+  "   pages         = {1226--1233},\n"
+  "   title         = {{A Method for Molecular Dynamics on Curved Surfaces}},\n"
+  "   volume        = {110},\n"
+  "   year          = {2016}\n"
   "}\n\n";
-
-
 /* ---------------------------------------------------------------------- */
 
 FixNVTManifoldRattle::FixNVTManifoldRattle(LAMMPS *lmp, int narg, char **arg,
diff --git a/src/USER-OMP/pair_airebo_omp.cpp b/src/USER-OMP/pair_airebo_omp.cpp
index 35a644007a..27089fd6b3 100644
--- a/src/USER-OMP/pair_airebo_omp.cpp
+++ b/src/USER-OMP/pair_airebo_omp.cpp
@@ -2022,6 +2022,8 @@ void PairAIREBOOMP::FLJ_thr(int ifrom, int ito, int evflag, int eflag,
       // if outside of 2-path cutoff but inside 4-path cutoff,
       //   best = 0.0, test 3-,4-paths
       // if inside 2-path cutoff, best = wij, only test 3-,4-paths if best < 1
+      npath = testpath = done = 0;
+      best = 0.0;
 
       if (rijsq >= cutljsq[itype][jtype]) continue;
       rij = sqrt(rijsq);
@@ -2038,7 +2040,6 @@ void PairAIREBOOMP::FLJ_thr(int ifrom, int ito, int evflag, int eflag,
         else testpath = 0;
       }
 
-      done = 0;
       if (testpath) {
 
         // test all 3-body paths = I-K-J
@@ -2059,7 +2060,7 @@ void PairAIREBOOMP::FLJ_thr(int ifrom, int ito, int evflag, int eflag,
           if (rsq < rcmaxsq[itype][ktype]) {
             rik = sqrt(rsq);
             wik = Sp(rik,rcmin[itype][ktype],rcmax[itype][ktype],dwik);
-          } else wik = 0.0;
+          } else { dwik = wik = 0.0; rikS = rik = 1.0; }
 
           if (wik > best) {
             deljk[0] = x[j][0] - x[k][0];
@@ -2109,7 +2110,7 @@ void PairAIREBOOMP::FLJ_thr(int ifrom, int ito, int evflag, int eflag,
               if (rsq < rcmaxsq[ktype][mtype]) {
                 rkm = sqrt(rsq);
                 wkm = Sp(rkm,rcmin[ktype][mtype],rcmax[ktype][mtype],dwkm);
-              } else wkm = 0.0;
+              } else { dwkm = wkm = 0.0; rkmS = rkm = 1.0; }
 
               if (wik*wkm > best) {
                 deljm[0] = x[j][0] - x[m][0];
diff --git a/src/USER-REAXC/pair_reax_c.cpp b/src/USER-REAXC/pair_reax_c.cpp
index 352a9468b3..59b39d4534 100644
--- a/src/USER-REAXC/pair_reax_c.cpp
+++ b/src/USER-REAXC/pair_reax_c.cpp
@@ -247,7 +247,7 @@ void PairReaxC::settings(int narg, char **arg)
       system->safezone = force->numeric(FLERR,arg[iarg+1]);
       if (system->safezone < 0.0)
 	error->all(FLERR,"Illegal pair_style reax/c safezone command");
-      system->saferzone = system->safezone + 0.2;
+      system->saferzone = system->safezone*1.2;
       iarg += 2;
     } else if (strcmp(arg[iarg],"mincap") == 0) {
       if (iarg+2 > narg) error->all(FLERR,"Illegal pair_style reax/c command");
diff --git a/src/VORONOI/compute_voronoi_atom.cpp b/src/VORONOI/compute_voronoi_atom.cpp
index f4f48fd5e8..5239a27ecc 100644
--- a/src/VORONOI/compute_voronoi_atom.cpp
+++ b/src/VORONOI/compute_voronoi_atom.cpp
@@ -61,6 +61,7 @@ ComputeVoronoi::ComputeVoronoi(LAMMPS *lmp, int narg, char **arg) :
   con_mono = NULL;
   con_poly = NULL;
   tags = NULL;
+  oldmaxtag = 0;
   occvec = sendocc = lroot = lnext = NULL;
   faces = NULL;
 
@@ -177,6 +178,8 @@ ComputeVoronoi::~ComputeVoronoi()
 
 void ComputeVoronoi::init()
 {
+  if (occupation && (atom->tag_enable == 0))
+    error->all(FLERR,"Compute voronoi/atom occupation requires atom IDs");
 }
 
 /* ----------------------------------------------------------------------
@@ -213,11 +216,14 @@ void ComputeVoronoi::compute_peratom()
       lnext = NULL;
       lmax = 0;
 
-      // build the occupation buffer
+      // build the occupation buffer.
+      // NOTE: we cannot make it of length oldnatoms, as tags may not be
+      // consecutive at this point, e.g. due to deleted or lost atoms.
       oldnatoms = atom->natoms;
-      memory->create(occvec,oldnatoms,"voronoi/atom:occvec");
+      oldmaxtag = atom->map_tag_max;
+      memory->create(occvec,oldmaxtag,"voronoi/atom:occvec");
 #ifdef NOTINPLACE
-      memory->create(sendocc,oldnatoms,"voronoi/atom:sendocc");
+      memory->create(sendocc,oldmaxtag,"voronoi/atom:sendocc");
 #endif
     }
 
@@ -440,7 +446,14 @@ void ComputeVoronoi::checkOccupation()
   // cherry pick currently owned atoms
   for (i=0; i<nlocal; i++) {
     // set the new atom count in the atom's first frame voronoi cell
-    voro[i][0] = occvec[atom->tag[i]-1];
+    // but take into account that new atoms might have been added to
+    // the system, so we can only look up occupancy for tags that are
+    // smaller or equal to the recorded largest tag.
+    tagint mytag = atom->tag[i];
+    if (mytag > oldmaxtag)
+      voro[i][0] = 0;
+    else
+      voro[i][0] = occvec[mytag-1];
   }
 }
 
diff --git a/src/VORONOI/compute_voronoi_atom.h b/src/VORONOI/compute_voronoi_atom.h
index 7d23a611bd..8fd187cea3 100644
--- a/src/VORONOI/compute_voronoi_atom.h
+++ b/src/VORONOI/compute_voronoi_atom.h
@@ -55,7 +55,7 @@ class ComputeVoronoi : public Compute {
   enum { VOROSURF_NONE, VOROSURF_ALL, VOROSURF_GROUP } surface;
   bool onlyGroup, occupation;
 
-  tagint *tags;
+  tagint *tags, oldmaxtag;
   int *occvec, *sendocc, *lroot, *lnext, lmax, oldnatoms, oldnall;
   int faces_flag, nfaces, nfacesmax;
   double **faces;
diff --git a/src/compute_com.cpp b/src/compute_com.cpp
index e086210322..27606cbc6f 100644
--- a/src/compute_com.cpp
+++ b/src/compute_com.cpp
@@ -51,6 +51,7 @@ void ComputeCOM::init()
 void ComputeCOM::compute_vector()
 {
   invoked_vector = update->ntimestep;
+  if (group->dynamic[igroup]) masstotal = group->mass(igroup);
 
   group->xcm(igroup,masstotal,vector);
 }
diff --git a/src/compute_gyration.cpp b/src/compute_gyration.cpp
index fab413e670..0c1c58a4fd 100644
--- a/src/compute_gyration.cpp
+++ b/src/compute_gyration.cpp
@@ -57,6 +57,7 @@ double ComputeGyration::compute_scalar()
   invoked_scalar = update->ntimestep;
 
   double xcm[3];
+  if (group->dynamic[igroup]) masstotal = group->mass(igroup);
   group->xcm(igroup,masstotal,xcm);
   scalar = group->gyration(igroup,masstotal,xcm);
   return scalar;
diff --git a/src/compute_msd.cpp b/src/compute_msd.cpp
index 32d93f9334..7f005d2759 100644
--- a/src/compute_msd.cpp
+++ b/src/compute_msd.cpp
@@ -34,6 +34,7 @@ ComputeMSD::ComputeMSD(LAMMPS *lmp, int narg, char **arg) :
   size_vector = 4;
   extvector = 0;
   create_attribute = 1;
+  dynamic_group_allow = 0;
 
   // optional args
 
diff --git a/src/info.cpp b/src/info.cpp
index 016a4c6526..459894ce48 100644
--- a/src/info.cpp
+++ b/src/info.cpp
@@ -327,9 +327,9 @@ void Info::command(int narg, char **arg)
               bstyles[domain->boundary[0][0]],bstyles[domain->boundary[0][1]],
               bstyles[domain->boundary[1][0]],bstyles[domain->boundary[1][1]],
               bstyles[domain->boundary[2][0]],bstyles[domain->boundary[2][1]]);
-      fprintf(out,"Xlo, zhi = %g, %g\n", domain->boxlo[0], domain->boxhi[0]);
-      fprintf(out,"Ylo, zhi = %g, %g\n", domain->boxlo[1], domain->boxhi[1]);
-      fprintf(out,"Zlo, zhi = %g, %g\n", domain->boxlo[2], domain->boxhi[2]);
+      fprintf(out,"xlo, xhi = %g, %g\n", domain->boxlo[0], domain->boxhi[0]);
+      fprintf(out,"ylo, yhi = %g, %g\n", domain->boxlo[1], domain->boxhi[1]);
+      fprintf(out,"zlo, zhi = %g, %g\n", domain->boxlo[2], domain->boxhi[2]);
       if (domain->triclinic)
           fprintf(out,"Xy, xz, yz = %g, %g, %g\n",
                   domain->xy, domain->xz, domain->yz);
-- 
GitLab