Merge pull request #248 from akohlmey/collected-small-fixes

collected small bugfixes and updates

doc/src/Section_howto.txt

@@ -2729,7 +2729,7 @@ production runs and is only required during equilibration. This way one
 is consistent with literature (based on the code packages DL_POLY or
 GULP for instance).
 
-The mentioned energy transfer will typically lead to a a small drift
+The mentioned energy transfer will typically lead to a small drift
 in total energy over time.  This internal energy can be monitored
 using the "compute chunk/atom"_compute_chunk_atom.html and "compute
 temp/chunk"_compute_temp_chunk.html commands.  The internal kinetic
@@ -2830,7 +2830,7 @@ temp/drude"_compute_temp_drude.html. This requires also to use the
 command {comm_modify vel yes}.
 
 Short-range damping of the induced dipole interactions can be achieved
-using Thole functions through the the "pair style
+using Thole functions through the "pair style
 thole"_pair_thole.html in "pair_style hybrid/overlay"_pair_hybrid.html
 with a Coulomb pair style. It may be useful to use {coul/long/cs} or
 similar from the CORESHELL package if the core and Drude particle come

doc/src/Section_start.txt

@@ -706,7 +706,7 @@ future changes to LAMMPS.
 User packages, such as user-atc or user-omp, have been contributed by
 users, and always begin with the user prefix.  If they are a single
 command (single file), they are typically in the user-misc package.
-Otherwise, they are a a set of files grouped together which add a
+Otherwise, they are a set of files grouped together which add a
 specific functionality to the code.
 
 User packages don't necessarily meet the requirements of the standard

doc/src/compute_centro_atom.txt

@@ -114,7 +114,7 @@ local defects surrounding the central atom, as described above.  For
 the {axes yes} case, the vector components are also unitless, since
 they represent spatial directions.
 
-Here are typical centro-symmetry values, from a a nanoindentation
+Here are typical centro-symmetry values, from a nanoindentation
 simulation into gold (FCC).  These were provided by Jon Zimmerman
 (Sandia):
 

doc/src/compute_chunk_atom.txt

@@ -536,7 +536,7 @@ For the {bin/cylinder} style the details are as follows.  If {discard}
 is set to {yes}, an out-of-domain atom will have its chunk ID set to
 0.  If {discard} is set to {no}, the atom will have its chunk ID set
 to the first or last bin in both the radial and axis dimensions.  If
-{discard} is set to {mixed}, which is the default, the the radial
+{discard} is set to {mixed}, which is the default, the radial
 dimension is treated the same as for {discard} = no.  But for the axis
 dimensinon, it will only have its chunk ID set to the first or last
 bin if bins extend to the simulation box boundary in the axis

doc/src/compute_stress_atom.txt

@@ -60,7 +60,7 @@ produced by a small set of atoms (e.g. 4 atoms in a dihedral or 3
 atoms in a Tersoff 3-body interaction) is assigned in equal portions
 to each atom in the set.  E.g. 1/4 of the dihedral virial to each of
 the 4 atoms, or 1/3 of the fix virial due to SHAKE constraints applied
-to atoms in a a water molecule via the "fix shake"_fix_shake.html
+to atoms in a water molecule via the "fix shake"_fix_shake.html
 command.
 
 If no extra keywords are listed, all of the terms in this formula are

doc/src/dihedral_table.txt

@@ -154,7 +154,7 @@ radians instead of degrees.  (Note: This changes the way the forces
 are scaled in the 4th column of the data file.)
 
 The optional "CHECKU" keyword is followed by a filename.  This allows
-the user to save all of the the {Ntable} different entries in the
+the user to save all of the {Ntable} different entries in the
 interpolated energy table to a file to make sure that the interpolated
 function agrees with the user's expectations.  (Note: You can
 temporarily increase the {Ntable} parameter to a high value for this

doc/src/dump_modify.txt

@@ -165,7 +165,7 @@ extra buffering.
 
 :line
 
