Merge pull request #593 from akohlmey/collected-bug-fixes-and-small-updates

Collected bug fixes and small updates

examples/README

@@ -58,6 +58,7 @@ These are the sample problems and their output in the various
 sub-directories:
 
 accelerate: use of all the various accelerator packages
+airebo:   example for using AIREBO and AIREBO-M
 balance:  dynamic load balancing, 2d system
 body:     body particles, 2d system
 cmap:     CMAP 5-body contributions to CHARMM force field

src/RIGID/fix_shake.cpp

@@ -1617,6 +1617,12 @@ void FixShake::shake3(int m)
 
     lamda01 = lamda01_new;
     lamda02 = lamda02_new;
+
+    // stop iterations before we have a floating point overflow
+    // max double is < 1.0e308, so 1e150 is a reasonable cutoff
+
+    if (fabs(lamda01) > 1e150 || fabs(lamda02) > 1e150) done = 1;
+
     niter++;
   }
 
@@ -1854,6 +1860,13 @@ void FixShake::shake4(int m)
     lamda01 = lamda01_new;
     lamda02 = lamda02_new;
     lamda03 = lamda03_new;
+
+    // stop iterations before we have a floating point overflow
+    // max double is < 1.0e308, so 1e150 is a reasonable cutoff
+
+    if (fabs(lamda01) > 1e150 || fabs(lamda02) > 1e150
+        || fabs(lamda03) > 1e150) done = 1;
+
     niter++;
   }
 
@@ -2097,6 +2110,13 @@ void FixShake::shake3angle(int m)
     lamda01 = lamda01_new;
     lamda02 = lamda02_new;
     lamda12 = lamda12_new;
+
+    // stop iterations before we have a floating point overflow
+    // max double is < 1.0e308, so 1e150 is a reasonable cutoff
+
+    if (fabs(lamda01) > 1e150 || fabs(lamda02) > 1e150
+        || fabs(lamda12) > 1e150) done = 1;
+
     niter++;
   }
 

src/USER-DRUDE/pair_lj_cut_thole_long.cpp

@@ -417,7 +417,7 @@ double PairLJCutTholeLong::init_one(int i, int j)
   lj3[i][j] = 4.0 * epsilon[i][j] * pow(sigma[i][j],12.0);
   lj4[i][j] = 4.0 * epsilon[i][j] * pow(sigma[i][j],6.0);
 
-  if (offset_flag) {
+  if (offset_flag && (cut_lj[i][j] > 0.0)) {
     double ratio = sigma[i][j] / cut_lj[i][j];
     offset[i][j] = 4.0 * epsilon[i][j] * (pow(ratio,12.0) - pow(ratio,6.0));
   } else offset[i][j] = 0.0;

src/USER-MISC/pair_morse_smooth_linear.cpp

@@ -296,7 +296,6 @@ void PairMorseSmoothLinear::read_restart(FILE *fp)
 void PairMorseSmoothLinear::write_restart_settings(FILE *fp)
 {
   fwrite(&cut_global,sizeof(double),1,fp);
-  // fwrite(&offset_flag,sizeof(int),1,fp);
   fwrite(&mix_flag,sizeof(int),1,fp);
 }
 
@@ -308,11 +307,9 @@ void PairMorseSmoothLinear::read_restart_settings(FILE *fp)
 {
   if (comm->me == 0) {
     fread(&cut_global,sizeof(double),1,fp);
-    // fread(&offset_flag,sizeof(int),1,fp);
     fread(&mix_flag,sizeof(int),1,fp);
   }
   MPI_Bcast(&cut_global,1,MPI_DOUBLE,0,world);
-  // MPI_Bcast(&offset_flag,1,MPI_INT,0,world);
   MPI_Bcast(&mix_flag,1,MPI_INT,0,world);
 }
 

src/USER-REAXC/reaxc_list.cpp

@@ -37,7 +37,7 @@ int Make_List(int n, int num_intrs, int type, reax_list *l, MPI_Comm comm)
   l->num_intrs = num_intrs;
 
   if (l->index) sfree(l->index, "list:index");
-  if (l->end_index) sfree(l->index, "list:end_index");
+  if (l->end_index) sfree(l->end_index, "list:end_index");
   l->index = (int*) smalloc( n * sizeof(int), "list:index", comm );
   l->end_index = (int*) smalloc( n * sizeof(int), "list:end_index", comm );
 

src/USER-TALLY/compute_force_tally.cpp

@@ -215,7 +215,7 @@ void ComputeForceTally::compute_peratom()
 
 double ComputeForceTally::memory_usage()
 {
-  double bytes = nmax*size_peratom_cols * sizeof(double);
+  double bytes = (nmax < 0) ? 0 : nmax*size_peratom_cols * sizeof(double);
   return bytes;
 }
 

src/USER-TALLY/compute_heat_flux_tally.cpp

@@ -275,7 +275,7 @@ void ComputeHeatFluxTally::compute_vector()
 
 double ComputeHeatFluxTally::memory_usage()
 {
-  double bytes = nmax*comm_reverse * sizeof(double);
+  double bytes = (nmax < 0) ? 0 : nmax*comm_reverse * sizeof(double);
   return bytes;
 }
 

src/USER-TALLY/compute_pe_tally.cpp

@@ -199,7 +199,7 @@ void ComputePETally::compute_peratom()
 
 double ComputePETally::memory_usage()
 {
-  double bytes = nmax*size_peratom_cols * sizeof(double);
+  double bytes = (nmax < 0) ? 0 : nmax*size_peratom_cols * sizeof(double);
   return bytes;
 }
 

src/USER-TALLY/compute_stress_tally.cpp

@@ -242,7 +242,7 @@ void ComputeStressTally::compute_peratom()
 
 double ComputeStressTally::memory_usage()
 {
-  double bytes = nmax*size_peratom_cols * sizeof(double);
+  double bytes = (nmax < 0) ? 0 : nmax*size_peratom_cols * sizeof(double);
   return bytes;
 }
 

src/dump_local.cpp

@@ -49,6 +49,9 @@ DumpLocal::DumpLocal(LAMMPS *lmp, int narg, char **arg) :
   nevery = force->inumeric(FLERR,arg[3]);
   if (nevery <= 0) error->all(FLERR,"Illegal dump local command");
 
+  if (binary)
+    error->all(FLERR,"Binary files are not supported with dump local");
+
   nfield = narg - 5;
 
   // expand args if any have wildcard character "*"

src/main.cpp

@@ -18,6 +18,10 @@
 #include <stdio.h>
 #include <stdlib.h>
 
+#if defined(LAMMPS_TRAP_FPE) && defined(_GNU_SOURCE)
+#include <fenv.h>
+#endif
+
 using namespace LAMMPS_NS;
 
 /* ----------------------------------------------------------------------
@@ -28,6 +32,18 @@ int main(int argc, char **argv)
 {
   MPI_Init(&argc,&argv);
 
+// enable trapping selected floating point exceptions.
+// this uses GNU extensions and is only tested on Linux
+// therefore we make it depend on -D_GNU_SOURCE, too.
+
+#if defined(LAMMPS_TRAP_FPE) && defined(_GNU_SOURCE)
+  fesetenv(FE_NOMASK_ENV);
+  fedisableexcept(FE_ALL_EXCEPT);
+  feenableexcept(FE_DIVBYZERO);
+  feenableexcept(FE_INVALID);
+  feenableexcept(FE_OVERFLOW);
+#endif
+
 #ifdef LAMMPS_EXCEPTIONS
   try {
     LAMMPS *lammps = new LAMMPS(argc,argv,MPI_COMM_WORLD);