sync with GH

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

lib/colvars/colvarmodule.cpp

@@ -137,6 +137,9 @@ int colvarmodule::parse_config(std::string &conf)
   cvm::log("Collective variables module (re)initialized.\n");
   cvm::log(cvm::line_marker);
 
+  // update any necessary proxy data
+  proxy->setup();
+
   if (cv_traj_os.is_open()) {
     // configuration might have changed, better redo the labels
     write_traj_label(cv_traj_os);

lib/colvars/colvarmodule.h

@@ -4,7 +4,7 @@
 #define COLVARMODULE_H
 
 #ifndef COLVARS_VERSION
-#define COLVARS_VERSION "2016-09-03"
+#define COLVARS_VERSION "2016-09-14"
 #endif
 
 #ifndef COLVARS_DEBUG

python/lammps.py

@@ -436,6 +436,42 @@ class Atom2D(Atom):
             self.lmp.eval("fy[%d]" % self.index))
 
 
+def get_thermo_data(output):
+    """ traverse output of runs and extract thermo data columns """
+    if isinstance(output, str):
+        lines = output.splitlines()
+    else:
+        lines = output
+
+    runs = []
+    columns = []
+    in_run = False
+
+    for line in lines:
+        if line.startswith("Memory usage per processor"):
+            in_run = True
+        elif in_run and len(columns) == 0:
+            # first line after memory usage are column names
+            columns = line.split()
+
+            current_run = {}
+
+            for col in columns:
+                current_run[col] = []
+
+        elif line.startswith("Loop time of "):
+            in_run = False
+            columns = None
+            thermo_data = namedtuple('ThermoData', list(current_run.keys()))(*list(current_run.values()))
+            r = {'thermo' : thermo_data }
+            runs.append(namedtuple('Run', list(r.keys()))(*list(r.values())))
+        elif in_run and len(columns) > 0:
+            values = [float(x) for x in line.split()]
+
+            for i, col in enumerate(columns):
+                current_run[col].append(values[i])
+    return runs
+
 class PyLammps(object):
   """
   More Python-like wrapper for LAMMPS (e.g., for iPython)
@@ -454,6 +490,8 @@ class PyLammps(object):
     else:
       self.lmp = lammps(name=name,cmdargs=cmdargs,ptr=None,comm=comm)
     print("LAMMPS output is captured by PyLammps wrapper")
+    self._cmd_history = []
+    self.runs = []
 
   def __del__(self):
     if self.lmp: self.lmp.close()
@@ -469,8 +507,26 @@ class PyLammps(object):
   def file(self,file):
     self.lmp.file(file)
 
+  def write_script(self,filename):
+    """ Write LAMMPS script file containing all commands executed up until now """
+    with open(filename, "w") as f:
+      for cmd in self._cmd_history:
+        f.write("%s\n" % cmd)
+
   def command(self,cmd):
     self.lmp.command(cmd)
+    self._cmd_history.append(cmd)
+
+  def run(self, *args, **kwargs):
+    output = self.__getattr__('run')(*args, **kwargs)
+    self.runs += get_thermo_data(output)
+    return output
+
+  @property
+  def last_run(self):
+    if len(self.runs) > 0:
+        return self.runs[-1]
+    return None
 
   @property
   def atoms(self):
@@ -643,7 +699,7 @@ class PyLammps(object):
       cmd_args = [name] + [str(x) for x in args]
 
       with OutputCapture() as capture:
-        self.lmp.command(' '.join(cmd_args))
+        self.command(' '.join(cmd_args))
         output = capture.output
 
       if 'verbose' in kwargs and kwargs['verbose']:

src/MEAM/pair_meam.cpp

@@ -35,12 +35,13 @@ using namespace LAMMPS_NS;
 #define MAXLINE 1024
 
 enum{FCC,BCC,HCP,DIM,DIAMOND,B1,C11,L12,B2};
-int nkeywords = 21;
-const char *keywords[] = {"Ec","alpha","rho0","delta","lattce",
-                          "attrac","repuls","nn2","Cmin","Cmax","rc","delr",
-                          "augt1","gsmooth_factor","re","ialloy",
-                          "mixture_ref_t","erose_form","zbl",
-                          "emb_lin_neg","bkgd_dyn"};
+static const int nkeywords = 21;
+static const char *keywords[] = {
+  "Ec","alpha","rho0","delta","lattce",
+  "attrac","repuls","nn2","Cmin","Cmax","rc","delr",
+  "augt1","gsmooth_factor","re","ialloy",
+  "mixture_ref_t","erose_form","zbl",
+  "emb_lin_neg","bkgd_dyn"};
 
 /* ---------------------------------------------------------------------- */
 

src/Make.py

@@ -26,7 +26,7 @@ libclasses = ("atc","awpmd","colvars","cuda","gpu","h5md",
 buildclasses = ("intel","kokkos")
 makeclasses = ("cc","flags","mpi","fft","jpg","png")
 
-setargs = ("gzip","#gzip","ffmpeg","#ffmpeg","smallbig","bigbig","smallsmall")
+setargs = ("gzip","#gzip","ffmpeg","#ffmpeg","smallbig","bigbig","smallsmall","exceptions","#exceptions")
 actionargs = ("lib-all","file","clean","exe")
 
 gpubuildflag = 0
@@ -459,6 +459,8 @@ class Actions(object):
             make.delvar("LMP_INC","-DLAMMPS_SMALLBIG")
             make.delvar("LMP_INC","-DLAMMPS_BIGBIG")
             make.addvar("LMP_INC","-DLAMMPS_SMALLSMALL")
+          elif one == "exceptions": make.addvar("LMP_INC","-DLAMMPS_EXCEPTIONS")
+          elif one == "#exception": make.delvar("LMP_INC","-DLAMMPS_EXCEPTIONS")
           
       # add FFT, JPG, PNG settings
 
@@ -809,7 +811,7 @@ class Packages(object):
 
     original = {}
     tmp = "cd %s; make ps" % dir.src
-    output = subprocess.check_output(tmp,stderr=subprocess.STDOUT,shell=True).decode()
+    output = subprocess.check_output(tmp,stderr=subprocess.STDOUT,shell=True).decode().split('\n')
     pattern = "Installed\s+(\w+): package (\S+)"
     for line in output:
       m = re.search(pattern,line)
@@ -966,13 +968,14 @@ class Settings(object):
   def help(self):
     return """
 -s set1 set2 ...