-The {element} keyword applies only to the the dump {cfg}, {xyz}, and
+The {element} keyword applies only to the dump {cfg}, {xyz}, and
 {image} styles.  It associates element names (e.g. H, C, Fe) with
 LAMMPS atom types.  See the list of element names at the bottom of
 this page.
@@ -574,7 +574,7 @@ e.g. its x-component of velocity if the atom-attribute "vx" was
 specified.
 
 The basic idea of a color map is that the atom-attribute will be
-within a range of values, and that range is associated with a a series
+within a range of values, and that range is associated with a series
 of colors (e.g. red, blue, green).  An atom's specific value (vx =
 -3.2) can then mapped to the series of colors (e.g. halfway between
 red and blue), and a specific color is determined via an interpolation

doc/src/fix_deform.txt

@@ -150,7 +150,7 @@ initial box length is 100 Angstroms, and V is 10 Angstroms/psec, then
 after 10 psec, the box length will have doubled.  After 20 psec, it
 will have tripled.
 
-The {erate} style changes a dimension of the the box at a "constant
+The {erate} style changes a dimension of the box at a "constant
 engineering strain rate".  The units of the specified strain rate are
 1/time.  See the "units"_units.html command for the time units
 associated with different choices of simulation units,

doc/src/fix_langevin.txt

@@ -237,7 +237,7 @@ described in the papers cited below, the purpose of this method is to
 enable longer timesteps to be used (up to the numerical stability
 limit of the integrator), while still producing the correct Boltzmann
 distribution of atom positions.  It is implemented within LAMMPS, by
-changing how the the random force is applied so that it is composed of
+changing how the random force is applied so that it is composed of
 the average of two random forces representing half-contributions from
 the previous and current time intervals.
 

doc/src/fix_lb_fluid.txt

@@ -233,7 +233,7 @@ present, the speed of sound squared is set equal to (1/3)*(dx/dt)^2.
 Setting a0 > (dx/dt)^2 is not allowed, as this may lead to
 instabilities.
 
-If the {noise} keyword is used, followed by a a positive temperature
+If the {noise} keyword is used, followed by a positive temperature
 value, and a positive integer random number seed, a thermal
 lattice-Boltzmann algorithm is used.  If {LBtype} is set equal to 1
 (i.e. the standard LB integrator is chosen), the thermal LB algorithm

doc/src/fix_rx.txt

@@ -87,7 +87,7 @@ end of each run. A positive value N means that the diagnostics are reported once
 per N time-steps.
 
 The diagnostics report the average # of integrator steps and RHS function evaluations
-and run-time per ODE as well as the the average/RMS/min/max per process. If the
+and run-time per ODE as well as the average/RMS/min/max per process. If the
 reporting frequency is 1, the RMS/min/max per ODE are also reported. The per ODE
 statistics can be used to adjust the tolerance and min/max step parameters. The
 statistics per MPI process can be useful to examine any load imbalance caused by the

doc/src/fix_smd_adjust_dt.txt

@@ -28,7 +28,7 @@ fix 1 all smd/adjust_dt 0.1 :pre
 The fix calculates a new stable time increment for use with the SMD time integrators.
 
 The stable time increment is based on multiple conditions. For the SPH pair styles, a
-CFL criterion (Courant, Friedrichs & Lewy, 1928) is evaluated, which determines the the speed of
+CFL criterion (Courant, Friedrichs & Lewy, 1928) is evaluated, which determines the speed of
 sound cannot propagate further than a typical spacing between particles within a single time step to ensure
 no information is lost. For the contact pair styles, a linear analysis of the pair potential determines a
 stable maximum time step.

doc/src/fix_srd.txt

@@ -101,7 +101,7 @@ particles move in the normal way via a time integration "fix"_fix.html
 with a short timestep dt.  SRD particles advect with a large timestep
 dt_SRD >= dt.
 
