ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Merge branch 'master' into more-fmtlib-use

Browse Source
This commit is contained in:
Axel Kohlmeyer
2020-06-29 13:04:22 -04:00
parent 6ace98b64e 5263e457cb
commit 3eda735d29
49 changed files with 470 additions and 514 deletions
Download Patch File Download Diff File Expand all files Collapse all files

4
doc/src/pair_eim.rst
View File

@ -48,7 +48,7 @@ and :math:`sigma_i` are calculated as
\sigma_i = & \sum_{j=i_1}^{i_N} q_j \cdot \psi_{ij} \left(r_{ij}\right) \\
E_i\left(q_i,\sigma_i\right) = & \frac{1}{2} \cdot q_i \cdot \sigma_i
where :math:`\eta_{ji} is a pairwise function describing electron flow from atom
where :math:`\eta_{ji}` is a pairwise function describing electron flow from atom
I to atom J, and :math:`\psi_{ij}` is another pairwise function. The multi-body
nature of the EIM potential is a result of the embedding energy term.
A complete list of all the pair functions used in EIM is summarized
@ -63,7 +63,7 @@ below
\right.\\
\eta_{ji} = & A_{\eta,ij}\left(\chi_j-\chi_i\right)f_c\left(r,r_{s,\eta,ij},r_{c,\eta,ij}\right) \\
\psi_{ij}\left(r\right) = & A_{\psi,ij}\exp\left(-\zeta_{ij}r\right)f_c\left(r,r_{s,\psi,ij},r_{c,\psi,ij}\right) \\
f_{c}\left(r,r_p,r_c\right) = & 0.510204 \mathrm{erfc}\left[\frac{1.64498\left(2r-r_p-r_c\right)}{r_c-r_p}\right] - 0.010204
f_{c}\left(r,r_p,r_c\right) = & 0.510204 \cdot \mathrm{erfc}\left[\frac{1.64498\left(2r-r_p-r_c\right)}{r_c-r_p}\right] - 0.010204
Here :math:`E_b, r_e, r_(c,\phi), \alpha, \beta, A_(\psi), \zeta, r_(s,\psi),
r_(c,\psi), A_(\eta), r_(s,\eta), r_(c,\eta), \chi,` and pair function type

18
examples/threebody/Si.sw
View File

@ -1,18 +0,0 @@
# DATE: 2007-06-11 CONTRIBUTOR: Aidan Thompson, athomps@sandia.gov CITATION: Stillinger and Weber, Phys Rev B, 31, 5262, (1985)
# Stillinger-Weber parameters for various elements and mixtures
# multiple entries can be added to this file, LAMMPS reads the ones it needs
# these entries are in LAMMPS "metal" units:
# epsilon = eV; sigma = Angstroms
# other quantities are unitless
# format of a single entry (one or more lines):
# element 1, element 2, element 3,
# epsilon, sigma, a, lambda, gamma, costheta0, A, B, p, q, tol
# Here are the original parameters in metal units, for Silicon from:
#
# Stillinger and Weber, Phys. Rev. B, v. 31, p. 5262, (1985)
#
Si Si Si 2.1683 2.0951 1.80 21.0 1.20 -0.333333333333
7.049556277 0.6022245584 4.0 0.0 0.0

1
examples/threebody/Si.sw Symbolic link
View File

@ -0,0 +1 @@
../../potentials/Si.sw

2
potentials/MOH.nb3b.harmonic
View File

@ -1,4 +1,4 @@
# DATE: 2016-07-29 UNITS: metal CONTRIBUTOR: Todd Zeitler, tzeitle@sandia.gov CITATION: none
# DATE: 2016-07-29 UNITS: real CONTRIBUTOR: Todd Zeitler, tzeitle@sandia.gov CITATION: none
# nb3b/harmonic (nonbonded 3-body harmonic) parameters for various elements
#
# multiple entries can be added to this file, LAMMPS reads the ones it needs

31
src/ASPHERE/fix_nph_asphere.cpp
View File

@ -13,6 +13,7 @@
#include "fix_nph_asphere.h"
#include <cstring>
#include <string>
#include "modify.h"
#include "error.h"
@ -35,35 +36,21 @@ FixNPHAsphere::FixNPHAsphere(LAMMPS *lmp, int narg, char **arg) :
// compute group = all since pressure is always global (group all)
// and thus its KE/temperature contribution should use group all
int n = strlen(id) + 6;
id_temp = new char[n];
strcpy(id_temp,id);
strcat(id_temp,"_temp");
std::string tcmd = id + std::string("_temp");
id_temp = new char[tcmd.size()+1];
strcpy(id_temp,tcmd.c_str());
char **newarg = new char*[3];
newarg[0] = id_temp;
newarg[1] = (char *) "all";
newarg[2] = (char *) "temp/asphere";
modify->add_compute(3,newarg);
delete [] newarg;
modify->add_compute(tcmd + " all temp/asphere");
tcomputeflag = 1;
// create a new compute pressure style
// id = fix-ID + press, compute group = all
// pass id_temp as 4th arg to pressure constructor
n = strlen(id) + 7;
id_press = new char[n];
strcpy(id_press,id);
strcat(id_press,"_press");
std::string pcmd = id + std::string("_press");
id_press = new char[pcmd.size()+1];
strcpy(id_press,pcmd.c_str());
newarg = new char*[4];
newarg[0] = id_press;
newarg[1] = (char *) "all";
newarg[2] = (char *) "pressure";
newarg[3] = id_temp;
modify->add_compute(4,newarg);
delete [] newarg;
modify->add_compute(pcmd + " all pressure " + std::string(id_temp));
pcomputeflag = 1;
}

31
src/ASPHERE/fix_npt_asphere.cpp
View File

@ -13,6 +13,7 @@
#include "fix_npt_asphere.h"
#include <cstring>
#include <string>
#include "modify.h"
#include "error.h"
@ -34,35 +35,21 @@ FixNPTAsphere::FixNPTAsphere(LAMMPS *lmp, int narg, char **arg) :
// compute group = all since pressure is always global (group all)
// and thus its KE/temperature contribution should use group all
int n = strlen(id) + 6;
id_temp = new char[n];
strcpy(id_temp,id);
strcat(id_temp,"_temp");
std::string tcmd = id + std::string("_temp");
id_temp = new char[tcmd.size()+1];
strcpy(id_temp,tcmd.c_str());
char **newarg = new char*[3];
newarg[0] = id_temp;
newarg[1] = (char *) "all";
newarg[2] = (char *) "temp/asphere";
modify->add_compute(3,newarg);
delete [] newarg;
modify->add_compute(tcmd + " all temp/asphere");
tcomputeflag = 1;
// create a new compute pressure style
// id = fix-ID + press, compute group = all
// pass id_temp as 4th arg to pressure constructor
n = strlen(id) + 7;
id_press = new char[n];
strcpy(id_press,id);
strcat(id_press,"_press");
std::string pcmd = id + std::string("_press");
id_press = new char[pcmd.size()+1];
strcpy(id_press,pcmd.c_str());
newarg = new char*[4];
newarg[0] = id_press;
newarg[1] = (char *) "all";
newarg[2] = (char *) "pressure";
newarg[3] = id_temp;
modify->add_compute(4,newarg);
delete [] newarg;
modify->add_compute(pcmd + " all pressure " + std::string(id_temp));
pcomputeflag = 1;
}

18
src/ASPHERE/fix_nvt_asphere.cpp
View File

@ -13,9 +13,11 @@
#include "fix_nvt_asphere.h"
#include <cstring>
#include <string>
#include "group.h"
#include "modify.h"
#include "error.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
using namespace FixConst;
@ -33,17 +35,11 @@ FixNVTAsphere::FixNVTAsphere(LAMMPS *lmp, int narg, char **arg) :
// create a new compute temp style
// id = fix-ID + temp
int n = strlen(id) + 6;
id_temp = new char[n];
strcpy(id_temp,id);
strcat(id_temp,"_temp");
std::string cmd = id + std::string("_temp");
id_temp = new char[cmd.size()+1];
strcpy(id_temp,cmd.c_str());
char **newarg = new char*[3];
newarg[0] = id_temp;
newarg[1] = group->names[igroup];
newarg[2] = (char *) "temp/asphere";
modify->add_compute(3,newarg);
delete [] newarg;
cmd += fmt::format(" {} temp/asphere",group->names[igroup]);
modify->add_compute(cmd);
tcomputeflag = 1;
}

31
src/BODY/fix_nph_body.cpp
View File

@ -17,6 +17,7 @@
#include "fix_nph_body.h"
#include <cstring>
#include <string>
#include "modify.h"
#include "error.h"
@ -38,35 +39,21 @@ FixNPHBody::FixNPHBody(LAMMPS *lmp, int narg, char **arg) :
// compute group = all since pressure is always global (group all)
// and thus its KE/temperature contribution should use group all
int n = strlen(id) + 6;
id_temp = new char[n];
strcpy(id_temp,id);
strcat(id_temp,"_temp");
std::string tcmd = id + std::string("_temp");
id_temp = new char[tcmd.size()+1];
strcpy(id_temp,tcmd.c_str());
char **newarg = new char*[3];
newarg[0] = id_temp;
newarg[1] = (char *) "all";
newarg[2] = (char *) "temp/body";
modify->add_compute(3,newarg);
delete [] newarg;
modify->add_compute(tcmd + " all temp/body");
tcomputeflag = 1;
// create a new compute pressure style
// id = fix-ID + press, compute group = all
// pass id_temp as 4th arg to pressure constructor
n = strlen(id) + 7;
id_press = new char[n];
strcpy(id_press,id);
strcat(id_press,"_press");
std::string pcmd = id + std::string("_press");
id_press = new char[pcmd.size()+1];
strcpy(id_press,pcmd.c_str());
newarg = new char*[4];
newarg[0] = id_press;
newarg[1] = (char *) "all";
newarg[2] = (char *) "pressure";
newarg[3] = id_temp;
modify->add_compute(4,newarg);
delete [] newarg;
modify->add_compute(pcmd + " all pressure " + std::string(id_temp));
pcomputeflag = 1;
}

31
src/BODY/fix_npt_body.cpp
View File

@ -17,6 +17,7 @@
#include "fix_npt_body.h"
#include <cstring>
#include <string>
#include "modify.h"
#include "error.h"
@ -38,35 +39,21 @@ FixNPTBody::FixNPTBody(LAMMPS *lmp, int narg, char **arg) :
// compute group = all since pressure is always global (group all)
// and thus its KE/temperature contribution should use group all
int n = strlen(id) + 6;
id_temp = new char[n];
strcpy(id_temp,id);
strcat(id_temp,"_temp");
std::string tcmd = id + std::string("_temp");
id_temp = new char[tcmd.size()+1];
strcpy(id_temp,tcmd.c_str());
char **newarg = new char*[3];
newarg[0] = id_temp;
newarg[1] = (char *) "all";
newarg[2] = (char *) "temp/body";
modify->add_compute(3,newarg);
delete [] newarg;
modify->add_compute(tcmd + " all temp/body");
tcomputeflag = 1;
// create a new compute pressure style
// id = fix-ID + press, compute group = all
// pass id_temp as 4th arg to pressure constructor
n = strlen(id) + 7;
id_press = new char[n];
strcpy(id_press,id);
strcat(id_press,"_press");
std::string pcmd = id + std::string("_press");
id_press = new char[pcmd.size()+1];
strcpy(id_press,pcmd.c_str());
newarg = new char*[4];
newarg[0] = id_press;
newarg[1] = (char *) "all";
newarg[2] = (char *) "pressure";
newarg[3] = id_temp;
modify->add_compute(4,newarg);
delete [] newarg;
modify->add_compute(pcmd + " all pressure " + std::string(id_temp));
pcomputeflag = 1;
}

16
src/BODY/fix_nvt_body.cpp
View File

@ -17,6 +17,7 @@
#include "fix_nvt_body.h"
#include <cstring>
#include <string>
#include "group.h"
#include "modify.h"
#include "error.h"
@ -37,17 +38,10 @@ FixNVTBody::FixNVTBody(LAMMPS *lmp, int narg, char **arg) :
// create a new compute temp style
// id = fix-ID + temp
int n = strlen(id) + 6;
id_temp = new char[n];
strcpy(id_temp,id);
strcat(id_temp,"_temp");
std::string tcmd = id + std::string("_temp");
id_temp = new char[tcmd.size()+1];
strcpy(id_temp,tcmd.c_str());
char **newarg = new char*[3];
newarg[0] = id_temp;
newarg[1] = group->names[igroup];
newarg[2] = (char *) "temp/body";
modify->add_compute(3,newarg);
delete [] newarg;
modify->add_compute(tcmd + " all temp/body");
tcomputeflag = 1;
}

20
src/CORESHELL/compute_temp_cs.cpp
View File

@ -19,6 +19,7 @@
#include "compute_temp_cs.h"
#include <mpi.h>
#include <cstring>
#include <string>
#include "atom.h"
#include "atom_vec.h"
#include "domain.h"
@ -30,6 +31,7 @@
#include "comm.h"
#include "memory.h"
#include "error.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
@ -66,21 +68,13 @@ ComputeTempCS::ComputeTempCS(LAMMPS *lmp, int narg, char **arg) :
// create a new fix STORE style
// id = compute-ID + COMPUTE_STORE, fix group = compute group
int n = strlen(id) + strlen("_COMPUTE_STORE") + 1;
id_fix = new char[n];
strcpy(id_fix,id);
strcat(id_fix,"_COMPUTE_STORE");
std::string fixcmd = id + std::string("_COMPUTE_STORE");
id_fix = new char[fixcmd.size()+1];
strcpy(id_fix,fixcmd.c_str());
char **newarg = new char*[6];
newarg[0] = id_fix;
newarg[1] = group->names[igroup];
newarg[2] = (char *) "STORE";
newarg[3] = (char *) "peratom";
newarg[4] = (char *) "0";
newarg[5] = (char *) "1";
modify->add_fix(6,newarg);
fixcmd += fmt::format(" {} STORE peratom 0 1", group->names[igroup]);
modify->add_fix(fixcmd);
fix = (FixStore *) modify->fix[modify->nfix-1];
delete [] newarg;
// set fix store values = 0 for now
// fill them in via setup() once Comm::borders() has been called

9
src/GRANULAR/pair_gran_hooke_history.cpp
View File

@ -19,6 +19,7 @@
#include <mpi.h>
#include <cmath>
#include <cstring>
#include <string>
#include "atom.h"
#include "force.h"
#include "update.h"
@ -65,13 +66,7 @@ PairGranHookeHistory::PairGranHookeHistory(LAMMPS *lmp) : Pair(lmp)
// this is so final order of Modify:fix will conform to input script
fix_history = NULL;
char **fixarg = new char*[3];
fixarg[0] = (char *) "NEIGH_HISTORY_HH_DUMMY";
fixarg[1] = (char *) "all";
fixarg[2] = (char *) "DUMMY";
modify->add_fix(3,fixarg,1);
delete [] fixarg;
modify->add_fix("NEIGH_HISTORY_HH_DUMMY all DUMMY");
fix_dummy = (FixDummy *) modify->fix[modify->nfix-1];
}

9
src/GRANULAR/pair_granular.cpp
View File

@ -21,6 +21,7 @@
#include <mpi.h>
#include <cmath>
#include <cstring>
#include <string>
#include "atom.h"
#include "force.h"
#include "update.h"
@ -96,13 +97,7 @@ PairGranular::PairGranular(LAMMPS *lmp) : Pair(lmp)
// this is so final order of Modify:fix will conform to input script
fix_history = NULL;
char **fixarg = new char*[3];
fixarg[0] = (char *) "NEIGH_HISTORY_GRANULAR_DUMMY";
fixarg[1] = (char *) "all";
fixarg[2] = (char *) "DUMMY";
modify->add_fix(3,fixarg,1);
delete [] fixarg;
modify->add_fix("NEIGH_HISTORY_GRANULAR_DUMMY all DUMMY");
fix_dummy = (FixDummy *) modify->fix[modify->nfix-1];
}

6
src/KIM/kim_init.cpp
View File

@ -291,11 +291,7 @@ void KimInit::do_init(char *model_name, char *user_units, char *model_units, KIM
int ifix = modify->find_fix("KIM_MODEL_STORE");
if (ifix >= 0) modify->delete_fix(ifix);
char *fixarg[3];
fixarg[0] = (char *)"KIM_MODEL_STORE";
fixarg[1] = (char *)"all";
fixarg[2] = (char *)"STORE/KIM";
modify->add_fix(3,fixarg);
modify->add_fix("KIM_MODEL_STORE all STORE/KIM");
ifix = modify->find_fix("KIM_MODEL_STORE");
FixStoreKIM *fix_store = (FixStoreKIM *) modify->fix[ifix];

15
src/KOKKOS/sna_kokkos_impl.h
View File

@ -552,10 +552,7 @@ void SNAKokkos<DeviceType>::compute_bi(const int& iatom_mod, const int& jjb, con
int idouble = 0;
for (int elem1 = 0; elem1 < nelements; elem1++) {
for (int elem2 = 0; elem2 < nelements; elem2++) {
const int jalloy = idouble;
for (int elem3 = 0; elem3 < nelements; elem3++) {
const int jjballoy = itriple;
double sumzu = 0.0;
double sumzu_temp = 0.0;
@ -566,7 +563,7 @@ void SNAKokkos<DeviceType>::compute_bi(const int& iatom_mod, const int& jjb, con
const int jjz_index = jjz+mb*(j+1)+ma;
if (2*mb == j) return; // I think we can remove this?
const auto utot = ulisttot_pack(iatom_mod, jju_index, elem3, iatom_div);
const auto zloc = zlist_pack(iatom_mod, jjz_index, jalloy, iatom_div);
const auto zloc = zlist_pack(iatom_mod, jjz_index, idouble, iatom_div);
sumzu_temp += utot.re * zloc.re + utot.im * zloc.im;
}
}
@ -582,7 +579,7 @@ void SNAKokkos<DeviceType>::compute_bi(const int& iatom_mod, const int& jjb, con
const int jjz_index = jjz+(mb-1)*(j+1)+(j+1)+ma;
const auto utot = ulisttot_pack(iatom_mod, jju_index, elem3, iatom_div);
const auto zloc = zlist_pack(iatom_mod, jjz_index, jalloy, iatom_div);
const auto zloc = zlist_pack(iatom_mod, jjz_index, idouble, iatom_div);
sumzu_temp += utot.re * zloc.re + utot.im * zloc.im;
}
@ -593,7 +590,7 @@ void SNAKokkos<DeviceType>::compute_bi(const int& iatom_mod, const int& jjb, con
const int jjz_index = jjz+(mb-1)*(j+1)+(j+1)+ma;
const auto utot = ulisttot_pack(iatom_mod, jju_index, elem3, iatom_div);
const auto zloc = zlist_pack(iatom_mod, jjz_index, jalloy, iatom_div);
const auto zloc = zlist_pack(iatom_mod, jjz_index, idouble, iatom_div);
sumzu += 0.5 * (utot.re * zloc.re + utot.im * zloc.im);
} // end if jeven
@ -607,7 +604,7 @@ void SNAKokkos<DeviceType>::compute_bi(const int& iatom_mod, const int& jjb, con
sumzu -= bzero[j];
}
}
blist_pack(iatom_mod, jjb, jjballoy, iatom_div) = sumzu;
blist_pack(iatom_mod, jjb, itriple, iatom_div) = sumzu;
//} // end loop over j
//} // end loop over j1, j2
itriple++;
@ -1078,7 +1075,7 @@ void SNAKokkos<DeviceType>::compute_bi_cpu(const typename Kokkos::TeamPolicy<Dev
int idouble = 0;
for (int elem1 = 0; elem1 < nelements; elem1++) {
for (int elem2 = 0; elem2 < nelements; elem2++) {
const auto jalloy = idouble;
auto jalloy = idouble; // must be non-const to work around gcc compiler bug
for (int elem3 = 0; elem3 < nelements; elem3++) {
Kokkos::parallel_for(Kokkos::TeamThreadRange(team,idxb_max),
[&] (const int& jjb) {
@ -1146,7 +1143,7 @@ void SNAKokkos<DeviceType>::compute_bi_cpu(const typename Kokkos::TeamPolicy<Dev
}
}
blist(jjb, jjballoy, iatom) = sumzu;
blist(jjb, itriple, iatom) = sumzu;
});
});
//} // end loop over j

43
src/MANYBODY/pair_eim.cpp
View File

@ -41,6 +41,7 @@ PairEIM::PairEIM(LAMMPS *lmp) : Pair(lmp)
restartinfo = 0;
one_coeff = 1;
manybody_flag = 1;
unit_convert_flag = utils::get_supported_conversions(utils::ENERGY);
setfl = NULL;
nmax = 0;
@ -477,7 +478,7 @@ void PairEIM::read_file(char *filename)
// read potential file
if( comm->me == 0) {
EIMPotentialFileReader reader(lmp, filename);
EIMPotentialFileReader reader(lmp, filename, unit_convert_flag);
reader.get_global(setfl);
@ -1050,14 +1051,18 @@ double PairEIM::memory_usage()
return bytes;
}
EIMPotentialFileReader::EIMPotentialFileReader(LAMMPS * lmp, const std::string & filename) :
EIMPotentialFileReader::EIMPotentialFileReader(LAMMPS *lmp,
const std::string &filename,
const int auto_convert) :
Pointers(lmp), filename(filename)
{
if (comm->me != 0) {
error->one(FLERR, "EIMPotentialFileReader should only be called by proc 0!");
}
FILE * fp = force->open_potential(filename.c_str());
int unit_convert = auto_convert;
FILE *fp = force->open_potential(filename.c_str(), &unit_convert);
conversion_factor = utils::get_conversion_factor(utils::ENERGY,unit_convert);
if (fp == NULL) {
error->one(FLERR, fmt::format("cannot open EIM potential file {}", filename));
@ -1186,7 +1191,7 @@ void EIMPotentialFileReader::parse(FILE * fp)
PairData data;
data.rcutphiA = values.next_double();
data.rcutphiR = values.next_double();
data.Eb = values.next_double();
data.Eb = values.next_double() * conversion_factor;
data.r0 = values.next_double();
data.alpha = values.next_double();
data.beta = values.next_double();
@ -1194,7 +1199,7 @@ void EIMPotentialFileReader::parse(FILE * fp)
data.Asigma = values.next_double();
data.rq = values.next_double();
data.rcutsigma = values.next_double();
data.Ac = values.next_double();
data.Ac = values.next_double() * conversion_factor;
data.zeta = values.next_double();
data.rs = values.next_double();
@ -1217,20 +1222,18 @@ void EIMPotentialFileReader::parse(FILE * fp)
}
}
void EIMPotentialFileReader::get_global(PairEIM::Setfl * setfl) {
void EIMPotentialFileReader::get_global(PairEIM::Setfl *setfl) {
setfl->division = division;
setfl->rbig = rbig;
setfl->rsmall = rsmall;
}
void EIMPotentialFileReader::get_element(PairEIM::Setfl * setfl, int i, const std::string & name) {
if (elements.find(name) == elements.end()) {
char str[128];
snprintf(str, 128, "Element %s not defined in EIM potential file", name.c_str());
error->one(FLERR, str);
}
void EIMPotentialFileReader::get_element(PairEIM::Setfl *setfl, int i,
const std::string &name) {
if (elements.find(name) == elements.end())
error->one(FLERR,"Element " + name + " not defined in EIM potential file");
ElementData & data = elements[name];
ElementData &data = elements[name];
setfl->ielement[i] = data.ielement;
setfl->mass[i] = data.mass;
setfl->negativity[i] = data.negativity;
@ -1240,16 +1243,16 @@ void EIMPotentialFileReader::get_element(PairEIM::Setfl * setfl, int i, const st
setfl->q0[i] = data.q0;
}
void EIMPotentialFileReader::get_pair(PairEIM::Setfl * setfl, int ij, const std::string & elemA, const std::string & elemB) {
void EIMPotentialFileReader::get_pair(PairEIM::Setfl *setfl, int ij,
const std::string &elemA,
const std::string &elemB) {
auto p = get_pair(elemA, elemB);
if (pairs.find(p) == pairs.end()) {
char str[128];
snprintf(str, 128, "Pair (%s, %s) not defined in EIM potential file", elemA.c_str(), elemB.c_str());
error->one(FLERR, str);
}
if (pairs.find(p) == pairs.end())
error->one(FLERR,"Element pair (" + elemA + ", " + elemB
+ ") is not defined in EIM potential file");
PairData & data = pairs[p];
PairData &data = pairs[p];
setfl->rcutphiA[ij] = data.rcutphiA;
setfl->rcutphiR[ij] = data.rcutphiR;
setfl->Eb[ij] = data.Eb;

18
src/MANYBODY/pair_eim.h
View File

@ -98,17 +98,21 @@ class EIMPotentialFileReader : protected Pointers {
std::string filename;
static const int MAXLINE = 1024;
char line[MAXLINE];
double conversion_factor;
void parse(FILE * fp);
char * next_line(FILE * fp);
std::pair<std::string, std::string> get_pair(const std::string & a, const std::string & b);
void parse(FILE *fp);
char *next_line(FILE *fp);
std::pair<std::string, std::string> get_pair(const std::string &a,
const std::string &b);
public:
EIMPotentialFileReader(class LAMMPS* lmp, const std::string & filename);
EIMPotentialFileReader(class LAMMPS* lmp, const std::string &filename,
const int auto_convert=0);
void get_global(PairEIM::Setfl * setfl);
void get_element(PairEIM::Setfl * setfl, int i, const std::string & name);
void get_pair(PairEIM::Setfl * setfl, int ij, const std::string & elemA, const std::string & elemB);
void get_global(PairEIM::Setfl *setfl);
void get_element(PairEIM::Setfl *setfl, int i, const std::string &name);
void get_pair(PairEIM::Setfl *setfl, int ij,
const std::string &elemA, const std::string &elemB);
private:
// potential parameters

13
src/MANYBODY/pair_gw.cpp
View File

@ -48,6 +48,7 @@ PairGW::PairGW(LAMMPS *lmp) : Pair(lmp)
restartinfo = 0;
one_coeff = 1;
manybody_flag = 1;
unit_convert_flag = utils::get_supported_conversions(utils::ENERGY);
nelements = 0;
elements = NULL;
@ -375,9 +376,14 @@ void PairGW::read_file(char *file)
// open file on proc 0
if (comm->me == 0) {
PotentialFileReader reader(lmp, file, "GW");
PotentialFileReader reader(lmp, file, "GW", unit_convert_flag);
char * line;
// transparently convert units for supported conversions
int unit_convert = reader.get_unit_convert();
double conversion_factor = utils::get_conversion_factor(utils::ENERGY,
unit_convert);
while((line = reader.next_line(NPARAMS_PER_LINE))) {
try {
ValueTokenizer values(line);
@ -427,6 +433,11 @@ void PairGW::read_file(char *file)
params[nparams].lam1 = values.next_double();
params[nparams].biga = values.next_double();
params[nparams].powermint = int(params[nparams].powerm);
if (unit_convert) {
params[nparams].biga *= conversion_factor;
params[nparams].bigb *= conversion_factor;
}
} catch (TokenizerException & e) {
error->one(FLERR, e.what());
}

12
src/MANYBODY/pair_gw_zbl.cpp
View File

@ -70,9 +70,14 @@ void PairGWZBL::read_file(char *file)
// open file on proc 0
if (comm->me == 0) {
PotentialFileReader reader(lmp, file, "GW/ZBL");
PotentialFileReader reader(lmp, file, "GW/ZBL", unit_convert_flag);
char * line;
// transparently convert units for supported conversions
int unit_convert = reader.get_unit_convert();
double conversion_factor = utils::get_conversion_factor(utils::ENERGY,
unit_convert);
while((line = reader.next_line(NPARAMS_PER_LINE))) {
try {
ValueTokenizer values(line);
@ -126,6 +131,11 @@ void PairGWZBL::read_file(char *file)
params[nparams].ZBLcut = values.next_double();
params[nparams].ZBLexpscale = values.next_double();
params[nparams].powermint = int(params[nparams].powerm);
if (unit_convert) {
params[nparams].biga *= conversion_factor;
params[nparams].bigb *= conversion_factor;
}
} catch (TokenizerException & e) {
error->one(FLERR, e.what());
}

10
src/MANYBODY/pair_nb3b_harmonic.cpp
View File

@ -47,6 +47,7 @@ PairNb3bHarmonic::PairNb3bHarmonic(LAMMPS *lmp) : Pair(lmp)
restartinfo = 0;
one_coeff = 1;
manybody_flag = 1;
unit_convert_flag = utils::get_supported_conversions(utils::ENERGY);
nelements = 0;
elements = NULL;
@ -291,9 +292,14 @@ void PairNb3bHarmonic::read_file(char *file)
// open file on proc 0
if (comm->me == 0) {
PotentialFileReader reader(lmp, file, "nb3b/harmonic");
PotentialFileReader reader(lmp, file, "nb3b/harmonic", unit_convert_flag);
char * line;
// transparently convert units for supported conversions
int unit_convert = reader.get_unit_convert();
double conversion_factor = utils::get_conversion_factor(utils::ENERGY,
unit_convert);
while((line = reader.next_line(NPARAMS_PER_LINE))) {
try {
ValueTokenizer values(line);
@ -331,6 +337,8 @@ void PairNb3bHarmonic::read_file(char *file)
params[nparams].k_theta = values.next_double();
params[nparams].theta0 = values.next_double();
params[nparams].cutoff = values.next_double();
if (unit_convert) params[nparams].k_theta *= conversion_factor;
} catch (TokenizerException & e) {
error->one(FLERR, e.what());
}

14
src/MC/fix_bond_swap.cpp
View File

@ -86,17 +86,11 @@ FixBondSwap::FixBondSwap(LAMMPS *lmp, int narg, char **arg) :
// create a new compute temp style
// id = fix-ID + temp, compute group = fix group
int n = strlen(id) + 6;
id_temp = new char[n];
strcpy(id_temp,id);
strcat(id_temp,"_temp");
std::string cmd = id + std::string("_temp");
id_temp = new char[cmd.size()+1];
strcpy(id_temp,cmd.c_str());
char **newarg = new char*[3];
newarg[0] = id_temp;
newarg[1] = (char *) "all";
newarg[2] = (char *) "temp";
modify->add_compute(3,newarg);
delete [] newarg;
modify->add_compute(cmd + " all temp");
tflag = 1;
// initialize atom list

10
src/PERI/pair_peri_eps.cpp
View File

@ -19,6 +19,7 @@
#include <mpi.h>
#include <cmath>
#include <cstring>
#include <string>
#include "atom.h"
#include "domain.h"
#include "lattice.h"
@ -515,14 +516,7 @@ void PairPeriEPS::init_style()
// if first init, create Fix needed for storing fixed neighbors
if (ifix_peri == -1) {
char **fixarg = new char*[3];
fixarg[0] = (char *) "PERI_NEIGH";
fixarg[1] = (char *) "all";
fixarg[2] = (char *) "PERI_NEIGH";
modify->add_fix(3,fixarg);
delete [] fixarg;
}
if (ifix_peri == -1) modify->add_fix("PERI_NEIGH all PERI_NEIGH");
// find associated PERI_NEIGH fix that must exist
// could have changed locations in fix list since created

10
src/PERI/pair_peri_lps.cpp
View File

@ -19,6 +19,7 @@
#include <mpi.h>
#include <cmath>
#include <cstring>
#include <string>
#include "atom.h"
#include "domain.h"
#include "lattice.h"
@ -439,14 +440,7 @@ void PairPeriLPS::init_style()
// if first init, create Fix needed for storing fixed neighbors
if (ifix_peri == -1) {
char **fixarg = new char*[3];
fixarg[0] = (char *) "PERI_NEIGH";
fixarg[1] = (char *) "all";
fixarg[2] = (char *) "PERI_NEIGH";
modify->add_fix(3,fixarg);
delete [] fixarg;
}
if (ifix_peri == -1) modify->add_fix("PERI_NEIGH all PERI_NEIGH");
// find associated PERI_NEIGH fix that must exist
// could have changed locations in fix list since created

10
src/PERI/pair_peri_pmb.cpp
View File

@ -20,6 +20,7 @@
#include <cmath>
#include <cfloat>
#include <cstring>
#include <string>
#include "atom.h"
#include "domain.h"
#include "lattice.h"
@ -369,14 +370,7 @@ void PairPeriPMB::init_style()
// if first init, create Fix needed for storing fixed neighbors
if (ifix_peri == -1) {
char **fixarg = new char*[3];
fixarg[0] = (char *) "PERI_NEIGH";
fixarg[1] = (char *) "all";
fixarg[2] = (char *) "PERI_NEIGH";
modify->add_fix(3,fixarg);
delete [] fixarg;
}
if (ifix_peri == -1) modify->add_fix("PERI_NEIGH all PERI_NEIGH");
// find associated PERI_NEIGH fix that must exist
// could have changed locations in fix list since created

10
src/PERI/pair_peri_ves.cpp
View File

@ -19,6 +19,7 @@
#include <mpi.h>
#include <cmath>
#include <cstring>
#include <string>
#include "atom.h"
#include "domain.h"
#include "lattice.h"
@ -495,14 +496,7 @@ void PairPeriVES::init_style()
// if first init, create Fix needed for storing fixed neighbors
if (ifix_peri == -1) {
char **fixarg = new char*[3];
fixarg[0] = (char *) "PERI_NEIGH";
fixarg[1] = (char *) "all";
fixarg[2] = (char *) "PERI_NEIGH";
modify->add_fix(3,fixarg);
delete [] fixarg;
}
if (ifix_peri == -1) modify->add_fix("PERI_NEIGH all PERI_NEIGH");
// find associated PERI_NEIGH fix that must exist
// could have changed locations in fix list since created

16
src/REPLICA/fix_neb.cpp
View File

@ -20,6 +20,7 @@
#include <mpi.h>
#include <cmath>
#include <cstring>
#include <string>
#include "universe.h"
#include "update.h"
#include "atom.h"
@ -148,17 +149,10 @@ FixNEB::FixNEB(LAMMPS *lmp, int narg, char **arg) :
// create a new compute pe style
// id = fix-ID + pe, compute group = all
int n = strlen(id) + 4;
id_pe = new char[n];
strcpy(id_pe,id);
strcat(id_pe,"_pe");
char **newarg = new char*[3];
newarg[0] = id_pe;
newarg[1] = (char *) "all";
newarg[2] = (char *) "pe";
modify->add_compute(3,newarg);
delete [] newarg;
std::string cmd = id + std::string("_pe");
id_pe = new char[cmd.size()+1];
strcpy(id_pe,cmd.c_str());
modify->add_compute(cmd + " all pe");
// initialize local storage

7
src/REPLICA/hyper.cpp
View File

@ -103,13 +103,8 @@ void Hyper::command(int narg, char **arg)
// create FixEventHyper class to store event and pre-quench states
char **args = new char*[3];
args[0] = (char *) "hyper_event";
args[1] = (char *) "all";
args[2] = (char *) "EVENT/HYPER";
modify->add_fix(3,args);
modify->add_fix("hyper_event all EVENT/HYPER");
fix_event = (FixEventHyper *) modify->fix[modify->nfix-1];
delete [] args;
// create Finish for timing output

14
src/REPLICA/prd.cpp
View File

@ -18,6 +18,7 @@
#include "prd.h"
#include <mpi.h>
#include <cstring>
#include <string>
#include "universe.h"
#include "update.h"
#include "atom.h"
@ -147,11 +148,7 @@ void PRD::command(int narg, char **arg)
// create ComputeTemp class to monitor temperature
char **args = new char*[3];
args[0] = (char *) "prd_temp";
args[1] = (char *) "all";
args[2] = (char *) "temp";
modify->add_compute(3,args);
modify->add_compute("prd_temp all temp");
temperature = modify->compute[modify->ncompute-1];
// create Velocity class for velocity creation in dephasing
@ -161,6 +158,7 @@ void PRD::command(int narg, char **arg)
velocity = new Velocity(lmp);
velocity->init_external("all");
char *args[2];
args[0] = (char *) "temp";
args[1] = (char *) "prd_temp";
velocity->options(2,args);
@ -173,10 +171,7 @@ void PRD::command(int narg, char **arg)
// create FixEventPRD class to store event and pre-quench states
args[0] = (char *) "prd_event";
args[1] = (char *) "all";
args[2] = (char *) "EVENT/PRD";
modify->add_fix(3,args);
modify->add_fix("prd_event all EVENT/PRD");
fix_event = (FixEventPRD *) modify->fix[modify->nfix-1];
// create Finish for timing output
@ -185,7 +180,6 @@ void PRD::command(int narg, char **arg)
// string clean-up
delete [] args;
delete [] loop_setting;
delete [] dist_setting;

7
src/USER-COLVARS/ndx_group.cpp
View File

@ -226,12 +226,7 @@ void Ndx2Group::create(char *name, bigint num, tagint *tags)
{
// wipe out all members if the group exists. gid==0 is group "all"
int gid = group->find(name);
if (gid > 0) {
char *cmd[2];
cmd[0] = name;
cmd[1] = (char *)"clear";
group->assign(2,cmd);
}
if (gid > 0) group->assign(std::string(name) + " clear");
// map from global to local
const int nlocal = atom->nlocal;

89
src/USER-REACTION/fix_bond_react.cpp
View File

@ -19,6 +19,7 @@ Contributing Author: Jacob Gissinger (jacob.gissinger@colorado.edu)
#include <mpi.h>
#include <cmath>
#include <cstring>
#include <string>
#include "update.h"
#include "modify.h"
#include "respa.h"
@ -41,6 +42,7 @@ Contributing Author: Jacob Gissinger (jacob.gissinger@colorado.edu)
#include "error.h"
#include "input.h"
#include "variable.h"
#include "fmt/format.h"
#include <algorithm>
@ -575,17 +577,11 @@ FixBondReact::~FixBondReact()
delete [] set;
if (group) {
char **newarg;
newarg = new char*[2];
newarg[0] = master_group;
newarg[1] = (char *) "delete";
group->assign(2,newarg);
group->assign(std::string(master_group) + " delete");
if (stabilization_flag == 1) {
newarg[0] = exclude_group;
group->assign(2,newarg);
group->assign(std::string(exclude_group) + " delete");
delete [] exclude_group;
}
delete [] newarg;
}
}
@ -608,59 +604,38 @@ it will have the name 'i_limit_tags' and will be intitialized to 0 (not in group
void FixBondReact::post_constructor()
{
int len;
// let's add the limit_tags per-atom property fix
int len = strlen("bond_react_props_internal") + 1;
id_fix2 = new char[len];
strcpy(id_fix2,"bond_react_props_internal");
std::string cmd = std::string("bond_react_props_internal");
id_fix2 = new char[cmd.size()+1];
strcpy(id_fix2,cmd.c_str());
int ifix = modify->find_fix(id_fix2);
if (ifix == -1) {
char **newarg = new char*[7];
newarg[0] = (char *) "bond_react_props_internal";
newarg[1] = (char *) "all"; // group ID is ignored
newarg[2] = (char *) "property/atom";
newarg[3] = (char *) "i_limit_tags";
newarg[4] = (char *) "i_react_tags";
newarg[5] = (char *) "ghost";
newarg[6] = (char *) "yes";
modify->add_fix(7,newarg);
delete [] newarg;
cmd += std::string(" all property/atom i_limit_tags i_react_tags ghost yes");
modify->add_fix(cmd);
}
// create master_group if not already existing
// NOTE: limit_tags and react_tags automaticaly intitialized to zero (unless read from restart)
group->find_or_create(master_group);
char **newarg;
newarg = new char*[5];
newarg[0] = master_group;
newarg[1] = (char *) "dynamic";
newarg[2] = (char *) "all";
newarg[3] = (char *) "property";
newarg[4] = (char *) "limit_tags";
group->assign(5,newarg);
delete [] newarg;
cmd = fmt::format("{} dynamic all property limit_tags",master_group);
group->assign(cmd);
if (stabilization_flag == 1) {
int igroup = group->find(exclude_group);
// create exclude_group if not already existing, or use as parent group if static
if (igroup == -1 || group->dynamic[igroup] == 0) {
// create stabilization per-atom property
len = strlen("bond_react_stabilization_internal") + 1;
id_fix3 = new char[len];
strcpy(id_fix3,"bond_react_stabilization_internal");
cmd = std::string("bond_react_stabilization_internal");
id_fix3 = new char[cmd.size()+1];
strcpy(id_fix3,cmd.c_str());
ifix = modify->find_fix(id_fix3);
if (ifix == -1) {
char **newarg = new char*[6];
newarg[0] = (char *) id_fix3;
newarg[1] = (char *) "all"; // group ID is ignored
newarg[2] = (char *) "property/atom";
newarg[3] = (char *) "i_statted_tags";
newarg[4] = (char *) "ghost";
newarg[5] = (char *) "yes";
modify->add_fix(6,newarg);
cmd += std::string(" all property/atom i_statted_tags ghost yes");
modify->add_fix(cmd);
fix3 = modify->fix[modify->nfix-1];
delete [] newarg;
}
len = strlen("statted_tags") + 1;
@ -680,16 +655,11 @@ void FixBondReact::post_constructor()
strcat(exclude_group,"_REACT");
group->find_or_create(exclude_group);
char **newarg;
newarg = new char*[5];
newarg[0] = exclude_group;
newarg[1] = (char *) "dynamic";
if (igroup == -1) newarg[2] = (char *) "all";
else newarg[2] = (char *) exclude_PARENT_group;
newarg[3] = (char *) "property";
newarg[4] = (char *) "statted_tags";
group->assign(5,newarg);
delete [] newarg;
if (igroup == -1)
cmd = fmt::format("{} dynamic all property statted_tags",exclude_group);
else
cmd = fmt::format("{} dynamic {} property statted_tags",exclude_group,exclude_PARENT_group);
group->assign(cmd);
delete [] exclude_PARENT_group;
// on to statted_tags (system-wide thermostat)
@ -737,21 +707,16 @@ void FixBondReact::post_constructor()
// let's create a new nve/limit fix to limit newly reacted atoms
len = strlen("bond_react_MASTER_nve_limit") + 1;
id_fix1 = new char[len];
strcpy(id_fix1,"bond_react_MASTER_nve_limit");
cmd = std::string("bond_react_MASTER_nve_limit");
id_fix1 = new char[cmd.size()+1];
strcpy(id_fix1,cmd.c_str());
ifix = modify->find_fix(id_fix1);
if (ifix == -1) {
char **newarg = new char*[4];
newarg[0] = id_fix1;
newarg[1] = master_group;
newarg[2] = (char *) "nve/limit";
newarg[3] = nve_limit_xmax;
modify->add_fix(4,newarg);
cmd += fmt::format(" {} nve/limit {}",master_group,nve_limit_xmax);
modify->add_fix(cmd);
fix1 = modify->fix[modify->nfix-1];
delete [] newarg;
}
}
}

2
src/balance.cpp
View File

@ -494,7 +494,7 @@ void Balance::weight_storage(char *prefix)
if (prefix) cmd = prefix;
cmd += "IMBALANCE_WEIGHTS";
int ifix = modify->find_fix(cmd.c_str());
int ifix = modify->find_fix(cmd);
if (ifix < 1) {
cmd += " all STORE peratom 0 1";
modify->add_fix(cmd);

20
src/compute_msd.cpp
View File

@ -14,6 +14,7 @@
#include "compute_msd.h"
#include <mpi.h>
#include <cstring>
#include <string>
#include "atom.h"
#include "update.h"
#include "group.h"
@ -21,6 +22,7 @@
#include "modify.h"
#include "fix_store.h"
#include "error.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
@ -63,21 +65,13 @@ ComputeMSD::ComputeMSD(LAMMPS *lmp, int narg, char **arg) :
// create a new fix STORE style for reference positions
// id = compute-ID + COMPUTE_STORE, fix group = compute group
int n = strlen(id) + strlen("_COMPUTE_STORE") + 1;
id_fix = new char[n];
strcpy(id_fix,id);
strcat(id_fix,"_COMPUTE_STORE");
std::string fixcmd = id + std::string("_COMPUTE_STORE");
id_fix = new char[fixcmd.size()+1];
strcpy(id_fix,fixcmd.c_str());
char **newarg = new char*[6];
newarg[0] = id_fix;
newarg[1] = group->names[igroup];
newarg[2] = (char *) "STORE";
newarg[3] = (char *) "peratom";
newarg[4] = (char *) "1";
newarg[5] = (char *) "3";
modify->add_fix(6,newarg);
fixcmd += fmt::format(" {} STORE peratom 1 3",group->names[igroup]);
modify->add_fix(fixcmd);
fix = (FixStore *) modify->fix[modify->nfix-1];
delete [] newarg;
// calculate xu,yu,zu for fix store array
// skip if reset from restart file

20
src/compute_msd_chunk.cpp
View File

@ -14,6 +14,7 @@
#include "compute_msd_chunk.h"
#include <mpi.h>
#include <cstring>
#include <string>
#include "atom.h"
#include "group.h"
#include "update.h"
@ -23,6 +24,7 @@
#include "fix_store.h"
#include "memory.h"
#include "error.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
@ -56,21 +58,13 @@ ComputeMSDChunk::ComputeMSDChunk(LAMMPS *lmp, int narg, char **arg) :
// potentially re-populate the fix array (and change it to correct size)
// otherwise size reset and init will be done in setup()
n = strlen(id) + strlen("_COMPUTE_STORE") + 1;
id_fix = new char[n];
strcpy(id_fix,id);
strcat(id_fix,"_COMPUTE_STORE");
std::string fixcmd = id + std::string("_COMPUTE_STORE");
id_fix = new char[fixcmd.size()+1];
strcpy(id_fix,fixcmd.c_str());
char **newarg = new char*[6];
newarg[0] = id_fix;
newarg[1] = group->names[igroup];
newarg[2] = (char *) "STORE";
newarg[3] = (char *) "global";
newarg[4] = (char *) "1";
newarg[5] = (char *) "1";
modify->add_fix(6,newarg);
fixcmd += fmt::format(" {} STORE global 1 1",group->names[igroup]);
modify->add_fix(fixcmd);
fix = (FixStore *) modify->fix[modify->nfix-1];
delete [] newarg;
}
/* ---------------------------------------------------------------------- */

21
src/compute_vacf.cpp
View File

@ -14,12 +14,14 @@
#include "compute_vacf.h"
#include <mpi.h>
#include <cstring>
#include <string>
#include "atom.h"
#include "update.h"
#include "group.h"
#include "modify.h"
#include "fix_store.h"
#include "error.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
@ -39,21 +41,12 @@ ComputeVACF::ComputeVACF(LAMMPS *lmp, int narg, char **arg) :
// create a new fix STORE style
// id = compute-ID + COMPUTE_STORE, fix group = compute group
int n = strlen(id) + strlen("_COMPUTE_STORE") + 1;
id_fix = new char[n];
strcpy(id_fix,id);
strcat(id_fix,"_COMPUTE_STORE");
char **newarg = new char*[6];
newarg[0] = id_fix;
newarg[1] = group->names[igroup];
newarg[2] = (char *) "STORE";
newarg[3] = (char *) "peratom";
newarg[4] = (char *) "1";
newarg[5] = (char *) "3";
modify->add_fix(6,newarg);
std::string fixcmd = id + std::string("_COMPUTE_STORE");
id_fix = new char[fixcmd.size()+1];
strcpy(id_fix,fixcmd.c_str());
fixcmd += fmt::format(" {} STORE peratom 1 3", group->names[igroup]);
modify->add_fix(fixcmd);
fix = (FixStore *) modify->fix[modify->nfix-1];
delete [] newarg;
// store current velocities in fix store array
// skip if reset from restart file

22
src/dump_custom.cpp
View File

@ -14,6 +14,7 @@
#include "dump_custom.h"
#include <cstdlib>
#include <cstring>
#include <string>
#include "atom.h"
#include "force.h"
#include "domain.h"
@ -28,6 +29,7 @@
#include "error.h"
#include "update.h"
#include "variable.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
@ -1981,22 +1983,12 @@ int DumpCustom::modify_param(int narg, char **arg)
"dump:thresh_fixID");
memory->grow(thresh_first,(nthreshlast+1),"dump:thresh_first");
int n = strlen(id) + strlen("_DUMP_STORE") + 8;
thresh_fixID[nthreshlast] = new char[n];
strcpy(thresh_fixID[nthreshlast],id);
sprintf(&thresh_fixID[nthreshlast][strlen(id)],"%d",nthreshlast);
strcat(thresh_fixID[nthreshlast],"_DUMP_STORE");
char **newarg = new char*[6];
newarg[0] = thresh_fixID[nthreshlast];
newarg[1] = group->names[igroup];
newarg[2] = (char *) "STORE";
newarg[3] = (char *) "peratom";
newarg[4] = (char *) "1";
newarg[5] = (char *) "1";
modify->add_fix(6,newarg);
std::string threshid = fmt::format("{}{}_DUMP_STORE",id,nthreshlast);
thresh_fixID[nthreshlast] = new char[threshid.size()+1];
strcpy(thresh_fixID[nthreshlast],threshid.c_str());
modify->add_fix(fmt::format("{} {} STORE peratom 1 1",threshid,
group->names[igroup]));
thresh_fix[nthreshlast] = (FixStore *) modify->fix[modify->nfix-1];
delete [] newarg;
thresh_last[nthreshlast] = nthreshlast;
thresh_first[nthreshlast] = 1;

32
src/fix_adapt.cpp
View File

@ -31,6 +31,7 @@
#include "memory.h"
#include "error.h"
#include "utils.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
using namespace FixConst;
@ -265,20 +266,12 @@ void FixAdapt::post_constructor()
id_fix_diam = NULL;
id_fix_chg = NULL;
char **newarg = new char*[6];
newarg[1] = group->names[igroup];
newarg[2] = (char *) "STORE";
newarg[3] = (char *) "peratom";
newarg[4] = (char *) "1";
newarg[5] = (char *) "1";
if (diamflag && atom->radius_flag) {
int n = strlen(id) + strlen("_FIX_STORE_DIAM") + 1;
id_fix_diam = new char[n];
strcpy(id_fix_diam,id);
strcat(id_fix_diam,"_FIX_STORE_DIAM");
newarg[0] = id_fix_diam;
modify->add_fix(6,newarg);
std::string fixcmd = id + std::string("_FIX_STORE_DIAM");
id_fix_diam = new char[fixcmd.size()+1];
strcpy(id_fix_diam,fixcmd.c_str());
fixcmd += fmt::format(" {} STORE peratom 1 1",group->names[igroup]);
modify->add_fix(fixcmd);
fix_diam = (FixStore *) modify->fix[modify->nfix-1];
if (fix_diam->restart_reset) fix_diam->restart_reset = 0;
@ -296,12 +289,11 @@ void FixAdapt::post_constructor()
}
if (chgflag && atom->q_flag) {
int n = strlen(id) + strlen("_FIX_STORE_CHG") + 1;
id_fix_chg = new char[n];
strcpy(id_fix_chg,id);
strcat(id_fix_chg,"_FIX_STORE_CHG");
newarg[0] = id_fix_chg;
modify->add_fix(6,newarg);
std::string fixcmd = id + std::string("_FIX_STORE_CHG");
id_fix_chg = new char[fixcmd.size()+1];
strcpy(id_fix_chg,fixcmd.c_str());
fixcmd += fmt::format(" {} STORE peratom 1 1",group->names[igroup]);
modify->add_fix(fixcmd);
fix_chg = (FixStore *) modify->fix[modify->nfix-1];
if (fix_chg->restart_reset) fix_chg->restart_reset = 0;
@ -317,8 +309,6 @@ void FixAdapt::post_constructor()
}
}
}
delete [] newarg;
}
/* ---------------------------------------------------------------------- */

6
src/fix_nph_sphere.cpp
View File

@ -39,8 +39,7 @@ FixNPHSphere::FixNPHSphere(LAMMPS *lmp, int narg, char **arg) :
id_temp = new char[tcmd.size()+1];
strcpy(id_temp,tcmd.c_str());
tcmd += " all temp/sphere";
modify->add_compute(tcmd);
modify->add_compute(tcmd + " all temp/sphere");
tcomputeflag = 1;
// create a new compute pressure style
@ -51,7 +50,6 @@ FixNPHSphere::FixNPHSphere(LAMMPS *lmp, int narg, char **arg) :
id_press = new char[pcmd.size()+1];
strcpy(id_press,pcmd.c_str());
pcmd += " all pressure " + std::string(id_temp);
modify->add_compute(pcmd);
modify->add_compute(pcmd + " all pressure " + std::string(id_temp));
pcomputeflag = 1;
}

6
src/fix_npt_sphere.cpp
View File

@ -39,8 +39,7 @@ FixNPTSphere::FixNPTSphere(LAMMPS *lmp, int narg, char **arg) :
id_temp = new char[tcmd.size()+1];
strcpy(id_temp,tcmd.c_str());
tcmd += " all temp/sphere";
modify->add_compute(tcmd);
modify->add_compute(tcmd + " all temp/sphere");
tcomputeflag = 1;
// create a new compute pressure style
@ -51,7 +50,6 @@ FixNPTSphere::FixNPTSphere(LAMMPS *lmp, int narg, char **arg) :
id_press = new char[pcmd.size()+1];
strcpy(id_press,pcmd.c_str());
pcmd += " all pressure " + std::string(id_temp);
modify->add_compute(pcmd);
modify->add_compute(pcmd + " all pressure " + std::string(id_temp));
pcomputeflag = 1;
}

1
src/fix_nvt_sphere.cpp
View File

@ -13,6 +13,7 @@
#include "fix_nvt_sphere.h"
#include <cstring>
#include <string>
#include "group.h"
#include "modify.h"
#include "error.h"

2
src/force.cpp
View File

@ -1033,7 +1033,7 @@ FILE *Force::open_potential(const char *name, int *auto_convert)
return nullptr;
}
} else {
if (units == unit_style) {
if (units.empty() || units == unit_style) {
*auto_convert = utils::NOCONVERT;
} else {
if ((units == "metal") && (unit_style == "real")

51
src/group.cpp
View File

@ -124,11 +124,8 @@ void Group::assign(int narg, char **arg)
int bits = inversemask[igroup];
for (i = 0; i < nlocal; i++) mask[i] &= bits;
if (dynamic[igroup]) {
std::string fixID = "GROUP_";
fixID += names[igroup];
modify->delete_fix(fixID.c_str());
}
if (dynamic[igroup])
modify->delete_fix(std::string("GROUP_") + names[igroup]);
delete [] names[igroup];
names[igroup] = NULL;
@ -491,24 +488,15 @@ void Group::assign(int narg, char **arg)
// if group is already dynamic, delete existing FixGroup
if (dynamic[igroup]) {
std::string fixID = "GROUP_";
fixID += names[igroup];
modify->delete_fix(fixID.c_str());
}
if (dynamic[igroup])
modify->delete_fix(std::string("GROUP_") + names[igroup]);
dynamic[igroup] = 1;
std::string fixID = "GROUP_";
fixID += names[igroup];
char **newarg = new char*[narg];
newarg[0] = (char *)fixID.c_str();
newarg[1] = arg[2];
newarg[2] = (char *) "GROUP";
for (int i = 3; i < narg; i++) newarg[i] = arg[i];
modify->add_fix(narg,newarg);
delete [] newarg;
std::string fixcmd = "GROUP_";
fixcmd += fmt::format("{} {} GROUP",names[igroup],arg[2]);
for (int i = 3; i < narg; i++) fixcmd += std::string(" ") + arg[i];
modify->add_fix(fixcmd);
// style = static
// remove dynamic FixGroup if necessary
@ -517,11 +505,8 @@ void Group::assign(int narg, char **arg)
if (narg != 2) error->all(FLERR,"Illegal group command");
if (dynamic[igroup]) {
std::string fixID = "GROUP_";
fixID += names[igroup];
modify->delete_fix(fixID.c_str());
}
if (dynamic[igroup])
modify->delete_fix(std::string("GROUP_") + names[igroup]);
dynamic[igroup] = 0;
@ -547,6 +532,22 @@ void Group::assign(int narg, char **arg)
}
}
/* ----------------------------------------------------------------------
convenience function to allow assigning to groups from a single string
------------------------------------------------------------------------- */
void Group::assign(const std::string &groupcmd)
{
std::vector<std::string> args = utils::split_words(groupcmd);
char **newarg = new char*[args.size()];
int i=0;
for (const auto &arg : args) {
newarg[i++] = (char *)arg.c_str();
}
assign(args.size(),newarg);
delete[] newarg;
}
/* ----------------------------------------------------------------------
add flagged atoms to a new or existing group
------------------------------------------------------------------------- */

2
src/group.h
View File

@ -16,6 +16,7 @@
#include "pointers.h"
#include <map>
#include <string>
namespace LAMMPS_NS {
@ -30,6 +31,7 @@ class Group : protected Pointers {
Group(class LAMMPS *);
~Group();
void assign(int, char **); // assign atoms to a group
void assign(const std::string &); // convenience function
void create(char *, int *); // add flagged atoms to a group
int find(const char *); // lookup name in list of groups
int find_or_create(const char *); // lookup name or create new group

11
src/min.cpp
View File

@ -113,18 +113,13 @@ Min::~Min()
void Min::init()
{
if (lmp->kokkos && !kokkosable)
error->all(FLERR,"Must use a Kokkos-enabled min style (e.g. min_style cg/kk) "
"with Kokkos minimize");
error->all(FLERR,"Must use a Kokkos-enabled min style "
"(e.g. min_style cg/kk) with Kokkos minimize");
// create fix needed for storing atom-based quantities
// will delete it at end of run
char **fixarg = new char*[3];
fixarg[0] = (char *) "MINIMIZE";
fixarg[1] = (char *) "all";
fixarg[2] = (char *) "MINIMIZE";
modify->add_fix(3,fixarg);
delete [] fixarg;
modify->add_fix("MINIMIZE all MINIMIZE");
fix_minimize = (FixMinimize *) modify->fix[modify->nfix-1];
// clear out extra global and per-atom dof

28
src/modify.cpp
View File

@ -1031,7 +1031,7 @@ void Modify::modify_fix(int narg, char **arg)
Atom class must update indices in its list of callbacks to fixes
------------------------------------------------------------------------- */
void Modify::delete_fix(const char *id)
void Modify::delete_fix(const std::string &id)
{
int ifix = find_fix(id);
if (ifix < 0) error->all(FLERR,"Could not find fix ID to delete");
@ -1055,14 +1055,12 @@ void Modify::delete_fix(int ifix)
return index of fix or -1 if not found
------------------------------------------------------------------------- */
int Modify::find_fix(const char *id)
int Modify::find_fix(const std::string &id)
{
if (id == NULL) return -1;
int ifix;
for (ifix = 0; ifix < nfix; ifix++)
if (strcmp(id,fix[ifix]->id) == 0) break;
if (ifix == nfix) return -1;
return ifix;
if (id.empty()) return -1;
for (int ifix = 0; ifix < nfix; ifix++)
if (id == fix[ifix]->id) return ifix;
return -1;
}
/* ----------------------------------------------------------------------
@ -1295,7 +1293,7 @@ void Modify::modify_compute(int narg, char **arg)
delete a Compute from list of Computes
------------------------------------------------------------------------- */
void Modify::delete_compute(const char *id)
void Modify::delete_compute(const std::string &id)
{
int icompute = find_compute(id);
if (icompute < 0) error->all(FLERR,"Could not find compute ID to delete");
@ -1312,14 +1310,12 @@ void Modify::delete_compute(const char *id)
return index of compute or -1 if not found
------------------------------------------------------------------------- */
int Modify::find_compute(const char *id)
int Modify::find_compute(const std::string &id)
{
if(id==NULL) return -1;
int icompute;
for (icompute = 0; icompute < ncompute; icompute++)
if (strcmp(id,compute[icompute]->id) == 0) break;
if (icompute == ncompute) return -1;
return icompute;
if(id.empty()) return -1;
for (int icompute = 0; icompute < ncompute; icompute++)
if (id == compute[icompute]->id) return icompute;
return -1;
}
/* ----------------------------------------------------------------------

8
src/modify.h
View File

@ -99,9 +99,9 @@ class Modify : protected Pointers {
void add_fix(const std::string &, int trysuffix=1);
void replace_fix(const char *, int, char **, int trysuffix=1);
void modify_fix(int, char **);
void delete_fix(const char *);
void delete_fix(const std::string &);
void delete_fix(int);
int find_fix(const char *);
int find_fix(const std::string &);
int find_fix_by_style(const char *);
int check_package(const char *);
int check_rigid_group_overlap(int);
@ -111,8 +111,8 @@ class Modify : protected Pointers {
void add_compute(int, char **, int trysuffix=1);
void add_compute(const std::string &, int trysuffix=1);
void modify_compute(int, char **);
void delete_compute(const char *);
int find_compute(const char *);
void delete_compute(const std::string &);
int find_compute(const std::string &);
void clearstep_compute();
void addstep_compute(bigint);

18
src/pair_table.cpp
View File

@ -1033,11 +1033,15 @@ void *PairTable::extract(const char *str, int &dim)
if (strcmp(str,"cut_coul") != 0) return NULL;
if (ntables == 0) error->all(FLERR,"All pair coeffs are not set");
double cut_coul = tables[0].cut;
for (int m = 1; m < ntables; m++)
if (tables[m].cut != cut_coul)
error->all(FLERR,
"Pair table cutoffs must all be equal to use with KSpace");
dim = 0;
return &tables[0].cut;
// only check for cutoff consistency if claiming to be KSpace compatible
if (ewaldflag || pppmflag || msmflag || dispersionflag || tip4pflag) {
double cut_coul = tables[0].cut;
for (int m = 1; m < ntables; m++)
if (tables[m].cut != cut_coul)
error->all(FLERR,
"Pair table cutoffs must all be equal to use with KSpace");
dim = 0;
return &tables[0].cut;
} else return NULL;
}

16
src/read_restart.cpp
View File

@ -423,19 +423,9 @@ void ReadRestart::command(int narg, char **arg)
// create a temporary fix to hold and migrate extra atom info
// necessary b/c irregular will migrate atoms
if (nextra) {
char cextra[8],fixextra[8];
sprintf(cextra,"%d",nextra);
sprintf(fixextra,"%d",modify->nfix_restart_peratom);
char **newarg = new char*[5];
newarg[0] = (char *) "_read_restart";
newarg[1] = (char *) "all";
newarg[2] = (char *) "READ_RESTART";
newarg[3] = cextra;
newarg[4] = fixextra;
modify->add_fix(5,newarg);
delete [] newarg;
}
if (nextra)
modify->add_fix(fmt::format("_read_restart all READ_RESTART {} {}",
nextra,modify->nfix_restart_peratom));
// move atoms to new processors via irregular()
// turn sorting on in migrate_atoms() to avoid non-reproducible restarts

19
src/respa.cpp
View File

@ -299,22 +299,11 @@ void Respa::init()
// create fix needed for storing atom-based respa level forces
// will delete it at end of run
char **fixarg = new char*[5];
fixarg[0] = (char *) "RESPA";
fixarg[1] = (char *) "all";
fixarg[2] = (char *) "RESPA";
fixarg[3] = new char[8];
sprintf(fixarg[3],"%d",nlevels);
// if supported, we also store torques on a per-level basis
if (atom->torque_flag) {
fixarg[4] = (char *) "torque";
modify->add_fix(5,fixarg);
} else {
modify->add_fix(4,fixarg);
}
delete [] fixarg[3];
delete [] fixarg;
std::string cmd = fmt::format("RESPA all RESPA {}",nlevels);
if (atom->torque_flag) modify->add_fix(cmd + " torque");
else modify->add_fix(cmd);
fix_respa = (FixRespa *) modify->fix[modify->nfix-1];
// insure respa inner/middle/outer is using Pair class that supports it

2
src/thermo.cpp
View File

@ -483,7 +483,7 @@ void Thermo::modify_params(int narg, char **arg)
icompute = modify->find_compute(id_compute[index_press_vector]);
if (icompute < 0) error->all(FLERR,
"Pressure ID for thermo does not exist");
} else icompute = modify->find_compute((char *) "thermo_press");
} else icompute = modify->find_compute("thermo_press");
modify->compute[icompute]->reset_extra_compute_fix(arg[iarg+1]);

172
unittest/formats/test_pair_unit_convert.cpp
View File

@ -353,6 +353,178 @@ TEST_F(PairUnitConvertTest, eam_cd)
EXPECT_NEAR(ev_convert * fold[i][j], f[i][j], fabs(f[i][j] * rel_error));
}
TEST_F(PairUnitConvertTest, eim)
{
// check if the prerequisite pair style is available
if (!info->has_style("pair", "eim")) GTEST_SKIP();
if (!verbose) ::testing::internal::CaptureStdout();
lmp->input->one("units metal");
lmp->input->one("read_data test_pair_unit_convert.data");
lmp->input->one("pair_style eim");
lmp->input->one("pair_coeff * * Na Cl ffield.eim Na Cl");
lmp->input->one("run 0 post no");
if (!verbose) ::testing::internal::GetCapturedStdout();
// copy pressure, energy, and force from first step
double pold;
lmp->output->thermo->evaluate_keyword("press", &pold);
double eold = lmp->force->pair->eng_vdwl + lmp->force->pair->eng_coul;
double **f = lmp->atom->f;
for (int i = 0; i < 4; ++i)
for (int j = 0; j < 3; ++j)
fold[i][j] = f[i][j];
if (!verbose) ::testing::internal::CaptureStdout();
lmp->input->one("clear");
lmp->input->one("units real");
lmp->input->one("read_data test_pair_unit_convert.data");
lmp->input->one("pair_style eim");
lmp->input->one("pair_coeff * * Na Cl ffield.eim Na Cl");
lmp->input->one("run 0 post no");
if (!verbose) ::testing::internal::GetCapturedStdout();
double pnew;
lmp->output->thermo->evaluate_keyword("press", &pnew);
EXPECT_NEAR(pold, p_convert * pnew, fabs(pnew * rel_error));
double enew = lmp->force->pair->eng_vdwl + lmp->force->pair->eng_coul;
EXPECT_NEAR(ev_convert * eold, enew, fabs(enew * rel_error));
f = lmp->atom->f;
for (int i = 0; i < 4; ++i)
for (int j = 0; j < 3; ++j)
EXPECT_NEAR(ev_convert * fold[i][j], f[i][j], fabs(f[i][j] * rel_error));
}
TEST_F(PairUnitConvertTest, gw)
{
// check if the prerequisite pair style is available
if (!info->has_style("pair", "gw")) GTEST_SKIP();
if (!verbose) ::testing::internal::CaptureStdout();
lmp->input->one("units metal");
lmp->input->one("read_data test_pair_unit_convert.data");
lmp->input->one("pair_style gw");
lmp->input->one("pair_coeff * * SiC.gw Si C");
lmp->input->one("run 0 post no");
if (!verbose) ::testing::internal::GetCapturedStdout();
// copy pressure, energy, and force from first step
double pold;
lmp->output->thermo->evaluate_keyword("press", &pold);
double eold = lmp->force->pair->eng_vdwl + lmp->force->pair->eng_coul;
double **f = lmp->atom->f;
for (int i = 0; i < 4; ++i)
for (int j = 0; j < 3; ++j)
fold[i][j] = f[i][j];
if (!verbose) ::testing::internal::CaptureStdout();
lmp->input->one("clear");
lmp->input->one("units real");
lmp->input->one("read_data test_pair_unit_convert.data");
lmp->input->one("pair_style gw");
lmp->input->one("pair_coeff * * SiC.gw Si C");
lmp->input->one("run 0 post no");
if (!verbose) ::testing::internal::GetCapturedStdout();
double pnew;
lmp->output->thermo->evaluate_keyword("press", &pnew);
EXPECT_NEAR(pold, p_convert * pnew, fabs(pnew * rel_error));
double enew = lmp->force->pair->eng_vdwl + lmp->force->pair->eng_coul;
EXPECT_NEAR(ev_convert * eold, enew, fabs(enew * rel_error));
f = lmp->atom->f;
for (int i = 0; i < 4; ++i)
for (int j = 0; j < 3; ++j)
EXPECT_NEAR(ev_convert * fold[i][j], f[i][j], fabs(f[i][j] * rel_error));
}
TEST_F(PairUnitConvertTest, gw_zbl)
{
// check if the prerequisite pair style is available
if (!info->has_style("pair", "gw/zbl")) GTEST_SKIP();
if (!verbose) ::testing::internal::CaptureStdout();
lmp->input->one("units metal");
lmp->input->one("read_data test_pair_unit_convert.data");
lmp->input->one("pair_style gw/zbl");
lmp->input->one("pair_coeff * * SiC.gw.zbl Si C");
lmp->input->one("run 0 post no");
if (!verbose) ::testing::internal::GetCapturedStdout();
// copy pressure, energy, and force from first step
double pold;
lmp->output->thermo->evaluate_keyword("press", &pold);
double eold = lmp->force->pair->eng_vdwl + lmp->force->pair->eng_coul;
double **f = lmp->atom->f;
for (int i = 0; i < 4; ++i)
for (int j = 0; j < 3; ++j)
fold[i][j] = f[i][j];
if (!verbose) ::testing::internal::CaptureStdout();
lmp->input->one("clear");
lmp->input->one("units real");
lmp->input->one("read_data test_pair_unit_convert.data");
lmp->input->one("pair_style gw/zbl");
lmp->input->one("pair_coeff * * SiC.gw.zbl Si C");
lmp->input->one("run 0 post no");
if (!verbose) ::testing::internal::GetCapturedStdout();
double pnew;
lmp->output->thermo->evaluate_keyword("press", &pnew);
EXPECT_NEAR(pold, p_convert * pnew, fabs(pnew * rel_error));
double enew = lmp->force->pair->eng_vdwl + lmp->force->pair->eng_coul;
EXPECT_NEAR(ev_convert * eold, enew, fabs(enew * rel_error));
f = lmp->atom->f;
for (int i = 0; i < 4; ++i)
for (int j = 0; j < 3; ++j)
EXPECT_NEAR(ev_convert * fold[i][j], f[i][j], fabs(f[i][j] * rel_error));
}
TEST_F(PairUnitConvertTest, nb3b_harmonic)
{
// check if the prerequisite pair style is available
if (!info->has_style("pair", "nb3b/harmonic")) GTEST_SKIP();
if (!verbose) ::testing::internal::CaptureStdout();
lmp->input->one("units metal");
lmp->input->one("read_data test_pair_unit_convert.data");
lmp->input->one("pair_style nb3b/harmonic");
lmp->input->one("pair_coeff * * MOH.nb3b.harmonic M O");
lmp->input->one("run 0 post no");
if (!verbose) ::testing::internal::GetCapturedStdout();
// copy pressure, energy, and force from first step
double pold;
lmp->output->thermo->evaluate_keyword("press", &pold);
double eold = lmp->force->pair->eng_vdwl + lmp->force->pair->eng_coul;
double **f = lmp->atom->f;
for (int i = 0; i < 4; ++i)
for (int j = 0; j < 3; ++j)
fold[i][j] = f[i][j];
if (!verbose) ::testing::internal::CaptureStdout();
lmp->input->one("clear");
lmp->input->one("units real");
lmp->input->one("read_data test_pair_unit_convert.data");
lmp->input->one("pair_style nb3b/harmonic");
lmp->input->one("pair_coeff * * MOH.nb3b.harmonic M O");
lmp->input->one("run 0 post no");
if (!verbose) ::testing::internal::GetCapturedStdout();
double pnew;
lmp->output->thermo->evaluate_keyword("press", &pnew);
EXPECT_NEAR(pold, p_convert * pnew, fabs(pnew * rel_error));
double enew = lmp->force->pair->eng_vdwl + lmp->force->pair->eng_coul;
EXPECT_NEAR(ev_convert * eold, enew, fabs(enew * rel_error));
f = lmp->atom->f;
for (int i = 0; i < 4; ++i)
for (int j = 0; j < 3; ++j)
EXPECT_NEAR(ev_convert * fold[i][j], f[i][j], fabs(f[i][j] * rel_error));
}
TEST_F(PairUnitConvertTest, sw)
{
// check if the prerequisite pair style is available