-  possible settings = gzip #gzip ffmpeg #ffmpeg smallbig bigbig smallsmall
+  possible settings = gzip #gzip ffmpeg #ffmpeg smallbig bigbig smallsmall exceptions #exceptions
   alter LAMMPS ifdef settings in Makefile.auto
     only happens if new Makefile.auto is created by use of "file" action
   gzip and #gzip turn on/off LAMMPS_GZIP setting
   ffmpeg and #ffmpeg turn on/off LAMMPS_FFMPEG setting
   smallbig, bigbig, smallsmall turn on LAMMPS_SMALLBIG, etc
     and turn off other two
+  exceptions and #exceptions turn on/off LAMMPS_EXCEPTIONS setting
 """
   
   def check(self):

src/USER-MISC/pair_lj_sf_dipole_sf.cpp

@@ -31,6 +31,8 @@
 
 using namespace LAMMPS_NS;
 
+static int warn_single = 0;
+
 /* ---------------------------------------------------------------------- */
 
 PairLJSFDipoleSF::PairLJSFDipoleSF(LAMMPS *lmp) : Pair(lmp)
@@ -87,7 +89,6 @@ void PairLJSFDipoleSF::compute(int eflag, int vflag)
   double **torque = atom->torque;
   int *type = atom->type;
   int nlocal = atom->nlocal;
-// The global scaling parameters aren't used anymore
   double *special_coul = force->special_coul;
   double *special_lj = force->special_lj;
   int newton_pair = force->newton_pair;
@@ -355,7 +356,7 @@ void PairLJSFDipoleSF::settings(int narg, char **arg)
 
 void PairLJSFDipoleSF::coeff(int narg, char **arg)
 {
-  if (narg < 4 || narg > 7)
+  if (narg < 4 || narg > 8)
     error->all(FLERR,"Incorrect args for pair coefficients");
   if (!allocated) allocate();
 
@@ -366,13 +367,27 @@ void PairLJSFDipoleSF::coeff(int narg, char **arg)
   double epsilon_one = force->numeric(FLERR,arg[2]);
   double sigma_one = force->numeric(FLERR,arg[3]);
 
-  double scale_one = 1.0;
-  if (narg >= 5) scale_one = force->numeric(FLERR,arg[4]);
-
   double cut_lj_one = cut_lj_global;
   double cut_coul_one = cut_coul_global;
-  if (narg >= 6) cut_coul_one = cut_lj_one = force->numeric(FLERR,arg[5]);
-  if (narg == 7) cut_coul_one = force->numeric(FLERR,arg[6]);
+  double scale_one = 1.0;
+  int iarg = 4;
+
+  if ((narg > iarg) && (strcmp(arg[iarg],"scale") != 0)) {
+    cut_coul_one = cut_lj_one = force->numeric(FLERR,arg[iarg]);
+    ++iarg;
+  }
+  if ((narg > iarg) && (strcmp(arg[iarg],"scale") != 0)) {
+    cut_coul_one = force->numeric(FLERR,arg[iarg]);
+    ++iarg;
+  }
+  if (narg > iarg) {
+    if (strcmp(arg[iarg],"scale") == 0) {
+      scale_one = force->numeric(FLERR,arg[iarg+1]);
+      iarg += 2;
+    } else error->all(FLERR,"Incorrect args for pair coefficients");
+  }
+  if (iarg != narg)
+    error->all(FLERR,"Incorrect args for pair coefficients");
 
   int count = 0;
   for (int i = ilo; i <= ihi; i++) {
@@ -532,6 +547,12 @@ double PairLJSFDipoleSF::single(int i, int j, int itype, int jtype, double rsq,
   double *q = atom->q;
   double **mu = atom->mu;
 
+  if (!warn_single) {
+    warn_single = 1;
+    if (comm->me == 0) {
+      error->warning(FLERR,"Single method for lj/sf/dipole/sf does not compute forces");
+    }
+  }
   qtmp = q[i];
   xtmp = x[i][0];
   ytmp = x[i][1];

src/USER-OMP/pair_lj_sf_dipole_sf_omp.cpp

@@ -247,7 +247,7 @@ void PairLJSFDipoleSFOMP::eval(int iifrom, int iito, ThrData * const thr)
 
         // total force
 
-        fq = factor_coul*qqrd2e;
+        fq = factor_coul*qqrd2e*scale[itype][jtype];
         fx = fq*forcecoulx + delx*forcelj;
         fy = fq*forcecouly + dely*forcelj;
         fz = fq*forcecoulz + delz*forcelj;
@@ -272,7 +272,7 @@ void PairLJSFDipoleSFOMP::eval(int iifrom, int iito, ThrData * const thr)
 
         if (EFLAG) {
           if (rsq < cut_coulsq[itype][jtype]) {
-            ecoul = (1.0-sqrt(rsq)/sqrt(cut_coulsq[itype][jtype]));
+            ecoul = (1.0-sqrt(rsq/cut_coulsq[itype][jtype]));
             ecoul *= ecoul;
             ecoul *= qtmp * q[j] * rinv;
             if (mu[i].w > 0.0 && mu[j].w > 0.0)
@@ -281,7 +281,7 @@ void PairLJSFDipoleSFOMP::eval(int iifrom, int iito, ThrData * const thr)
               ecoul += -q[j] * r3inv * pqfac * pidotr;
             if (mu[j].w > 0.0 && qtmp != 0.0)
               ecoul += qtmp * r3inv * qpfac * pjdotr;
-            ecoul *= factor_coul*qqrd2e;
+            ecoul *= factor_coul*qqrd2e*scale[itype][jtype];
           } else ecoul = 0.0;
 
           if (rsq < cut_ljsq[itype][jtype]) {

src/fix_temp_csvr.cpp

@@ -49,7 +49,12 @@ double FixTempCSVR::gamdev(const int ia)
     x=1.0;
     for (j=1; j<=ia; j++)
       x *= random->uniform();
-    x = -log(x);
+
+    // make certain, that -log() doesn't overflow.
+    if (x < 2.2250759805e-308)
+      x = 708.4;
+    else
+      x = -log(x);
   } else {
   restart:
     do {