-If the {lamda} keyword is not specified, the the SRD temperature
+If the {lamda} keyword is not specified, the SRD temperature
 {Tsrd} is used in the above formula to compute lamda.  If the {lamda}
 keyword is specified, then the {Tsrd} setting is ignored and the above
 equation is used to compute the SRD temperature.

doc/src/fix_ttm.txt

@@ -107,7 +107,7 @@ specified as parameters to the fix.  The other quantities are derived.
 The form of the heat diffusion equation used here is almost the same
 as that in equation 6 of "(Duffy)"_#Duffy, with the exception that the
 electronic density is explicitly reprensented, rather than being part
-of the the specific heat parameter.
+of the specific heat parameter.
 
 Currently, fix ttm assumes that none of the user-supplied parameters
 will vary with temperature. Note that "(Duffy)"_#Duffy used a tanh()

doc/src/improper_style.txt

@@ -27,7 +27,7 @@ between quadruplets of atoms, which remain in force for the duration
 of the simulation.  The list of improper quadruplets is read in by a
 "read_data"_read_data.html or "read_restart"_read_restart.html command
 from a data or restart file.  Note that the ordering of the 4 atoms in
-an improper quadruplet determines the the definition of the improper
+an improper quadruplet determines the definition of the improper
 angle used in the formula for each style.  See the doc pages of
 individual styles for details.
 

doc/src/neb.txt

@@ -284,7 +284,7 @@ ID2 x2 y2 z2
 ...
 IDN xN yN zN :pre
 
-The fields are the the atom ID, followed by the x,y,z coordinates.
+The fields are the atom ID, followed by the x,y,z coordinates.
 The lines can be listed in any order.  Additional trailing information
 on the line is OK, such as a comment.
 

doc/src/next.txt

@@ -44,7 +44,7 @@ one value from their respective list of values.  A {file}-style
 variable reads the next line from its associated file.  An
 {atomfile}-style variable reads the next set of lines (one per atom)
 from its associated file.  {String-} or {atom}- or {equal}- or
-{world}-style variables cannot be used with the the next command,
+{world}-style variables cannot be used with the next command,
 since they only store a single value.
 
 When any of the variables in the next command has no more values, a

doc/src/pair_srp.txt

@@ -99,7 +99,7 @@ The optional {exclude} keyword determines if forces are computed
 between first neighbor (directly connected) bonds.  For a setting of
 {no}, first neighbor forces are computed; for {yes} they are not
 computed. A setting of {no} cannot be used with the {min} option for
-distance calculation because the the minimum distance between directly
+distance calculation because the minimum distance between directly
 connected bonds is zero.
 
 Pair style {srp} turns off normalization of thermodynamic properties

doc/src/prd.txt

@@ -214,7 +214,7 @@ when a correlated event occurs during the third stage of the loop
 listed above, i.e. when only one replica is running dynamics.
 
 When more than one replica detects an event at the end of the same
-event check (every {t_event} steps) during the the second stage, then
+event check (every {t_event} steps) during the second stage, then
 one of them is chosen at random.  The number of coincident events is
 the number of replicas that detected an event.  Normally, this value
 should be 1.  If it is often greater than 1, then either the number of
@@ -241,7 +241,7 @@ time was spent in each stage (dephasing, dynamics, quenching, etc).
 
 Any "dump files"_dump.html defined in the input script, will be
 written to during a PRD run at timesteps corresponding to both
-uncorrelated and correlated events.  This means the the requested dump
+uncorrelated and correlated events.  This means the requested dump
 frequency in the "dump"_dump.html command is ignored.  There will be
 one dump file (per dump command) created for all partitions.
 

doc/src/python.txt

@@ -188,7 +188,7 @@ is assumed to have been previously loaded by another python command.
 
 Note that the Python code that is loaded and run must contain a
 function with the specified {func} name.  To operate properly when
-later invoked, the the function code must match the {input} and
+later invoked, the function code must match the {input} and
 {return} and {format} keywords specified by the python command.
 Otherwise Python will generate an error.
 

