From 09a005414bf530cb5540be52bef14007f48b63d5 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sat, 6 Jun 2020 17:27:10 -0400 Subject: [PATCH 01/19] recover compilation of the KIM package --- src/KIM/kim_init.cpp | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/KIM/kim_init.cpp b/src/KIM/kim_init.cpp index 5328769d22..9f16a2bcdb 100644 --- a/src/KIM/kim_init.cpp +++ b/src/KIM/kim_init.cpp @@ -57,6 +57,8 @@ ------------------------------------------------------------------------- */ #include "kim_init.h" +#include "fix_store_kim.h" +#include "kim_units.h" #include #include #include @@ -71,8 +73,7 @@ #include "input.h" #include "variable.h" #include "citeme.h" -#include "fix_store_kim.h" -#include "kim_units.h" +#include "utils.h" extern "C" { #include "KIM_SimulatorHeaders.h" From 586f2c00b026945edbb5173dcb9ce9cff737e045 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 10:13:35 -0400 Subject: [PATCH 02/19] add missing linefeed character --- src/force.cpp | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/force.cpp b/src/force.cpp index 2827c43268..840d7e31ca 100644 --- a/src/force.cpp +++ b/src/force.cpp @@ -1019,7 +1019,8 @@ FILE *Force::open_potential(const char *name) date = utils::get_potential_date(filepath, "potential"); if(!date.empty()) { - utils::logmesg(lmp, fmt::format("Reading potential file {} with DATE: {}", name, date)); + utils::logmesg(lmp, fmt::format("Reading potential file {} " + "with DATE: {}\n", name, date)); } return fopen(filepath.c_str(), "r"); } From edc7237f1527a93a7f2b1035ff3c7a391114cfb2 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 10:46:35 -0400 Subject: [PATCH 03/19] more thorough checking if the reaxff force field is consistent. --- src/USER-REAXC/reaxc_ffield.cpp | 26 +++++++++++++++++++------- 1 file changed, 19 insertions(+), 7 deletions(-) diff --git a/src/USER-REAXC/reaxc_ffield.cpp b/src/USER-REAXC/reaxc_ffield.cpp index c8e097eb1c..b865b48ad3 100644 --- a/src/USER-REAXC/reaxc_ffield.cpp +++ b/src/USER-REAXC/reaxc_ffield.cpp @@ -337,6 +337,8 @@ char Read_Force_Field( FILE *fp, reax_interaction *reax, j = atoi(tmp[0]) - 1; k = atoi(tmp[1]) - 1; + if ((c < 10) || (j < 0) || (k < 0)) + control->error_ptr->all(FLERR, "Inconsistent force field file"); if (j < reax->num_atom_types && k < reax->num_atom_types) { @@ -361,6 +363,8 @@ char Read_Force_Field( FILE *fp, reax_interaction *reax, /* line 2 */ fgets(s,MAX_LINE,fp); c=Tokenize(s,&tmp); + if (c < 8) + control->error_ptr->all(FLERR, "Inconsistent force field file"); val = atof(tmp[0]); reax->tbp[j][k].p_be2 = val; reax->tbp[k][j].p_be2 = val; @@ -491,6 +495,9 @@ char Read_Force_Field( FILE *fp, reax_interaction *reax, j = atoi(tmp[0]) - 1; k = atoi(tmp[1]) - 1; + if ((c < (lgflag ? 9 : 8)) || (j < 0) || (k < 0)) + control->error_ptr->all(FLERR, "Inconsistent force field file"); + if (j < reax->num_atom_types && k < reax->num_atom_types) { val = atof(tmp[2]); if (val > 0.0) { @@ -528,10 +535,12 @@ char Read_Force_Field( FILE *fp, reax_interaction *reax, reax->tbp[k][j].r_pp = val; } - val = atof(tmp[8]); - if (val >= 0.0) { - reax->tbp[j][k].lgcij = val; - reax->tbp[k][j].lgcij = val; + if (lgflag) { + val = atof(tmp[8]); + if (val >= 0.0) { + reax->tbp[j][k].lgcij = val; + reax->tbp[k][j].lgcij = val; + } } } } @@ -552,6 +561,8 @@ char Read_Force_Field( FILE *fp, reax_interaction *reax, j = atoi(tmp[0]) - 1; k = atoi(tmp[1]) - 1; m = atoi(tmp[2]) - 1; + if ((c < 10) || (j < 0) || (k < 0) || (m < 0)) + control->error_ptr->all(FLERR, "Inconsistent force field file"); if (j < reax->num_atom_types && k < reax->num_atom_types && m < reax->num_atom_types) { @@ -611,6 +622,8 @@ char Read_Force_Field( FILE *fp, reax_interaction *reax, k = atoi(tmp[1]) - 1; m = atoi(tmp[2]) - 1; n = atoi(tmp[3]) - 1; + if ((c < 9) || (k < 0) || (m < 0)) + control->error_ptr->all(FLERR, "Inconsistent force field file"); if (j >= 0 && n >= 0) { // this means the entry is not in compact form if (j < reax->num_atom_types && k < reax->num_atom_types && @@ -667,8 +680,6 @@ char Read_Force_Field( FILE *fp, reax_interaction *reax, } } - - /* next line is number of hydrogen bond params and some comments */ fgets( s, MAX_LINE, fp ); c = Tokenize( s, &tmp ); @@ -686,7 +697,8 @@ char Read_Force_Field( FILE *fp, reax_interaction *reax, j = atoi(tmp[0]) - 1; k = atoi(tmp[1]) - 1; m = atoi(tmp[2]) - 1; - + if ((c < 7) || (j < 0) || (k < 0) || (m < 0)) + control->error_ptr->all(FLERR, "Inconsistent force field file"); if (j < reax->num_atom_types && m < reax->num_atom_types) { val = atof(tmp[3]); From cee7cd5fe9e7390caa7fee5bba62711f70f67459 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 11:44:42 -0400 Subject: [PATCH 04/19] consolidate enumerator for per-atom array data types --- src/atom.cpp | 2 - src/atom.h | 1 + src/atom_vec.cpp | 146 +++++++++++++++++++++++------------------------ 3 files changed, 73 insertions(+), 76 deletions(-) diff --git a/src/atom.cpp b/src/atom.cpp index f47906f203..f6379213ea 100644 --- a/src/atom.cpp +++ b/src/atom.cpp @@ -48,8 +48,6 @@ using namespace MathConst; #define DELTA_PERATOM 64 #define EPSILON 1.0e-6 -enum{DOUBLE,INT,BIGINT}; - /* ---------------------------------------------------------------------- */ Atom::Atom(LAMMPS *lmp) : Pointers(lmp) diff --git a/src/atom.h b/src/atom.h index fca682aff6..e20ee0596d 100644 --- a/src/atom.h +++ b/src/atom.h @@ -24,6 +24,7 @@ class Atom : protected Pointers { public: char *atom_style; class AtomVec *avec; + enum{DOUBLE,INT,BIGINT}; // atom counts diff --git a/src/atom_vec.cpp b/src/atom_vec.cpp index b9a8d6b09d..29f2526751 100644 --- a/src/atom_vec.cpp +++ b/src/atom_vec.cpp @@ -33,8 +33,6 @@ using namespace MathConst; #define DELTA 16384 #define DELTA_BONUS 8192 -enum{DOUBLE,INT,BIGINT}; - /* ---------------------------------------------------------------------- */ AtomVec::AtomVec(LAMMPS *lmp) : Pointers(lmp) @@ -214,7 +212,7 @@ void AtomVec::grow(int n) datatype = mgrow.datatype[i]; cols = mgrow.cols[i]; const int nthreads = threads[i] ? comm->nthreads : 1; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) memory->grow(*((double **) pdata),nmax*nthreads,"atom:dvec"); else if (cols > 0) @@ -223,7 +221,7 @@ void AtomVec::grow(int n) maxcols = *(mgrow.maxcols[i]); memory->grow(*((double ***) pdata),nmax*nthreads,maxcols,"atom:darray"); } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) memory->grow(*((int **) pdata),nmax*nthreads,"atom:ivec"); else if (cols > 0) @@ -232,7 +230,7 @@ void AtomVec::grow(int n) maxcols = *(mgrow.maxcols[i]); memory->grow(*((int ***) pdata),nmax*nthreads,maxcols,"atom:iarray"); } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) memory->grow(*((bigint **) pdata),nmax*nthreads,"atom:bvec"); else if (cols > 0) @@ -275,7 +273,7 @@ void AtomVec::copy(int i, int j, int delflag) pdata = mcopy.pdata[n]; datatype = mcopy.datatype[n]; cols = mcopy.cols[n]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); vec[j] = vec[i]; @@ -292,7 +290,7 @@ void AtomVec::copy(int i, int j, int delflag) for (m = 0; m < ncols; m++) array[j][m] = array[i][m]; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); vec[j] = vec[i]; @@ -309,7 +307,7 @@ void AtomVec::copy(int i, int j, int delflag) for (m = 0; m < ncols; m++) array[j][m] = array[i][m]; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); vec[j] = vec[i]; @@ -377,7 +375,7 @@ int AtomVec::pack_comm(int n, int *list, double *buf, pdata = mcomm.pdata[nn]; datatype = mcomm.datatype[nn]; cols = mcomm.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); for (i = 0; i < n; i++) { @@ -392,7 +390,7 @@ int AtomVec::pack_comm(int n, int *list, double *buf, buf[m++] = array[j][mm]; } } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); for (i = 0; i < n; i++) { @@ -407,7 +405,7 @@ int AtomVec::pack_comm(int n, int *list, double *buf, buf[m++] = ubuf(array[j][mm]).d; } } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); for (i = 0; i < n; i++) { @@ -498,7 +496,7 @@ int AtomVec::pack_comm_vel(int n, int *list, double *buf, pdata = mcomm_vel.pdata[nn]; datatype = mcomm_vel.datatype[nn]; cols = mcomm_vel.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); for (i = 0; i < n; i++) { @@ -513,7 +511,7 @@ int AtomVec::pack_comm_vel(int n, int *list, double *buf, buf[m++] = array[j][mm]; } } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); for (i = 0; i < n; i++) { @@ -528,7 +526,7 @@ int AtomVec::pack_comm_vel(int n, int *list, double *buf, buf[m++] = ubuf(array[j][mm]).d; } } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); for (i = 0; i < n; i++) { @@ -572,7 +570,7 @@ void AtomVec::unpack_comm(int n, int first, double *buf) pdata = mcomm.pdata[nn]; datatype = mcomm.datatype[nn]; cols = mcomm.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); for (i = first; i < last; i++) @@ -583,7 +581,7 @@ void AtomVec::unpack_comm(int n, int first, double *buf) for (mm = 0; mm < cols; mm++) array[i][mm] = buf[m++]; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); for (i = first; i < last; i++) @@ -594,7 +592,7 @@ void AtomVec::unpack_comm(int n, int first, double *buf) for (mm = 0; mm < cols; mm++) array[i][mm] = (int) ubuf(buf[m++]).i; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); for (i = first; i < last; i++) @@ -635,7 +633,7 @@ void AtomVec::unpack_comm_vel(int n, int first, double *buf) pdata = mcomm_vel.pdata[nn]; datatype = mcomm_vel.datatype[nn]; cols = mcomm_vel.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); for (i = first; i < last; i++) @@ -646,7 +644,7 @@ void AtomVec::unpack_comm_vel(int n, int first, double *buf) for (mm = 0; mm < cols; mm++) array[i][mm] = buf[m++]; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); for (i = first; i < last; i++) @@ -657,7 +655,7 @@ void AtomVec::unpack_comm_vel(int n, int first, double *buf) for (mm = 0; mm < cols; mm++) array[i][mm] = (int) ubuf(buf[m++]).i; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); for (i = first; i < last; i++) @@ -695,7 +693,7 @@ int AtomVec::pack_reverse(int n, int first, double *buf) pdata = mreverse.pdata[nn]; datatype = mreverse.datatype[nn]; cols = mreverse.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); for (i = first; i < last; i++) { @@ -708,7 +706,7 @@ int AtomVec::pack_reverse(int n, int first, double *buf) buf[m++] = array[i][mm]; } } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); for (i = first; i < last; i++) { @@ -721,7 +719,7 @@ int AtomVec::pack_reverse(int n, int first, double *buf) buf[m++] = ubuf(array[i][mm]).d; } } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); for (i = first; i < last; i++) { @@ -761,7 +759,7 @@ void AtomVec::unpack_reverse(int n, int *list, double *buf) pdata = mreverse.pdata[nn]; datatype = mreverse.datatype[nn]; cols = mreverse.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); for (i = 0; i < n; i++) { @@ -776,7 +774,7 @@ void AtomVec::unpack_reverse(int n, int *list, double *buf) array[j][mm] += buf[m++]; } } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); for (i = 0; i < n; i++) { @@ -791,7 +789,7 @@ void AtomVec::unpack_reverse(int n, int *list, double *buf) array[j][mm] += buf[m++]; } } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); for (i = 0; i < n; i++) { @@ -856,7 +854,7 @@ int AtomVec::pack_border(int n, int *list, double *buf, int pbc_flag, int *pbc) pdata = mborder.pdata[nn]; datatype = mborder.datatype[nn]; cols = mborder.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); for (i = 0; i < n; i++) { @@ -871,7 +869,7 @@ int AtomVec::pack_border(int n, int *list, double *buf, int pbc_flag, int *pbc) buf[m++] = array[j][mm]; } } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); for (i = 0; i < n; i++) { @@ -886,7 +884,7 @@ int AtomVec::pack_border(int n, int *list, double *buf, int pbc_flag, int *pbc) buf[m++] = ubuf(array[j][mm]).d; } } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); for (i = 0; i < n; i++) { @@ -990,7 +988,7 @@ int AtomVec::pack_border_vel(int n, int *list, double *buf, pdata = mborder_vel.pdata[nn]; datatype = mborder_vel.datatype[nn]; cols = mborder_vel.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); for (i = 0; i < n; i++) { @@ -1005,7 +1003,7 @@ int AtomVec::pack_border_vel(int n, int *list, double *buf, buf[m++] = array[j][mm]; } } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); for (i = 0; i < n; i++) { @@ -1020,7 +1018,7 @@ int AtomVec::pack_border_vel(int n, int *list, double *buf, buf[m++] = ubuf(array[j][mm]).d; } } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); for (i = 0; i < n; i++) { @@ -1072,7 +1070,7 @@ void AtomVec::unpack_border(int n, int first, double *buf) pdata = mborder.pdata[nn]; datatype = mborder.datatype[nn]; cols = mborder.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); for (i = first; i < last; i++) @@ -1083,7 +1081,7 @@ void AtomVec::unpack_border(int n, int first, double *buf) for (mm = 0; mm < cols; mm++) array[i][mm] = buf[m++]; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); for (i = first; i < last; i++) @@ -1094,7 +1092,7 @@ void AtomVec::unpack_border(int n, int first, double *buf) for (mm = 0; mm < cols; mm++) array[i][mm] = (int) ubuf(buf[m++]).i; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); for (i = first; i < last; i++) @@ -1144,7 +1142,7 @@ void AtomVec::unpack_border_vel(int n, int first, double *buf) pdata = mborder_vel.pdata[nn]; datatype = mborder_vel.datatype[nn]; cols = mborder_vel.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); for (i = first; i < last; i++) @@ -1155,7 +1153,7 @@ void AtomVec::unpack_border_vel(int n, int first, double *buf) for (mm = 0; mm < cols; mm++) array[i][mm] = buf[m++]; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); for (i = first; i < last; i++) @@ -1166,7 +1164,7 @@ void AtomVec::unpack_border_vel(int n, int first, double *buf) for (mm = 0; mm < cols; mm++) array[i][mm] = (int) ubuf(buf[m++]).i; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); for (i = first; i < last; i++) @@ -1216,7 +1214,7 @@ int AtomVec::pack_exchange(int i, double *buf) pdata = mexchange.pdata[nn]; datatype = mexchange.datatype[nn]; cols = mexchange.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); buf[m++] = vec[i]; @@ -1233,7 +1231,7 @@ int AtomVec::pack_exchange(int i, double *buf) for (mm = 0; mm < ncols; mm++) buf[m++] = array[i][mm]; } - } if (datatype == INT) { + } if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); buf[m++] = ubuf(vec[i]).d; @@ -1250,7 +1248,7 @@ int AtomVec::pack_exchange(int i, double *buf) for (mm = 0; mm < ncols; mm++) buf[m++] = ubuf(array[i][mm]).d; } - } if (datatype == BIGINT) { + } if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); buf[m++] = ubuf(vec[i]).d; @@ -1308,7 +1306,7 @@ int AtomVec::unpack_exchange(double *buf) pdata = mexchange.pdata[nn]; datatype = mexchange.datatype[nn]; cols = mexchange.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); vec[nlocal] = buf[m++]; @@ -1325,7 +1323,7 @@ int AtomVec::unpack_exchange(double *buf) for (mm = 0; mm < ncols; mm++) array[nlocal][mm] = buf[m++]; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); vec[nlocal] = ubuf(buf[m++]).i; @@ -1342,7 +1340,7 @@ int AtomVec::unpack_exchange(double *buf) for (mm = 0; mm < ncols; mm++) array[nlocal][mm] = (int) ubuf(buf[m++]).i; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); vec[nlocal] = (bigint) ubuf(buf[m++]).i; @@ -1450,7 +1448,7 @@ int AtomVec::pack_restart(int i, double *buf) pdata = mrestart.pdata[nn]; datatype = mrestart.datatype[nn]; cols = mrestart.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); buf[m++] = vec[i]; @@ -1467,7 +1465,7 @@ int AtomVec::pack_restart(int i, double *buf) for (mm = 0; mm < ncols; mm++) buf[m++] = array[i][mm]; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); buf[m++] = ubuf(vec[i]).d; @@ -1484,7 +1482,7 @@ int AtomVec::pack_restart(int i, double *buf) for (mm = 0; mm < ncols; mm++) buf[m++] = ubuf(array[i][mm]).d; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); buf[m++] = ubuf(vec[i]).d; @@ -1551,7 +1549,7 @@ int AtomVec::unpack_restart(double *buf) pdata = mrestart.pdata[nn]; datatype = mrestart.datatype[nn]; cols = mrestart.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); vec[nlocal] = buf[m++]; @@ -1568,7 +1566,7 @@ int AtomVec::unpack_restart(double *buf) for (mm = 0; mm < ncols; mm++) array[nlocal][mm] = buf[m++]; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); vec[nlocal] = ubuf(buf[m++]).i; @@ -1585,7 +1583,7 @@ int AtomVec::unpack_restart(double *buf) for (mm = 0; mm < ncols; mm++) array[nlocal][mm] = (int) ubuf(buf[m++]).i; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); vec[nlocal] = (bigint) ubuf(buf[m++]).i; @@ -1654,7 +1652,7 @@ void AtomVec::create_atom(int itype, double *coord) pdata = mcreate.pdata[n]; datatype = mcreate.datatype[n]; cols = mcreate.cols[n]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); vec[nlocal] = 0.0; @@ -1663,7 +1661,7 @@ void AtomVec::create_atom(int itype, double *coord) for (m = 0; m < cols; m++) array[nlocal][m] = 0.0; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); vec[nlocal] = 0; @@ -1672,7 +1670,7 @@ void AtomVec::create_atom(int itype, double *coord) for (m = 0; m < cols; m++) array[nlocal][m] = 0; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); vec[nlocal] = 0; @@ -1718,7 +1716,7 @@ void AtomVec::data_atom(double *coord, imageint imagetmp, char **values) pdata = mdata_atom.pdata[n]; datatype = mdata_atom.datatype[n]; cols = mdata_atom.cols[n]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); vec[nlocal] = utils::numeric(FLERR,values[ivalue++],true,lmp); @@ -1731,7 +1729,7 @@ void AtomVec::data_atom(double *coord, imageint imagetmp, char **values) for (m = 0; m < cols; m++) array[nlocal][m] = utils::numeric(FLERR,values[ivalue++],true,lmp); } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); vec[nlocal] = utils::inumeric(FLERR,values[ivalue++],true,lmp); @@ -1740,7 +1738,7 @@ void AtomVec::data_atom(double *coord, imageint imagetmp, char **values) for (m = 0; m < cols; m++) array[nlocal][m] = utils::inumeric(FLERR,values[ivalue++],true,lmp); } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); vec[nlocal] = utils::bnumeric(FLERR,values[ivalue++],true,lmp); @@ -1788,7 +1786,7 @@ void AtomVec::pack_data(double **buf) pdata = mdata_atom.pdata[n]; datatype = mdata_atom.datatype[n]; cols = mdata_atom.cols[n]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); buf[i][j++] = vec[i]; @@ -1797,7 +1795,7 @@ void AtomVec::pack_data(double **buf) for (m = 0; m < cols; m++) buf[i][j++] = array[i][m]; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); buf[i][j++] = ubuf(vec[i]).d; @@ -1806,7 +1804,7 @@ void AtomVec::pack_data(double **buf) for (m = 0; m < cols; m++) buf[i][j++] = ubuf(array[i][m]).d; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); buf[i][j++] = ubuf(vec[i]).d; @@ -1844,21 +1842,21 @@ void AtomVec::write_data(FILE *fp, int n, double **buf) for (nn = 1; nn < ndata_atom; nn++) { datatype = mdata_atom.datatype[nn]; cols = mdata_atom.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { fprintf(fp," %-1.16e",buf[i][j++]); } else { for (m = 0; m < cols; m++) fprintf(fp," %-1.16e",buf[i][j++]); } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { fprintf(fp," %d",(int) ubuf(buf[i][j++]).i); } else { for (m = 0; m < cols; m++) fprintf(fp," %d",(int) ubuf(buf[i][j++]).i); } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { fmt::print(fp," {}",(bigint) ubuf(buf[i][j++]).i); } else { @@ -1895,7 +1893,7 @@ void AtomVec::data_vel(int ilocal, char **values) pdata = mdata_vel.pdata[n]; datatype = mdata_vel.datatype[n]; cols = mdata_vel.cols[n]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); vec[ilocal] = utils::numeric(FLERR,values[ivalue++],true,lmp); @@ -1904,7 +1902,7 @@ void AtomVec::data_vel(int ilocal, char **values) for (m = 0; m < cols; m++) array[ilocal][m] = utils::numeric(FLERR,values[ivalue++],true,lmp); } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); vec[ilocal] = utils::inumeric(FLERR,values[ivalue++],true,lmp); @@ -1913,7 +1911,7 @@ void AtomVec::data_vel(int ilocal, char **values) for (m = 0; m < cols; m++) array[ilocal][m] = utils::inumeric(FLERR,values[ivalue++],true,lmp); } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); vec[ilocal] = utils::bnumeric(FLERR,values[ivalue++],true,lmp); @@ -1944,7 +1942,7 @@ void AtomVec::pack_vel(double **buf) pdata = mdata_vel.pdata[n]; datatype = mdata_vel.datatype[n]; cols = mdata_vel.cols[n]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { double *vec = *((double **) pdata); buf[i][j++] = vec[i]; @@ -1953,7 +1951,7 @@ void AtomVec::pack_vel(double **buf) for (m = 0; m < cols; m++) buf[i][j++] = array[i][m]; } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); buf[i][j++] = ubuf(vec[i]).d; @@ -1962,7 +1960,7 @@ void AtomVec::pack_vel(double **buf) for (m = 0; m < cols; m++) buf[i][j++] = ubuf(array[i][m]).d; } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bigint *vec = *((bigint **) pdata); buf[i][j++] = ubuf(vec[i]).d; @@ -1992,21 +1990,21 @@ void AtomVec::write_vel(FILE *fp, int n, double **buf) for (nn = 1; nn < ndata_vel; nn++) { datatype = mdata_vel.datatype[nn]; cols = mdata_vel.cols[nn]; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { fprintf(fp," %-1.16e",buf[i][j++]); } else { for (m = 0; m < cols; m++) fprintf(fp," %-1.16e",buf[i][j++]); } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { fprintf(fp," %d",(int) ubuf(buf[i][j++]).i); } else { for (m = 0; m < cols; m++) fprintf(fp," %d",(int) ubuf(buf[i][j++]).i); } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { fmt::print(fp," {}",(bigint) ubuf(buf[i][j++]).i); } else { @@ -2290,7 +2288,7 @@ bigint AtomVec::memory_usage() datatype = mgrow.datatype[i]; cols = mgrow.cols[i]; const int nthreads = threads[i] ? comm->nthreads : 1; - if (datatype == DOUBLE) { + if (datatype == Atom::DOUBLE) { if (cols == 0) { bytes += memory->usage(*((double **) pdata),nmax*nthreads); } else if (cols > 0) { @@ -2299,7 +2297,7 @@ bigint AtomVec::memory_usage() maxcols = *(mgrow.maxcols[i]); bytes += memory->usage(*((double ***) pdata),nmax*nthreads,maxcols); } - } else if (datatype == INT) { + } else if (datatype == Atom::INT) { if (cols == 0) { bytes += memory->usage(*((int **) pdata),nmax*nthreads); } else if (cols > 0) { @@ -2308,7 +2306,7 @@ bigint AtomVec::memory_usage() maxcols = *(mgrow.maxcols[i]); bytes += memory->usage(*((int ***) pdata),nmax*nthreads,maxcols); } - } else if (datatype == BIGINT) { + } else if (datatype == Atom::BIGINT) { if (cols == 0) { bytes += memory->usage(*((bigint **) pdata),nmax*nthreads); } else if (cols > 0) { From db9543ede2b558767f4d83ad4a1577bba8632586 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 14:36:41 -0400 Subject: [PATCH 05/19] use more common coding patterns, set maxexchange, use correct argument conversion functions --- src/USER-QTB/fix_qtb.cpp | 43 ++++++++++++++++++++++------------------ src/USER-QTB/fix_qtb.h | 27 +++++++++++++------------ 2 files changed, 38 insertions(+), 32 deletions(-) diff --git a/src/USER-QTB/fix_qtb.cpp b/src/USER-QTB/fix_qtb.cpp index fa15385859..95f876a529 100644 --- a/src/USER-QTB/fix_qtb.cpp +++ b/src/USER-QTB/fix_qtb.cpp @@ -62,27 +62,34 @@ FixQTB::FixQTB(LAMMPS *lmp, int narg, char **arg) : while (iarg < narg) { if (strcmp(arg[iarg],"temp") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qtb command"); - t_target = atof(arg[iarg+1]); if (t_target < 0.0) error->all(FLERR,"Fix qtb temp must be >= 0.0"); + t_target = force->numeric(FLERR,arg[iarg+1]); + if (t_target < 0.0) error->all(FLERR,"Fix qtb temp must be >= 0.0"); iarg += 2; } else if (strcmp(arg[iarg],"damp") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qtb command"); - t_period = atof(arg[iarg+1]); if (t_period <= 0.0) error->all(FLERR,"Fix qtb damp must be > 0.0"); fric_coef = 1/t_period; + t_period = force->numeric(FLERR,arg[iarg+1]); + if (t_period <= 0.0) error->all(FLERR,"Fix qtb damp must be > 0.0"); fric_coef = 1/t_period; iarg += 2; } else if (strcmp(arg[iarg],"seed") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qtb command"); - seed = atof(arg[iarg+1]); if (seed <= 0) error->all(FLERR,"Illegal fix qtb command"); + seed = force->inumeric(FLERR,arg[iarg+1]); + if (seed <= 0) error->all(FLERR,"Illegal fix qtb command"); iarg += 2; } else if (strcmp(arg[iarg],"f_max") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qtb command"); - f_max = atof(arg[iarg+1]); if (f_max <= 0) error->all(FLERR,"Illegal fix qtb command"); + f_max = force->numeric(FLERR,arg[iarg+1]); + if (f_max <= 0) error->all(FLERR,"Illegal fix qtb command"); iarg += 2; } else if (strcmp(arg[iarg],"N_f") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qtb command"); - N_f = atof(arg[iarg+1]); if (N_f <= 0) error->all(FLERR,"Illegal fix qtb command"); + N_f = force->inumeric(FLERR,arg[iarg+1]); + if (N_f <= 0) error->all(FLERR,"Illegal fix qtb command"); iarg += 2; } else error->all(FLERR,"Illegal fix qtb command"); } + maxexchange = 6*N_f+3; + // allocate qtb gfactor1 = NULL; gfactor3 = NULL; @@ -105,10 +112,6 @@ FixQTB::FixQTB(LAMMPS *lmp, int narg, char **arg) : // allocate random-arrays and fran grow_arrays(atom->nmax); atom->add_callback(0); - // memory->create(random_array_0,atom->nmax+300,2*N_f,"qtb:random_array_0"); - // memory->create(random_array_1,atom->nmax+300,2*N_f,"qtb:random_array_1"); - // memory->create(random_array_2,atom->nmax+300,2*N_f,"qtb:random_array_2"); - // memory->create(fran,atom->nmax+300,3,"qtb:fran"); // allocate omega_H and time_H memory->create(omega_H,2*N_f,"qtb:omega_H"); @@ -404,11 +407,12 @@ void FixQTB::copy_arrays(int i, int j, int /*delflag*/) ------------------------------------------------------------------------- */ int FixQTB::pack_exchange(int i, double *buf) { - for (int m = 0; m < 2*N_f; m++) buf[m] = random_array_0[i][m]; - for (int m = 0; m < 2*N_f; m++) buf[m+2*N_f] = random_array_1[i][m]; - for (int m = 0; m < 2*N_f; m++) buf[m+4*N_f] = random_array_2[i][m]; - for (int m = 0; m < 3; m++) buf[m+6*N_f] = fran[i][m]; - return 6*N_f+3; + int n = 0; + for (int m = 0; m < 2*N_f; m++) buf[n++] = random_array_0[i][m]; + for (int m = 0; m < 2*N_f; m++) buf[n++] = random_array_1[i][m]; + for (int m = 0; m < 2*N_f; m++) buf[n++] = random_array_2[i][m]; + for (int m = 0; m < 3; m++) buf[n++] = fran[i][m]; + return n; } /* ---------------------------------------------------------------------- @@ -416,9 +420,10 @@ int FixQTB::pack_exchange(int i, double *buf) ------------------------------------------------------------------------- */ int FixQTB::unpack_exchange(int nlocal, double *buf) { - for (int m = 0; m < 2*N_f; m++) random_array_0[nlocal][m] = buf[m]; - for (int m = 0; m < 2*N_f; m++) random_array_1[nlocal][m] = buf[m+2*N_f]; - for (int m = 0; m < 2*N_f; m++) random_array_2[nlocal][m] = buf[m+4*N_f]; - for (int m = 0; m < 3; m++) fran[nlocal][m] = buf[m+6*N_f]; - return 6*N_f+3; + int n = 0; + for (int m = 0; m < 2*N_f; m++) random_array_0[nlocal][m] = buf[n++]; + for (int m = 0; m < 2*N_f; m++) random_array_1[nlocal][m] = buf[n++]; + for (int m = 0; m < 2*N_f; m++) random_array_2[nlocal][m] = buf[n++]; + for (int m = 0; m < 3; m++) fran[nlocal][m] = buf[n++]; + return n; } diff --git a/src/USER-QTB/fix_qtb.h b/src/USER-QTB/fix_qtb.h index 8fd14f5c5b..86a91bc9f6 100644 --- a/src/USER-QTB/fix_qtb.h +++ b/src/USER-QTB/fix_qtb.h @@ -47,22 +47,23 @@ class FixQTB : public Fix { private: // qtb parameters - int counter_mu; // counter l and mu - double t_period, fric_coef; // friction coefficient - int seed; // seed for the random number generator - double f_max; // frequency cutoff - int N_f; // number of frequency grid - double t_target; // target qtb temperature + int counter_mu; // counter l and mu + double t_period, fric_coef; // friction coefficient + int seed; // seed for the random number generator + double f_max; // frequency cutoff + int N_f; // number of frequency grid + double t_target; // target qtb temperature char *id_temp; class Compute *temperature; - double h_timestep; // time step to update the random forces - int alpha; // number of time steps to update the random forces - class RanMars *random; // random number generator - double *gfactor1,*gfactor3; // factors of frictions and random forces - double *omega_H,*time_H; // H gives the desired power spectrum - double **random_array_0, **random_array_1, **random_array_2; // random number arrays give independence between atoms and directions + double h_timestep; // time step to update the random forces + int alpha; // number of time steps to update the random forces + class RanMars *random; // random number generator + double *gfactor1,*gfactor3; // factors of frictions and random forces + double *omega_H,*time_H; // H gives the desired power spectrum + // random number arrays give independence between atoms and directions + double **random_array_0, **random_array_1, **random_array_2; int nlevels_respa; - double **fran, fsum[3], fsumall[3]; // random forces and their sums + double **fran, fsum[3], fsumall[3]; // random forces and their sums }; } From 5c1236084c494403c3c94f8d418b6ec4625823c5 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 14:37:20 -0400 Subject: [PATCH 06/19] correctly compute the size of the maxexchange buffer. we put all fixes into one buffer. --- src/comm.cpp | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/src/comm.cpp b/src/comm.cpp index 3df3cc4d44..9102e2991d 100644 --- a/src/comm.cpp +++ b/src/comm.cpp @@ -258,10 +258,8 @@ void Comm::init_exchange() int onefix; maxexchange_fix = 0; - for (int i = 0; i < nfix; i++) { - onefix = fix[i]->maxexchange; - maxexchange_fix = MAX(maxexchange_fix,onefix); - } + for (int i = 0; i < nfix; i++) + maxexchange_fix += fix[i]->maxexchange; maxexchange = maxexchange_atom + maxexchange_fix; bufextra = maxexchange + BUFEXTRA; From fc1be785fcf0cb84c4ae5b6fb09bd605ca08ba63 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 15:08:10 -0400 Subject: [PATCH 07/19] port bugfix from fix qtb to fix qbmsst and some output simplification --- src/USER-QTB/fix_qbmsst.cpp | 163 +++++++++++++++++------------------- src/USER-QTB/fix_qtb.cpp | 9 +- 2 files changed, 84 insertions(+), 88 deletions(-) diff --git a/src/USER-QTB/fix_qbmsst.cpp b/src/USER-QTB/fix_qbmsst.cpp index 5ac611b13e..34350c9b05 100644 --- a/src/USER-QTB/fix_qbmsst.cpp +++ b/src/USER-QTB/fix_qbmsst.cpp @@ -21,6 +21,7 @@ #include #include #include +#include #include "atom.h" #include "force.h" #include "update.h" @@ -34,6 +35,7 @@ #include "kspace.h" #include "math_const.h" #include "utils.h" +#include "fmt/format.h" using namespace LAMMPS_NS; using namespace FixConst; @@ -98,58 +100,71 @@ FixQBMSST::FixQBMSST(LAMMPS *lmp, int narg, char **arg) : while (iarg < narg) { if (strcmp(arg[iarg],"q") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - qmass = atof(arg[iarg+1]); if (qmass < 0.0) error->all(FLERR,"Fix qbmsst qmass must be >= 0.0"); + qmass = force->numeric(FLERR,arg[iarg+1]); + if (qmass < 0.0) error->all(FLERR,"Fix qbmsst qmass must be >= 0.0"); iarg += 2; } else if (strcmp(arg[iarg],"mu") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - mu = atof(arg[iarg+1]); if (mu < 0.0) error->all(FLERR,"Fix qbmsst mu must be >= 0.0"); + mu = force->numeric(FLERR,arg[iarg+1]); + if (mu < 0.0) error->all(FLERR,"Fix qbmsst mu must be >= 0.0"); iarg += 2; } else if (strcmp(arg[iarg],"p0") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - p0 = atof(arg[iarg+1]); if (p0 < 0.0) error->all(FLERR,"Fix qbmsst p0 must be >= 0.0"); + p0 = force->numeric(FLERR,arg[iarg+1]); + if (p0 < 0.0) error->all(FLERR,"Fix qbmsst p0 must be >= 0.0"); p0_set = 1; iarg += 2; } else if (strcmp(arg[iarg],"v0") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - v0 = atof(arg[iarg+1]); if (v0 < 0.0) error->all(FLERR,"Fix qbmsst v0 must be >= 0.0"); + v0 = force->numeric(FLERR,arg[iarg+1]); + if (v0 < 0.0) error->all(FLERR,"Fix qbmsst v0 must be >= 0.0"); v0_set = 1; iarg += 2; } else if (strcmp(arg[iarg],"e0") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - e0 = atof(arg[iarg+1]); + e0 = force->numeric(FLERR,arg[iarg+1]); e0_set = 1; iarg += 2; } else if (strcmp(arg[iarg],"tscale") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - tscale = atof(arg[iarg+1]); if (tscale < 0.0 || tscale > 1.0) error->all(FLERR,"Fix qbmsst tscale must satisfy 0 <= tscale < 1"); + tscale = force->numeric(FLERR,arg[iarg+1]); + if (tscale < 0.0 || tscale > 1.0) error->all(FLERR,"Fix qbmsst tscale must satisfy 0 <= tscale < 1"); iarg += 2; } else if (strcmp(arg[iarg],"damp") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - t_period = atof(arg[iarg+1]); if (t_period <= 0.0) error->all(FLERR,"Fix qbmsst damp must be > 0.0"); fric_coef = 1/t_period; + t_period = force->numeric(FLERR,arg[iarg+1]); + if (t_period <= 0.0) error->all(FLERR,"Fix qbmsst damp must be > 0.0"); + fric_coef = 1/t_period; iarg += 2; } else if (strcmp(arg[iarg],"seed") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - seed = atof(arg[iarg+1]); if (seed <= 0) error->all(FLERR,"Fix qbmsst seed must be a positive integer"); + seed = force->inumeric(FLERR,arg[iarg+1]); + if (seed <= 0) error->all(FLERR,"Fix qbmsst seed must be a positive integer"); iarg += 2; } else if (strcmp(arg[iarg],"f_max") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - f_max = atof(arg[iarg+1]); if (f_max <= 0) error->all(FLERR,"Fix qbmsst f_max must be > 0.0"); + f_max = force->numeric(FLERR,arg[iarg+1]); + if (f_max <= 0) error->all(FLERR,"Fix qbmsst f_max must be > 0.0"); iarg += 2; } else if (strcmp(arg[iarg],"N_f") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - N_f = atof(arg[iarg+1]); if (N_f <= 0) error->all(FLERR,"Fix qbmsst N_f must be a positive integer"); + N_f = force->inumeric(FLERR,arg[iarg+1]); + if (N_f <= 0) error->all(FLERR,"Fix qbmsst N_f must be a positive integer"); iarg += 2; } else if (strcmp(arg[iarg],"eta") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - eta = atof(arg[iarg+1]); if (eta <= 0) error->all(FLERR,"Fix qbmsst eta must be >= 0.0"); + eta = force->numeric(FLERR,arg[iarg+1]); + if (eta <= 0) error->all(FLERR,"Fix qbmsst eta must be >= 0.0"); iarg += 2; } else if (strcmp(arg[iarg],"beta") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - beta = atof(arg[iarg+1]); if (beta <= 0) error->all(FLERR,"Fix qbmsst beta must be a positive integer"); + beta = force->inumeric(FLERR,arg[iarg+1]); + if (beta <= 0) error->all(FLERR,"Fix qbmsst beta must be a positive integer"); iarg += 2; } else if (strcmp(arg[iarg],"T_init") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix qbmsst command"); - t_init = atof(arg[iarg+1]); if (t_init <= 0) error->all(FLERR,"Fix qbmsst T_init must be >= 0.0"); + t_init = force->numeric(FLERR,arg[iarg+1]); + if (t_init <= 0) error->all(FLERR,"Fix qbmsst T_init must be >= 0.0"); iarg += 2; } else error->all(FLERR,"Illegal fix qbmsst command"); } @@ -157,48 +172,32 @@ FixQBMSST::FixQBMSST(LAMMPS *lmp, int narg, char **arg) : // check for periodicity in controlled dimensions if (domain->nonperiodic) error->all(FLERR,"Fix qbmsst requires a periodic box"); + maxexchange = 6*N_f+3; + // comments if (comm->me == 0) { - if (screen) { - fprintf(screen,"QBMSST parameters:\n"); - if (direction == 0) fprintf(screen," Shock in x direction\n"); - else if (direction == 1) fprintf(screen," Shock in y direction\n"); - else if (direction == 2) fprintf(screen," Shock in z direction\n"); + std::string msg = "QBMSST parameters:\n"; + + if (direction == 0) msg += " Shock in x direction\n"; + else if (direction == 1) msg += " Shock in y direction\n"; + else if (direction == 2) msg += " Shock in z direction\n"; - fprintf(screen," Cell mass-like parameter qmass (units of mass^2/length^4) = %12.5e\n", qmass); - fprintf(screen," Shock velocity = %12.5e\n", velocity); - fprintf(screen," Artificial viscosity (units of mass/length/time) = %12.5e\n", mu); + msg += fmt::format(" Cell mass-like parameter qmass " + "(units of mass^2/length^4) = {:12.5e}\n", qmass); + msg += fmt::format(" Shock velocity = {:12.5e}\n", velocity); + msg += fmt::format(" Artificial viscosity (units of " + "mass/length/time) = {:12.5e}\n", mu); - if (p0_set) - fprintf(screen," Initial pressure specified to be %12.5e\n", p0); - else fprintf(screen," Initial pressure calculated on first step\n"); - if (v0_set) - fprintf(screen," Initial volume specified to be %12.5e\n", v0); - else fprintf(screen," Initial volume calculated on first step\n"); - if (e0_set) - fprintf(screen," Initial energy specified to be %12.5e\n", e0); - else fprintf(screen," Initial energy calculated on first step\n"); - } - if (logfile) { - fprintf(logfile,"QBMSST parameters:\n"); - if (direction == 0) fprintf(logfile," Shock in x direction\n"); - else if (direction == 1) fprintf(logfile," Shock in y direction\n"); - else if (direction == 2) fprintf(logfile," Shock in z direction\n"); - - fprintf(logfile," Cell mass-like parameter qmass (units of mass^2/length^4) = %12.5e\n", qmass); - fprintf(logfile," Shock velocity = %12.5e\n", velocity); - fprintf(logfile," Artificial viscosity (units of mass/length/time) = %12.5e\n", mu); - - if (p0_set) - fprintf(logfile," Initial pressure specified to be %12.5e\n", p0); - else fprintf(logfile," Initial pressure calculated on first step\n"); - if (v0_set) - fprintf(logfile," Initial volume specified to be %12.5e\n", v0); - else fprintf(logfile," Initial volume calculated on first step\n"); - if (e0_set) - fprintf(logfile," Initial energy specified to be %12.5e\n", e0); - else fprintf(logfile," Initial energy calculated on first step\n"); - } + if (p0_set) + msg += fmt::format(" Initial pressure specified to be {:12.5e}\n", p0); + else msg += " Initial pressure calculated on first step\n"; + if (v0_set) + msg += fmt::format(" Initial volume specified to be {:12.5e}\n", v0); + else msg += " Initial volume calculated on first step\n"; + if (e0_set) + msg += fmt::format(" Initial energy specified to be {:12.5e}\n", e0); + else msg += " Initial energy calculated on first step\n"; + utils::logmesg(lmp,msg); } // create a new compute temp style @@ -363,14 +362,17 @@ void FixQBMSST::init() //set up the h time step for updating the random force \delta{}h=\frac{\pi}{\Omega_{max}} if (int(1.0/(2*f_max*dtv)) == 0) { - if (comm->me == 0) printf ("Warning: Either f_max is too high or the time step is too big, setting f_max to be 1/timestep!\n"); + if (comm->me == 0) error->warning(FLERR,"Either f_max is too high or the time step " + "is too big, setting f_max to be 1/timestep!\n"); h_timestep=dtv; alpha=1; } else { alpha=int(1.0/(2*f_max*dtv)); h_timestep=alpha*dtv; } - if (comm->me == 0) printf ("The effective maximum frequncy is now %f inverse time unit with alpha value as %d!\n", 0.5/h_timestep, alpha); + if (comm->me == 0 && screen) + fmt::print(screen,"The effective maximum frequency is now {} inverse time unit " + "with alpha value as {}!\n", 0.5/h_timestep, alpha); //gfactor is the random force \sqrt{\frac{2\gamma{}m_{i}}{\alpha*\delta{}t}}, \sqrt{12} makes the random array variance equal to unit. for (int i = 1; i <= atom->ntypes; i++) { @@ -432,20 +434,16 @@ void FixQBMSST::setup(int /*vflag*/) if ( v0_set == 0 ) { v0 = compute_vol(); v0_set = 1; - if (comm->me == 0) { - if ( screen ) fprintf(screen,"Fix QBMSST v0 = %12.5e\n", v0); - if ( logfile ) fprintf(logfile,"Fix QBMSST v0 = %12.5e\n", v0); - } + if (comm->me == 0) + utils::logmesg(lmp,fmt::format("Fix QBMSST v0 = {:12.5e}\n", v0)); } if ( p0_set == 0 ) { p0 = p_current[direction]; p0_set = 1; - if ( comm->me == 0 ) { - if ( screen ) fprintf(screen,"Fix QBMSST p0 = %12.5e\n", p0); - if ( logfile ) fprintf(logfile,"Fix QBMSST p0 = %12.5e\n", p0); - } + if ( comm->me == 0 ) + utils::logmesg(lmp,fmt::format("Fix QBMSST p0 = {:12.5e}\n", p0)); } if ( e0_set == 0 ) { @@ -453,11 +451,8 @@ void FixQBMSST::setup(int /*vflag*/) e0_set = 1; old_eavg = e0; - if ( comm->me == 0 ) { - if ( screen ) fprintf(screen,"Fix QBMSST e0 = to be %12.5e\n",e0); - if ( logfile ) fprintf(logfile,"Fix QBMSST e0 = to be %12.5e\n",e0); - } - + if ( comm->me == 0 ) + utils::logmesg(lmp,fmt::format("Fix QBMSST e0 = to be {:12.5e}\n",e0)); } temperature->compute_vector(); @@ -476,16 +471,10 @@ void FixQBMSST::setup(int /*vflag*/) omega[direction]=-1*sqrt(fac1); double fac2 = omega[direction]/v0; - if ( comm->me == 0 && tscale != 1.0) { - if ( screen ) - fprintf(screen,"Fix QBMSST initial strain rate of %12.5e established " - "by reducing temperature by factor of %12.5e\n", - fac2,tscale); - if ( logfile ) - fprintf(logfile,"Fix QBMSST initial strain rate of %12.5e established " - "by reducing temperature by factor of %12.5e\n", - fac2,tscale); - } + if ( comm->me == 0 && tscale != 1.0) + utils::logmesg(lmp,fmt::format("Fix QBMSST initial strain rate of {:12.5e} " + "established by reducing temperature by " + "factor of {:12.5e}\n",fac2,tscale)); for (int i = 0; i < atom->nlocal; i++) { if (mask[i] & groupbit) { for (int k = 0; k < 3; k++ ) { @@ -1170,11 +1159,12 @@ void FixQBMSST::copy_arrays(int i, int j, int /*delflag*/) ------------------------------------------------------------------------- */ int FixQBMSST::pack_exchange(int i, double *buf) { - for (int m = 0; m < 2*N_f; m++) buf[m] = random_array_0[i][m]; - for (int m = 0; m < 2*N_f; m++) buf[m+2*N_f] = random_array_1[i][m]; - for (int m = 0; m < 2*N_f; m++) buf[m+4*N_f] = random_array_2[i][m]; - for (int m = 0; m < 3; m++) buf[m+6*N_f] = fran[i][m]; - return 6*N_f+3; + int n = 0; + for (int m = 0; m < 2*N_f; m++) buf[n++] = random_array_0[i][m]; + for (int m = 0; m < 2*N_f; m++) buf[n++] = random_array_1[i][m]; + for (int m = 0; m < 2*N_f; m++) buf[n++] = random_array_2[i][m]; + for (int m = 0; m < 3; m++) buf[n++] = fran[i][m]; + return n; } /* ---------------------------------------------------------------------- @@ -1182,9 +1172,10 @@ int FixQBMSST::pack_exchange(int i, double *buf) ------------------------------------------------------------------------- */ int FixQBMSST::unpack_exchange(int nlocal, double *buf) { - for (int m = 0; m < 2*N_f; m++) random_array_0[nlocal][m] = buf[m]; - for (int m = 0; m < 2*N_f; m++) random_array_1[nlocal][m] = buf[m+2*N_f]; - for (int m = 0; m < 2*N_f; m++) random_array_2[nlocal][m] = buf[m+4*N_f]; - for (int m = 0; m < 3; m++) fran[nlocal][m] = buf[m+6*N_f]; - return 6*N_f+3; + int n = 0; + for (int m = 0; m < 2*N_f; m++) random_array_0[nlocal][m] = buf[n++]; + for (int m = 0; m < 2*N_f; m++) random_array_1[nlocal][m] = buf[n++]; + for (int m = 0; m < 2*N_f; m++) random_array_2[nlocal][m] = buf[n++]; + for (int m = 0; m < 3; m++) fran[nlocal][m] = buf[n++]; + return n; } diff --git a/src/USER-QTB/fix_qtb.cpp b/src/USER-QTB/fix_qtb.cpp index 95f876a529..646cf06ce2 100644 --- a/src/USER-QTB/fix_qtb.cpp +++ b/src/USER-QTB/fix_qtb.cpp @@ -32,6 +32,8 @@ #include "math_const.h" #include "memory.h" #include "error.h" +#include "utils.h" +#include "fmt/format.h" using namespace LAMMPS_NS; using namespace FixConst; @@ -163,14 +165,17 @@ void FixQTB::init() //set up the h time step for updating the random force \delta{}h=\frac{\pi}{\Omega_{max}} if (int(1.0/(2*f_max*dtv)) == 0) { - if (comm->me == 0) printf ("Warning: Either f_max is too high or the time step is too big, setting f_max to be 1/timestep!\n"); + if (comm->me == 0) error->warning(FLERR,"Either f_max is too high or the time step " + "is too big, setting f_max to be 1/timestep!\n"); h_timestep=dtv; alpha=1; } else { alpha=int(1.0/(2*f_max*dtv)); h_timestep=alpha*dtv; } - if (comm->me == 0) printf ("The effective maximum frequncy is now %f inverse time unit with alpha value as %d!\n", 0.5/h_timestep, alpha); + if (comm->me == 0 && screen) + fmt::print(screen,"The effective maximum frequency is now {} inverse time unit " + "with alpha value as {}!\n", 0.5/h_timestep, alpha); // set force prefactors if (!atom->rmass) { From 4da8ff3d219c6b63f8ad0b34ff97ca3f3553318e Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 15:47:46 -0400 Subject: [PATCH 08/19] avoid division by zero --- src/USER-QTB/fix_qbmsst.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/USER-QTB/fix_qbmsst.cpp b/src/USER-QTB/fix_qbmsst.cpp index 34350c9b05..bb396201cb 100644 --- a/src/USER-QTB/fix_qbmsst.cpp +++ b/src/USER-QTB/fix_qbmsst.cpp @@ -520,8 +520,8 @@ void FixQBMSST::initial_integrate(int /*vflag*/) // decide if the qtb temperature need to be updated or not if (counter_l == 0) { t_current -= dtv*fric_coef*eta*beta*(old_eavg-e0)/(3*ntotal*boltz); + if (t_current <= 0.0) break; old_eavg = 0;//clear old energy average - if (t_current < 0.0) t_current=0; // load omega_H with calculated spectrum at a specific temperature (corrected spectrum), omega_H is the Fourier transformation of time_H for (int k = 0; k < 2*N_f; k++) { @@ -529,7 +529,7 @@ void FixQBMSST::initial_integrate(int /*vflag*/) if(k == N_f) { omega_H[k]=sqrt(force->boltz * t_current); } else { - double energy_k= force->hplanck * fabs(f_k); + double energy_k= force->hplanck * fabs(f_k); omega_H[k]=sqrt( energy_k * (0.5+1.0/( exp(energy_k/(force->boltz * t_current)) - 1.0 )) ); omega_H[k]*=alpha*sin((k-N_f)*MY_PI/(2*alpha*N_f))/sin((k-N_f)*MY_PI/(2*N_f)); } From bbb6f408beede86b7bdfcfb993029849e765bfaa Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 15:56:18 -0400 Subject: [PATCH 09/19] fix syntax issue --- src/USER-QTB/fix_qbmsst.cpp | 43 +++++++++++++++++++------------------ 1 file changed, 22 insertions(+), 21 deletions(-) diff --git a/src/USER-QTB/fix_qbmsst.cpp b/src/USER-QTB/fix_qbmsst.cpp index bb396201cb..17eb4e9cea 100644 --- a/src/USER-QTB/fix_qbmsst.cpp +++ b/src/USER-QTB/fix_qbmsst.cpp @@ -520,30 +520,31 @@ void FixQBMSST::initial_integrate(int /*vflag*/) // decide if the qtb temperature need to be updated or not if (counter_l == 0) { t_current -= dtv*fric_coef*eta*beta*(old_eavg-e0)/(3*ntotal*boltz); - if (t_current <= 0.0) break; - old_eavg = 0;//clear old energy average + if (t_current > 0.0) { + old_eavg = 0;//clear old energy average - // load omega_H with calculated spectrum at a specific temperature (corrected spectrum), omega_H is the Fourier transformation of time_H - for (int k = 0; k < 2*N_f; k++) { - double f_k=(k-N_f)/(2*N_f*h_timestep); //\omega_k=\frac{2\pi}{\delta{}h}\frac{k}{2N_f} for k from -N_f to N_f-1 - if(k == N_f) { - omega_H[k]=sqrt(force->boltz * t_current); - } else { - double energy_k= force->hplanck * fabs(f_k); - omega_H[k]=sqrt( energy_k * (0.5+1.0/( exp(energy_k/(force->boltz * t_current)) - 1.0 )) ); - omega_H[k]*=alpha*sin((k-N_f)*MY_PI/(2*alpha*N_f))/sin((k-N_f)*MY_PI/(2*N_f)); - } - } - - // construct the signal filter H, filter has the unit of of sqrt(energy) \sqrt{2N_f}^{-1}H\left(t_n\right) - for (int n = 0; n < 2*N_f; n++) { - time_H[n] = 0; - double t_n=(n-N_f); + // load omega_H with calculated spectrum at a specific temperature (corrected spectrum), omega_H is the Fourier transformation of time_H for (int k = 0; k < 2*N_f; k++) { - double omega_k=(k-N_f)*MY_PI/N_f; - time_H[n] += omega_H[k]*(cos(omega_k*t_n)); + double f_k=(k-N_f)/(2*N_f*h_timestep); //\omega_k=\frac{2\pi}{\delta{}h}\frac{k}{2N_f} for k from -N_f to N_f-1 + if(k == N_f) { + omega_H[k]=sqrt(force->boltz * t_current); + } else { + double energy_k= force->hplanck * fabs(f_k); + omega_H[k]=sqrt( energy_k * (0.5+1.0/( exp(energy_k/(force->boltz * t_current)) - 1.0 )) ); + omega_H[k]*=alpha*sin((k-N_f)*MY_PI/(2*alpha*N_f))/sin((k-N_f)*MY_PI/(2*N_f)); + } + } + + // construct the signal filter H, filter has the unit of of sqrt(energy) \sqrt{2N_f}^{-1}H\left(t_n\right) + for (int n = 0; n < 2*N_f; n++) { + time_H[n] = 0; + double t_n=(n-N_f); + for (int k = 0; k < 2*N_f; k++) { + double omega_k=(k-N_f)*MY_PI/N_f; + time_H[n] += omega_H[k]*(cos(omega_k*t_n)); + } + time_H[n]/=(2.0*N_f); } - time_H[n]/=(2.0*N_f); } } From dd746f76956d635f39c40cfa8bc5d594353fad36 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 15:59:50 -0400 Subject: [PATCH 10/19] whitespace --- src/USER-QTB/fix_qbmsst.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/USER-QTB/fix_qbmsst.cpp b/src/USER-QTB/fix_qbmsst.cpp index 17eb4e9cea..d58c917308 100644 --- a/src/USER-QTB/fix_qbmsst.cpp +++ b/src/USER-QTB/fix_qbmsst.cpp @@ -173,11 +173,11 @@ FixQBMSST::FixQBMSST(LAMMPS *lmp, int narg, char **arg) : if (domain->nonperiodic) error->all(FLERR,"Fix qbmsst requires a periodic box"); maxexchange = 6*N_f+3; - + // comments if (comm->me == 0) { std::string msg = "QBMSST parameters:\n"; - + if (direction == 0) msg += " Shock in x direction\n"; else if (direction == 1) msg += " Shock in y direction\n"; else if (direction == 2) msg += " Shock in z direction\n"; From a9d1932032f548837ae98cf3892a6d5accc764fc Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 16:11:26 -0400 Subject: [PATCH 11/19] update singularity definition files to include kim-api packages --- tools/singularity/centos7.def | 2 +- tools/singularity/fedora32_mingw.def | 2 +- tools/singularity/ubuntu18.04.def | 5 ++++- tools/singularity/ubuntu20.04.def | 7 +++++-- 4 files changed, 11 insertions(+), 5 deletions(-) diff --git a/tools/singularity/centos7.def b/tools/singularity/centos7.def index 7803f3f386..6826ede8f6 100644 --- a/tools/singularity/centos7.def +++ b/tools/singularity/centos7.def @@ -11,7 +11,7 @@ From: centos:7 hdf5-devel python36-virtualenv python36-pip python-pip \ netcdf-devel netcdf-cxx-devel netcdf-mpich-devel netcdf-openmpi-devel \ python-virtualenv fftw-devel voro++-devel eigen3-devel gsl-devel openblas-devel enchant \ - blas-devel lapack-devel libyaml-devel + blas-devel lapack-devel libyaml-devel openkim-models kim-api-devel yum clean all # set custom prompt indicating the container name diff --git a/tools/singularity/fedora32_mingw.def b/tools/singularity/fedora32_mingw.def index 7c6c463243..860611e741 100644 --- a/tools/singularity/fedora32_mingw.def +++ b/tools/singularity/fedora32_mingw.def @@ -35,7 +35,7 @@ From: fedora:32 texlive-latex-bin texlive-lualatex-math texlive-fncychap texlive-tabulary \ texlive-framed texlive-wrapfig texlive-upquote texlive-capt-of \ texlive-needspace texlive-titlesec texlive-anysize texlive-dvipng \ - blas-devel lapack-devel libyaml-devel + blas-devel lapack-devel libyaml-devel openkim-models kim-api-devel dnf clean all # set custom prompt indicating the container name diff --git a/tools/singularity/ubuntu18.04.def b/tools/singularity/ubuntu18.04.def index 48b64d815b..2805a4a788 100644 --- a/tools/singularity/ubuntu18.04.def +++ b/tools/singularity/ubuntu18.04.def @@ -3,6 +3,7 @@ From: ubuntu:18.04 %post export DEBIAN_FRONTEND=noninteractive + add-apt-repository ppa:openkim/latest apt-get update apt-get upgrade --no-install-recommends -y apt-get install --no-install-recommends -y \ @@ -66,7 +67,9 @@ From: ubuntu:18.04 wget \ xxd \ valgrind \ - gdb + gdb \ + libkim-api-dev \ + openkim-models # clean cache rm -rf /var/lib/apt/lists/* diff --git a/tools/singularity/ubuntu20.04.def b/tools/singularity/ubuntu20.04.def index 0922bc5e86..7d5f29bdb9 100644 --- a/tools/singularity/ubuntu20.04.def +++ b/tools/singularity/ubuntu20.04.def @@ -3,6 +3,7 @@ From: ubuntu:20.04 %post export DEBIAN_FRONTEND=noninteractive + add-apt-repository ppa:openkim/latest apt-get update apt-get upgrade --no-install-recommends -y apt-get install --no-install-recommends -y \ @@ -62,7 +63,9 @@ From: ubuntu:20.04 wget \ xxd \ valgrind \ - gdb + gdb \ + libkim-api-dev \ + openkim-models # clean cache rm -rf /var/lib/apt/lists/* @@ -75,10 +78,10 @@ PS1="[ubuntu20.04:\u@\h] \W> " EOF chmod 755 $CUSTOM_PROMPT_ENV - %environment LC_ALL=C export LC_ALL + export PATH=/usr/lib/ccache:$PATH # restrict OpenMPI to shared memory comm by default OMPI_MCA_btl="tcp,self" # do not warn about unused components as this messes up testing From 864103f93eba9271ae1170f719ab51ab3685ed2b Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 21:35:06 -0400 Subject: [PATCH 12/19] fix cut-n-paste bug that crashed LAMMPS on reading molecular data files with -DLAMMPS_BIGBIG --- src/atom_vec.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/atom_vec.cpp b/src/atom_vec.cpp index 29f2526751..00c1abacb6 100644 --- a/src/atom_vec.cpp +++ b/src/atom_vec.cpp @@ -237,7 +237,7 @@ void AtomVec::grow(int n) memory->grow(*((bigint ***) pdata),nmax*nthreads,cols,"atom:barray"); else { maxcols = *(mgrow.maxcols[i]); - memory->grow(*((int ***) pdata),nmax*nthreads,maxcols,"atom:barray"); + memory->grow(*((bigint ***) pdata),nmax*nthreads,maxcols,"atom:barray"); } } } From e5f937388aa6b416b82975f4bf7a5e5fa85421a2 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 22:14:51 -0400 Subject: [PATCH 13/19] fix bugs with reading restart files when using -DLAMMPS_BIGBIG --- src/atom_vec.cpp | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/src/atom_vec.cpp b/src/atom_vec.cpp index 00c1abacb6..8cd7db4fd4 100644 --- a/src/atom_vec.cpp +++ b/src/atom_vec.cpp @@ -648,7 +648,7 @@ void AtomVec::unpack_comm_vel(int n, int first, double *buf) if (cols == 0) { int *vec = *((int **) pdata); for (i = first; i < last; i++) - vec[i] = ubuf(buf[m++]).i; + vec[i] = (int) ubuf(buf[m++]).i; } else { int **array = *((int ***) pdata); for (i = first; i < last; i++) @@ -779,14 +779,14 @@ void AtomVec::unpack_reverse(int n, int *list, double *buf) int *vec = *((int **) pdata); for (i = 0; i < n; i++) { j = list[i]; - vec[j] += buf[m++]; + vec[j] += (int) ubuf(buf[m++]).i; } } else { int **array = *((int ***) pdata); for (i = 0; i < n; i++) { j = list[i]; for (mm = 0; mm < cols; mm++) - array[j][mm] += buf[m++]; + array[j][mm] += (int) ubuf(buf[m++]).i; } } } else if (datatype == Atom::BIGINT) { @@ -794,14 +794,14 @@ void AtomVec::unpack_reverse(int n, int *list, double *buf) bigint *vec = *((bigint **) pdata); for (i = 0; i < n; i++) { j = list[i]; - vec[j] += buf[m++]; + vec[j] += (bigint) ubuf(buf[m++]).i; } } else { bigint **array = *((bigint ***) pdata); for (i = 0; i < n; i++) { j = list[i]; for (mm = 0; mm < cols; mm++) - array[j][mm] += buf[m++]; + array[j][mm] += (bigint) ubuf(buf[m++]).i; } } } @@ -1157,7 +1157,7 @@ void AtomVec::unpack_border_vel(int n, int first, double *buf) if (cols == 0) { int *vec = *((int **) pdata); for (i = first; i < last; i++) - vec[i] = ubuf(buf[m++]).i; + vec[i] = (int) ubuf(buf[m++]).i; } else { int **array = *((int ***) pdata); for (i = first; i < last; i++) @@ -1326,7 +1326,7 @@ int AtomVec::unpack_exchange(double *buf) } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); - vec[nlocal] = ubuf(buf[m++]).i; + vec[nlocal] = (int) ubuf(buf[m++]).i; } else if (cols > 0) { int **array = *((int ***) pdata); for (mm = 0; mm < cols; mm++) @@ -1569,7 +1569,7 @@ int AtomVec::unpack_restart(double *buf) } else if (datatype == Atom::INT) { if (cols == 0) { int *vec = *((int **) pdata); - vec[nlocal] = ubuf(buf[m++]).i; + vec[nlocal] = (int) ubuf(buf[m++]).i; } else if (cols > 0) { int **array = *((int ***) pdata); for (mm = 0; mm < cols; mm++) @@ -1592,7 +1592,7 @@ int AtomVec::unpack_restart(double *buf) for (mm = 0; mm < cols; mm++) array[nlocal][mm] = (bigint) ubuf(buf[m++]).i; } else { - int **array = *((int ***) pdata); + bigint **array = *((bigint ***) pdata); collength = mexchange.collength[nn]; plength = mexchange.plength[nn]; if (collength) ncols = (*((int ***) plength))[nlocal][collength-1]; From 448bccd1380f8007db3da893a712105a69dd1b52 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 7 Jun 2020 22:20:46 -0400 Subject: [PATCH 14/19] tweak test epsilons so they pass with -DLAMMPS_BIGBIG --- unittest/force-styles/tests/bond-gromos.yaml | 2 +- unittest/force-styles/tests/bond-nonlinear.yaml | 2 +- unittest/force-styles/tests/mol-pair-coul_cut.yaml | 2 +- unittest/force-styles/tests/mol-pair-coul_dsf.yaml | 2 +- unittest/force-styles/tests/mol-pair-coul_wolf.yaml | 2 +- .../tests/mol-pair-lj_charmm_coul_charmm_implicit.yaml | 2 +- 6 files changed, 6 insertions(+), 6 deletions(-) diff --git a/unittest/force-styles/tests/bond-gromos.yaml b/unittest/force-styles/tests/bond-gromos.yaml index a9f5dba6b4..112fcb9748 100644 --- a/unittest/force-styles/tests/bond-gromos.yaml +++ b/unittest/force-styles/tests/bond-gromos.yaml @@ -1,7 +1,7 @@ --- lammps_version: 5 May 2020 date_generated: Sat May 30 17:49:16 202 -epsilon: 2.5e-13 +epsilon: 5e-13 prerequisites: ! | atom full bond gromos diff --git a/unittest/force-styles/tests/bond-nonlinear.yaml b/unittest/force-styles/tests/bond-nonlinear.yaml index a8b6c20151..0d4935d0bb 100644 --- a/unittest/force-styles/tests/bond-nonlinear.yaml +++ b/unittest/force-styles/tests/bond-nonlinear.yaml @@ -1,7 +1,7 @@ --- lammps_version: 5 May 2020 date_generated: Sat May 30 17:49:16 202 -epsilon: 2.5e-13 +epsilon: 5e-13 prerequisites: ! | atom full bond nonlinear diff --git a/unittest/force-styles/tests/mol-pair-coul_cut.yaml b/unittest/force-styles/tests/mol-pair-coul_cut.yaml index 43d9dfe67e..bbae7a8456 100644 --- a/unittest/force-styles/tests/mol-pair-coul_cut.yaml +++ b/unittest/force-styles/tests/mol-pair-coul_cut.yaml @@ -1,7 +1,7 @@ --- lammps_version: 5 May 2020 date_generated: Sun May 31 08:39:49 202 -epsilon: 5e-14 +epsilon: 1e-13 prerequisites: ! | atom full pair coul/cut diff --git a/unittest/force-styles/tests/mol-pair-coul_dsf.yaml b/unittest/force-styles/tests/mol-pair-coul_dsf.yaml index e4f9983aae..98f2a64d3c 100644 --- a/unittest/force-styles/tests/mol-pair-coul_dsf.yaml +++ b/unittest/force-styles/tests/mol-pair-coul_dsf.yaml @@ -1,7 +1,7 @@ --- lammps_version: 5 May 2020 date_generated: Sun May 31 08:50:20 202 -epsilon: 1e-13 +epsilon: 2e-13 prerequisites: ! | atom full pair coul/dsf diff --git a/unittest/force-styles/tests/mol-pair-coul_wolf.yaml b/unittest/force-styles/tests/mol-pair-coul_wolf.yaml index 4da58e6ea6..bc8c243e44 100644 --- a/unittest/force-styles/tests/mol-pair-coul_wolf.yaml +++ b/unittest/force-styles/tests/mol-pair-coul_wolf.yaml @@ -1,7 +1,7 @@ --- lammps_version: 5 May 2020 date_generated: Sun May 31 09:23:33 202 -epsilon: 1e-13 +epsilon: 5e-13 prerequisites: ! | atom full pair coul/wolf diff --git a/unittest/force-styles/tests/mol-pair-lj_charmm_coul_charmm_implicit.yaml b/unittest/force-styles/tests/mol-pair-lj_charmm_coul_charmm_implicit.yaml index 1a1db6784e..2b20a41890 100644 --- a/unittest/force-styles/tests/mol-pair-lj_charmm_coul_charmm_implicit.yaml +++ b/unittest/force-styles/tests/mol-pair-lj_charmm_coul_charmm_implicit.yaml @@ -1,7 +1,7 @@ --- lammps_version: 5 May 2020 date_generated: Sun May 31 07:05:48 202 -epsilon: 5e-14 +epsilon: 2e-13 prerequisites: ! | atom full pair lj/charmm/coul/charmm/implicit From d106ed8981210b747743b62467ce302f490b8b2b Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Mon, 8 Jun 2020 01:49:16 -0400 Subject: [PATCH 15/19] add missing pair_style command to commands overview --- doc/src/Commands_all.rst | 1 + 1 file changed, 1 insertion(+) diff --git a/doc/src/Commands_all.rst b/doc/src/Commands_all.rst index 92c0ba107a..458238cca6 100644 --- a/doc/src/Commands_all.rst +++ b/doc/src/Commands_all.rst @@ -94,6 +94,7 @@ An alphabetic list of all general LAMMPS commands. * :doc:`package ` * :doc:`pair_coeff ` * :doc:`pair_modify ` + * :doc:`pair_style ` * :doc:`pair_write ` * :doc:`partition ` * :doc:`prd ` From 8b154cfbf623cb9b641f8c58b49f3e11e98610f6 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Mon, 8 Jun 2020 10:38:20 -0400 Subject: [PATCH 16/19] must not set *any* communicators in plumed lib when using MPI_STUBS --- src/USER-PLUMED/fix_plumed.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/USER-PLUMED/fix_plumed.cpp b/src/USER-PLUMED/fix_plumed.cpp index e1f0ea0bfe..3a52bbaaa1 100644 --- a/src/USER-PLUMED/fix_plumed.cpp +++ b/src/USER-PLUMED/fix_plumed.cpp @@ -82,6 +82,7 @@ FixPlumed::FixPlumed(LAMMPS *lmp, int narg, char **arg) : error->all(FLERR,"Incompatible API version for PLUMED in fix plumed. " "Only Plumed 2.4.x, 2.5.x, and 2.6.x are tested and supported."); +#if !defined(MPI_STUBS) // If the -partition option is activated then enable // inter-partition communication @@ -108,7 +109,6 @@ FixPlumed::FixPlumed(LAMMPS *lmp, int narg, char **arg) : // whereas if partitions are not defined then world is equal to // MPI_COMM_WORLD. -#if !defined(MPI_STUBS) // plumed does not know about LAMMPS using the MPI STUBS library and will // fail if this is called under these circumstances p->cmd("setMPIComm",&world); From 2668e3deb6eeb853be904a420b55954d76a469f3 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Mon, 8 Jun 2020 14:39:34 -0400 Subject: [PATCH 17/19] add manual installation of plumed library --- tools/singularity/centos7.def | 34 ++++++++++++++++++++++++++ tools/singularity/centos8.def | 36 +++++++++++++++++++++++++++- tools/singularity/fedora32_mingw.def | 34 ++++++++++++++++++++++++++ tools/singularity/ubuntu18.04.def | 15 ++++++++++++ tools/singularity/ubuntu20.04.def | 15 ++++++++++++ 5 files changed, 133 insertions(+), 1 deletion(-) diff --git a/tools/singularity/centos7.def b/tools/singularity/centos7.def index 6826ede8f6..5369c6977b 100644 --- a/tools/singularity/centos7.def +++ b/tools/singularity/centos7.def @@ -14,6 +14,40 @@ From: centos:7 blas-devel lapack-devel libyaml-devel openkim-models kim-api-devel yum clean all + # we need to reset any module variables + # inherited from the host. + unset __LMOD_REF_COUNT__LMFILES_ + unset __LMOD_REF_COUNT_PATH + unset __LMOD_REF_COUNT_LD_LIBRARY_PATH + unset __LMOD_REF_COUNT_MANPATH + unset __LMOD_REF_COUNT_MODULEPATH + unset __LMOD_REF_COUNT_LOADEDMODULES + unset _LMFILES_ + unset MODULEPATH + unset MODULESHOME + unset MODULEPATH_ROOT + unset LOADEDMODULES + unset LMOD_SYSTEM_DEFAULT_MODULES + + # load MPI by default + . /etc/profile + module load mpi + + # manually install Plumed + mkdir plumed + cd plumed + version=2.6.0 + curl -L -o plumed.tar.gz https://github.com/plumed/plumed2/releases/download/v${version}/plumed-src-${version}.tgz + tar -xzf plumed.tar.gz + cd plumed-${version} + ./configure --disable-doc --prefix=/usr + make + make install + # fix up installation for CentOS and Fedora + mv -v /usr/lib/pkgconfig/plumed* /usr/share/pkgconfig/ + cd ../../ + rm -rvf plumed + # set custom prompt indicating the container name CUSTOM_PROMPT_ENV=/.singularity.d/env/99-zz_custom_prompt.sh cat >$CUSTOM_PROMPT_ENV <$CUSTOM_PROMPT_ENV <$CUSTOM_PROMPT_ENV <$CUSTOM_PROMPT_ENV <$CUSTOM_PROMPT_ENV < Date: Mon, 8 Jun 2020 16:08:03 -0400 Subject: [PATCH 18/19] need one more step with ubuntu to allow enabling PPA repos --- tools/singularity/ubuntu18.04.def | 2 ++ tools/singularity/ubuntu20.04.def | 2 ++ 2 files changed, 4 insertions(+) diff --git a/tools/singularity/ubuntu18.04.def b/tools/singularity/ubuntu18.04.def index 93cc4e8939..bfc598b157 100644 --- a/tools/singularity/ubuntu18.04.def +++ b/tools/singularity/ubuntu18.04.def @@ -3,6 +3,8 @@ From: ubuntu:18.04 %post export DEBIAN_FRONTEND=noninteractive + apt-get update + apt-get install --no-install-recommends -y software-properties-common add-apt-repository ppa:openkim/latest apt-get update apt-get upgrade --no-install-recommends -y diff --git a/tools/singularity/ubuntu20.04.def b/tools/singularity/ubuntu20.04.def index cbc691b4b7..ed7a67bbe0 100644 --- a/tools/singularity/ubuntu20.04.def +++ b/tools/singularity/ubuntu20.04.def @@ -3,6 +3,8 @@ From: ubuntu:20.04 %post export DEBIAN_FRONTEND=noninteractive + apt-get update + apt-get install --no-install-recommends -y software-properties-common add-apt-repository ppa:openkim/latest apt-get update apt-get upgrade --no-install-recommends -y From 016629252d468ed66e02b2e0165abde07f237613 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Mon, 8 Jun 2020 19:18:01 -0400 Subject: [PATCH 19/19] a few more tweaks needed for centos8 and fedora32 --- tools/singularity/centos8.def | 2 +- tools/singularity/fedora32_mingw.def | 22 ++++------------------ 2 files changed, 5 insertions(+), 19 deletions(-) diff --git a/tools/singularity/centos8.def b/tools/singularity/centos8.def index f64df09d1f..d2eba155ec 100644 --- a/tools/singularity/centos8.def +++ b/tools/singularity/centos8.def @@ -49,7 +49,7 @@ From: centos:8 make make install # fix up installation for CentOS and Fedora - mv -v /usr/lib/pkgconfig/plumed* /usr/share/pkgconfig/ + mv -v /usr/lib64/pkgconfig/plumed* /usr/share/pkgconfig/ cd ../../ rm -rvf plumed diff --git a/tools/singularity/fedora32_mingw.def b/tools/singularity/fedora32_mingw.def index 79dca4dd23..6eb1944ec2 100644 --- a/tools/singularity/fedora32_mingw.def +++ b/tools/singularity/fedora32_mingw.def @@ -38,23 +38,9 @@ From: fedora:32 blas-devel lapack-devel libyaml-devel openkim-models kim-api-devel dnf clean all - # we need to reset any module variables - # inherited from the host. - unset __LMOD_REF_COUNT__LMFILES_ - unset __LMOD_REF_COUNT_PATH - unset __LMOD_REF_COUNT_LD_LIBRARY_PATH - unset __LMOD_REF_COUNT_MANPATH - unset __LMOD_REF_COUNT_MODULEPATH - unset __LMOD_REF_COUNT_LOADEDMODULES - unset _LMFILES_ - unset MODULEPATH - unset MODULESHOME - unset MODULEPATH_ROOT - unset LOADEDMODULES - unset LMOD_SYSTEM_DEFAULT_MODULES - - # load MPI by default - . /etc/profile + # enable Lmod and load MPI + source /usr/share/lmod/lmod/init/profile + module purge module load mpi # manually install Plumed @@ -68,7 +54,7 @@ From: fedora:32 make make install # fix up installation for CentOS and Fedora - mv -v /usr/lib/pkgconfig/plumed* /usr/share/pkgconfig/ + mv -v /usr/lib64/pkgconfig/plumed* /usr/share/pkgconfig/ cd ../../ rm -rvf plumed