doc/src/read_dump.txt

@@ -185,7 +185,7 @@ For dump files in {xyz} format, only the {x}, {y}, and {z} fields are
 supported.  The dump file does not store atom IDs, so these are
 assigned consecutively to the atoms as they appear in the dump file,
 starting from 1.  Thus you should insure that order of atoms is
-consistent from snapshot to snapshot in the the XYZ dump file.  See
+consistent from snapshot to snapshot in the XYZ dump file.  See
 the "dump_modify sort"_dump_modify.html command if the XYZ dump file
 was written by LAMMPS.
 
@@ -195,7 +195,7 @@ velocities, or their respective plugins may not support reading of
 velocities.  The molfile dump files do not store atom IDs, so these
 are assigned consecutively to the atoms as they appear in the dump
 file, starting from 1.  Thus you should insure that order of atoms are
-consistent from snapshot to snapshot in the the molfile dump file.
+consistent from snapshot to snapshot in the molfile dump file.
 See the "dump_modify sort"_dump_modify.html command if the dump file
 was written by LAMMPS.
 

doc/src/tad.txt

@@ -231,7 +231,7 @@ time was spent in each stage (NEB, dynamics, quenching, etc).
 
 Any "dump files"_dump.html defined in the input script will be written
 to during a TAD run at timesteps when an event is executed.  This
-means the the requested dump frequency in the "dump"_dump.html command
+means the requested dump frequency in the "dump"_dump.html command
 is ignored.  There will be one dump file (per dump command) created
 for all partitions.  The atom coordinates of the dump snapshot are
 those of the minimum energy configuration resulting from quenching

doc/src/tutorial_drude.txt

@@ -33,7 +33,7 @@ created by the surrounding particles.  Drude oscillators represent
 these dipoles by two fixed charges: the core (DC) and the Drude
 particle (DP) bound by a harmonic potential. The Drude particle can be
 thought of as the electron cloud whose center can be displaced from
-the position of the the corresponding nucleus.
+the position of the corresponding nucleus.
 
 The sum of the masses of a core-Drude pair should be the mass of the
 initial (unsplit) atom, \(m_C + m_D = m\).  The sum of their charges

src/KSPACE/ewald_disp.cpp

@@ -141,6 +141,7 @@ void EwaldDisp::init()
   init_coeffs();
   init_coeff_sums();
   if (function[0]) qsum_qsq();
+  else qsqsum = qsum = 0.0;
   natoms_original = atom->natoms;
 
   // turn off coulombic if no charge
@@ -152,6 +153,7 @@ void EwaldDisp::init()
   }
 
   double bsbsum = 0.0;
+  M2 = 0.0;
   if (function[1]) bsbsum = sum[1].x2;
   if (function[2]) bsbsum = sum[2].x2;
 
@@ -489,8 +491,9 @@ void EwaldDisp::init_coeffs()
     double **b = (double **) force->pair->extract("B",tmp);
     delete [] B;
     B = new double[n+1];
+    B[0] = 0.0;
     bytes += (n+1)*sizeof(double);
-    for (int i=0; i<=n; ++i) B[i] = sqrt(fabs(b[i][i]));
+    for (int i=1; i<=n; ++i) B[i] = sqrt(fabs(b[i][i]));
   }
   if (function[2]) {                                        // arithmetic 1/r^6
     double **epsilon = (double **) force->pair->extract("epsilon",tmp);
@@ -502,7 +505,9 @@ void EwaldDisp::init_coeffs()
     if (!(epsilon&&sigma))
       error->all(
           FLERR,"Epsilon or sigma reference not set by pair style in ewald/n");
-    for (int i=0; i<=n; ++i) {
+    for (int j=0; j<7; ++j)
+      *(bi++) = 0.0;
+    for (int i=1; i<=n; ++i) {
       eps_i = sqrt(epsilon[i][i]);
       sigma_i = sigma[i][i];
       sigma_n = 1.0;
@@ -523,6 +528,7 @@ void EwaldDisp::init_coeff_sums()
   Sum sum_local[EWALD_MAX_NSUMS];
 
   memset(sum_local, 0, EWALD_MAX_NSUMS*sizeof(Sum));
+  memset(sum,       0, EWALD_MAX_NSUMS*sizeof(Sum));
 
   // now perform qsum and qsq via parent qsum_qsq()
 

src/KSPACE/pair_buck_long_coul_long.cpp

@@ -51,8 +51,6 @@ PairBuckLongCoulLong::PairBuckLongCoulLong(LAMMPS *lmp) : Pair(lmp)
   dispersionflag = ewaldflag = pppmflag = 1;
   respa_enable = 1;
   writedata = 1;
-  ftable = NULL;
-  fdisptable = NULL;
 }
 
 /* ----------------------------------------------------------------------
@@ -134,6 +132,7 @@ PairBuckLongCoulLong::~PairBuckLongCoulLong()
     memory->destroy(offset);
   }
   if (ftable) free_tables();
+  if (fdisptable) free_disp_tables();
 }
 
 /* ----------------------------------------------------------------------
@@ -288,8 +287,8 @@ void PairBuckLongCoulLong::init_style()
   if (ewald_order&(1<<6)) g_ewald_6 = force->kspace->g_ewald_6;
   // setup force tables
 
-  if (ncoultablebits) init_tables(cut_coul,cut_respa);
-  if (ndisptablebits) init_tables_disp(cut_buck_global);
+  if (ncoultablebits && (ewald_order&(1<<1))) init_tables(cut_coul,cut_respa);
+  if (ndisptablebits && (ewald_order&(1<<6))) init_tables_disp(cut_buck_global);
 }
 
 /* ----------------------------------------------------------------------

src/KSPACE/pair_lj_long_coul_long.cpp

@@ -130,6 +130,7 @@ PairLJLongCoulLong::~PairLJLongCoulLong()
     memory->destroy(offset);
   }
   if (ftable) free_tables();
+  if (fdisptable) free_disp_tables();
 }
 
 /* ----------------------------------------------------------------------
@@ -282,8 +283,8 @@ void PairLJLongCoulLong::init_style()
 
   // setup force tables
 
-  if (ncoultablebits) init_tables(cut_coul,cut_respa);
-  if (ndisptablebits) init_tables_disp(cut_lj_global);
+  if (ncoultablebits && (ewald_order&(1<<1))) init_tables(cut_coul,cut_respa);
+  if (ndisptablebits && (ewald_order&(1<<6))) init_tables_disp(cut_lj_global);
 
 }
 

src/KSPACE/pppm_disp.cpp

@@ -1357,7 +1357,8 @@ void PPPMDisp::init_coeffs()				// local pair coeffs
   if (function[1]) {					// geometric 1/r^6
     double **b = (double **) force->pair->extract("B",tmp);
     B = new double[n+1];
-    for (int i=0; i<=n; ++i) B[i] = sqrt(fabs(b[i][i]));
+    B[0] = 0.0;
+    for (int i=1; i<=n; ++i) B[i] = sqrt(fabs(b[i][i]));
   }
   if (function[2]) {					// arithmetic 1/r^6
     //cannot use epsilon, because this has not been set yet

src/USER-OMP/pppm_disp_omp.cpp

@@ -50,6 +50,13 @@ PPPMDispOMP::PPPMDispOMP(LAMMPS *lmp, int narg, char **arg) :
   suffix_flag |= Suffix::OMP;
 }
 
+/* ---------------------------------------------------------------------- */
+
+PPPMDispOMP::~PPPMDispOMP()
+{
+  deallocate();
+}
+
 /* ----------------------------------------------------------------------
    allocate memory that depends on # of K-vectors and order
 ------------------------------------------------------------------------- */

src/USER-OMP/pppm_disp_omp.h

@@ -28,7 +28,7 @@ namespace LAMMPS_NS {
   class PPPMDispOMP : public PPPMDisp, public ThrOMP {
  public:
   PPPMDispOMP(class LAMMPS *, int, char **);
-  virtual ~PPPMDispOMP () {};
+  virtual ~PPPMDispOMP ();
   virtual void compute(int, int);
 
  protected:

src/output.cpp

@@ -749,6 +749,7 @@ void Output::create_restart(int narg, char **arg)
     } else restart_every_single = every;
 
     int n = strlen(arg[1]) + 3;
+    delete [] restart1;
     restart1 = new char[n];
     strcpy(restart1,arg[1]);
     if (strchr(restart1,'*') == NULL) strcat(restart1,".*");
@@ -765,6 +766,8 @@ void Output::create_restart(int narg, char **arg)
       restart_every_double = 0;
     } else restart_every_double = every;
 
+    delete [] restart2a;
+    delete [] restart2b;
     restart_toggle = 0;
     int n = strlen(arg[1]) + 3;
     restart2a = new char[n];

src/pair.cpp

@@ -87,6 +87,8 @@ Pair::Pair(LAMMPS *lmp) : Pointers(lmp)
   ndisptablebits = 12;
   tabinner = sqrt(2.0);
   tabinner_disp = sqrt(2.0);
+  ftable = NULL;
+  fdisptable = NULL;
 
   allocated = 0;
   suffix_flag = Suffix::NONE;

src/region_cylinder.cpp

@@ -36,8 +36,7 @@ RegCylinder::RegCylinder(LAMMPS *lmp, int narg, char **arg) :
 
   // check open face settings
 
-  if (openflag && (open_faces[2] || open_faces[3] ||
-                   open_faces[4] || open_faces[5]))
+  if (openflag && (open_faces[3] || open_faces[4] || open_faces[5]))
     error->all(FLERR,"Invalid region cylinder open setting");
 
   if (strcmp(arg[2],"x") && strcmp(arg[2],"y") && strcmp(arg[2],"z"))

src/universe.cpp

@@ -163,17 +163,49 @@ void Universe::add_world(char *str)
   int n,nper;
   char *ptr;
 
-  if (str == NULL) {
-    n = 1;
-    nper = nprocs;
-  } else if ((ptr = strchr(str,'x')) != NULL) {
-    *ptr = '\0';
-    n = atoi(str);
-    nper = atoi(ptr+1);
-  } else {
-    n = 1;
-    nper = atoi(str);
-  }
+  n = 1;
+  if (str != NULL) {
+
+    // check for valid partition argument
+
+    bool valid = true;
+
+    // str may not be empty and may only consist of digits or 'x'
+
+    int len = strlen(str);
+    if (len < 1) valid = false;
+    for (int i=0; i < len; ++i)
+      if (isdigit(str[i]) || str[i] == 'x') continue;
+      else valid = false;
+
+    if (valid) {
+      if ((ptr = strchr(str,'x')) != NULL) {
+
+        // 'x' may not be the first or last character
+
+        if (ptr == str) {
+          valid = false;
+        } else if (strlen(str) == len-1) {
+          valid = false;
+        } else {
+          *ptr = '\0';
+          n = atoi(str);
+          nper = atoi(ptr+1);
+          *ptr = 'x';
+        }
+      } else nper = atoi(str);
+    }
+
+    // require minimum of 1 partition with 1 processor
+
+    if (n < 1 || nper < 1) valid = false;
+
+    if (!valid) {
+      char msg[128];
+      sprintf(msg,"Invalid partition string '%s'",str);
+      error->universe_all(FLERR,msg);
+    }
+  } else nper = nprocs;
 
   memory->grow(procs_per_world,nworlds+n,"universe:procs_per_world");
   memory->grow(root_proc,(nworlds+n),"universe:root_proc");