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

Merge pull request #2638 from akohlmey/parse-molecule-refactor

Browse Source 
Further refactor parsing of molecule files
This commit is contained in:
Axel Kohlmeyer
2021-03-08 20:04:07 -05:00
committed by GitHub
parent f17ac30875 f455869de3
commit 37afa0e6d5
7 changed files with 653 additions and 220 deletions
Download Patch File Download Diff File Expand all files Collapse all files

468
src/molecule.cpp
View File

@ -55,14 +55,10 @@ Molecule::Molecule(LAMMPS *lmp, int narg, char **arg, int &index) :
if (index >= narg) error->all(FLERR,"Illegal molecule command"); if (index >= narg) error->all(FLERR,"Illegal molecule command");
int n = strlen(arg[0]) + 1; id = utils::strdup(arg[0]);
id = new char[n]; if (!utils::is_id(id))
strcpy(id,arg[0]); error->all(FLERR,"Molecule template ID must have only "
"alphanumeric or underscore characters");
for (int i = 0; i < n-1; i++)
if (!isalnum(id[i]) && id[i] != '_')
error->all(FLERR,"Molecule template ID must be "
"alphanumeric or underscore characters");
// parse args until reach unknown arg (next file) // parse args until reach unknown arg (next file)
@ -129,16 +125,21 @@ Molecule::Molecule(LAMMPS *lmp, int narg, char **arg, int &index) :
initialize(); initialize();
// scan file for sizes of all fields and allocate them // scan file for sizes of all fields and allocate storage for them
if (me == 0) open(arg[ifile]); if (me == 0) {
fp = fopen(arg[ifile],"r");
if (fp == nullptr)
error->one(FLERR,fmt::format("Cannot open molecule file {}: {}",
arg[ifile], utils::getsyserror()));
}
read(0); read(0);
if (me == 0) fclose(fp); if (me == 0) fclose(fp);
allocate(); allocate();
// read file again to populate all fields // read file again to populate all fields
if (me == 0) open(arg[ifile]); if (me == 0) fp = fopen(arg[ifile],"r");
read(1); read(1);
if (me == 0) fclose(fp); if (me == 0) fclose(fp);
@ -397,8 +398,8 @@ void Molecule::compute_inertia()
void Molecule::read(int flag) void Molecule::read(int flag)
{ {
char line[MAXLINE],keyword[MAXLINE]; char line[MAXLINE];
char *eof,*ptr; char *eof;
// skip 1st line of file // skip 1st line of file
@ -415,15 +416,14 @@ void Molecule::read(int flag)
readline(line); readline(line);
// trim anything from '#' onward // trim comments. if line is blank, continue
// if line is blank, continue
if ((ptr = strchr(line,'#'))) *ptr = '\0'; auto text = utils::trim(utils::trim_comment(line));
if (strspn(line," \t\n\r") == strlen(line)) continue; if (text.empty()) continue;
// search line for header keywords and set corresponding variable // search line for header keywords and set corresponding variable
try { try {
ValueTokenizer values(line); ValueTokenizer values(text);
int nmatch = values.count(); int nmatch = values.count();
int nwant = 0; int nwant = 0;
@ -485,10 +485,16 @@ void Molecule::read(int flag)
nibody = values.next_int(); nibody = values.next_int();
ndbody = values.next_int(); ndbody = values.next_int();
nwant = 3; nwant = 3;
} else break; } else {
// unknown header keyword
if (utils::strmatch(text,"^\\d+\\s+\\S+")) {
values.next_int();
auto keyword = values.next_string();
error->one(FLERR,fmt::format("Invalid header keyword: {}",keyword));
} else break;
}
if (nmatch != nwant) if (nmatch != nwant)
error->one(FLERR,"Invalid header in molecule file"); error->one(FLERR,"Invalid header line format in molecule file");
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid header in molecule file\n" error->one(FLERR, fmt::format("Invalid header in molecule file\n"
"{}", e.what())); "{}", e.what()));
@ -498,7 +504,7 @@ void Molecule::read(int flag)
// error checks // error checks
if (natoms < 1) if (natoms < 1)
error->all(FLERR,"No count or invalid atom count in molecule file"); error->all(FLERR,"No atoms or invalid atom count in molecule file");
if (nbonds < 0) error->all(FLERR,"Invalid bond count in molecule file"); if (nbonds < 0) error->all(FLERR,"Invalid bond count in molecule file");
if (nangles < 0) error->all(FLERR,"Invalid angle count in molecule file"); if (nangles < 0) error->all(FLERR,"Invalid angle count in molecule file");
if (ndihedrals < 0) if (ndihedrals < 0)
@ -512,108 +518,113 @@ void Molecule::read(int flag)
// grab keyword and skip next line // grab keyword and skip next line
parse_keyword(0,line,keyword); std::string keyword = parse_keyword(0,line);
readline(line); readline(line);
// loop over sections of molecule file // loop over sections of molecule file
while (strlen(keyword) > 0) { while (!keyword.empty()) {
if (strcmp(keyword,"Coords") == 0) { if (keyword == "Coords") {
xflag = 1; xflag = 1;
if (flag) coords(line); if (flag) coords(line);
else skip_lines(natoms,line); else skip_lines(natoms,line,keyword);
} else if (strcmp(keyword,"Types") == 0) { } else if (keyword == "Types") {
typeflag = 1; typeflag = 1;
if (flag) types(line); if (flag) types(line);
else skip_lines(natoms,line); else skip_lines(natoms,line,keyword);
} else if (strcmp(keyword,"Molecules") == 0) { } else if (keyword == "Molecules") {
moleculeflag = 1; moleculeflag = 1;
if (flag) molecules(line); if (flag) molecules(line);
else skip_lines(natoms,line); else skip_lines(natoms,line,keyword);
} else if (strcmp(keyword,"Fragments") == 0) { } else if (keyword == "Fragments") {
if (nfragments == 0) if (nfragments == 0)
error->all(FLERR,"Molecule file has fragments but no nfragments setting"); error->all(FLERR,"Found Fragments section but no nfragments setting in header");
fragmentflag = 1; fragmentflag = 1;
if (flag) fragments(line); if (flag) fragments(line);
else skip_lines(nfragments,line); else skip_lines(nfragments,line,keyword);
} else if (strcmp(keyword,"Charges") == 0) { } else if (keyword == "Charges") {
qflag = 1; qflag = 1;
if (flag) charges(line); if (flag) charges(line);
else skip_lines(natoms,line); else skip_lines(natoms,line,keyword);
} else if (strcmp(keyword,"Diameters") == 0) { } else if (keyword == "Diameters") {
radiusflag = 1; radiusflag = 1;
if (flag) diameters(line); if (flag) diameters(line);
else skip_lines(natoms,line); else skip_lines(natoms,line,keyword);
} else if (strcmp(keyword,"Masses") == 0) { } else if (keyword == "Masses") {
rmassflag = 1; rmassflag = 1;
if (flag) masses(line); if (flag) masses(line);
else skip_lines(natoms,line); else skip_lines(natoms,line,keyword);
} else if (strcmp(keyword,"Bonds") == 0) { } else if (keyword == "Bonds") {
if (nbonds == 0) if (nbonds == 0)
error->all(FLERR,"Molecule file has bonds but no nbonds setting"); error->all(FLERR,"Found Bonds section but no nbonds setting in header");
bondflag = tag_require = 1; bondflag = tag_require = 1;
bonds(flag,line); bonds(flag,line);
} else if (strcmp(keyword,"Angles") == 0) { } else if (keyword == "Angles") {
if (nangles == 0) if (nangles == 0)
error->all(FLERR,"Molecule file has angles but no nangles setting"); error->all(FLERR,"Found Angles section but no nangles setting in header");
angleflag = tag_require = 1; angleflag = tag_require = 1;
angles(flag,line); angles(flag,line);
} else if (strcmp(keyword,"Dihedrals") == 0) { } else if (keyword == "Dihedrals") {
if (ndihedrals == 0) error->all(FLERR,"Molecule file has dihedrals " if (ndihedrals == 0) error->all(FLERR,"Found Dihedrals section"
"but no ndihedrals setting"); "but no ndihedrals setting in header");
dihedralflag = tag_require = 1; dihedralflag = tag_require = 1;
dihedrals(flag,line); dihedrals(flag,line);
} else if (strcmp(keyword,"Impropers") == 0) { } else if (keyword == "Impropers") {
if (nimpropers == 0) error->all(FLERR,"Molecule file has impropers " if (nimpropers == 0) error->all(FLERR,"Found Impropers section"
"but no nimpropers setting"); "but no nimpropers setting in header");
improperflag = tag_require = 1; improperflag = tag_require = 1;
impropers(flag,line); impropers(flag,line);
} else if (strcmp(keyword,"Special Bond Counts") == 0) { } else if (keyword == "Special Bond Counts") {
nspecialflag = 1; nspecialflag = 1;
nspecial_read(flag,line); nspecial_read(flag,line);
} else if (strcmp(keyword,"Special Bonds") == 0) { } else if (keyword == "Special Bonds") {
specialflag = tag_require = 1; specialflag = tag_require = 1;
if (flag) special_read(line); if (flag) special_read(line);
else skip_lines(natoms,line); else skip_lines(natoms,line,keyword);
} else if (strcmp(keyword,"Shake Flags") == 0) { } else if (keyword == "Shake Flags") {
shakeflagflag = 1; shakeflagflag = 1;
if (flag) shakeflag_read(line); if (flag) shakeflag_read(line);
else skip_lines(natoms,line); else skip_lines(natoms,line,keyword);
} else if (strcmp(keyword,"Shake Atoms") == 0) { } else if (keyword == "Shake Atoms") {
shakeatomflag = tag_require = 1; shakeatomflag = tag_require = 1;
if (shaketypeflag) shakeflag = 1; if (shaketypeflag) shakeflag = 1;
if (!shakeflagflag) if (!shakeflagflag)
error->all(FLERR,"Molecule file shake flags not before shake atoms"); error->all(FLERR,"Shake Flags section must come before Shake Atoms section");
if (flag) shakeatom_read(line); if (flag) shakeatom_read(line);
else skip_lines(natoms,line); else skip_lines(natoms,line,keyword);
} else if (strcmp(keyword,"Shake Bond Types") == 0) { } else if (keyword == "Shake Bond Types") {
shaketypeflag = 1; shaketypeflag = 1;
if (shakeatomflag) shakeflag = 1; if (shakeatomflag) shakeflag = 1;
if (!shakeflagflag) if (!shakeflagflag)
error->all(FLERR,"Molecule file shake flags not before shake bonds"); error->all(FLERR,"Shake Flags section must come before Shake Bonds section");
if (flag) shaketype_read(line); if (flag) shaketype_read(line);
else skip_lines(natoms,line); else skip_lines(natoms,line,keyword);
} else if (strcmp(keyword,"Body Integers") == 0) { } else if (keyword == "Body Integers") {
if (bodyflag == 0 || nibody == 0) if (bodyflag == 0 || nibody == 0)
error->all(FLERR,"Molecule file has body params " error->all(FLERR,"Found Body Integers section but no setting in header");
"but no setting for them");
ibodyflag = 1; ibodyflag = 1;
body(flag,0,line); body(flag,0,line);
} else if (strcmp(keyword,"Body Doubles") == 0) { } else if (keyword == "Body Doubles") {
if (bodyflag == 0 || ndbody == 0) if (bodyflag == 0 || ndbody == 0)
error->all(FLERR,"Molecule file has body params " error->all(FLERR,"Found Body Doubles section but no setting in header");
"but no setting for them");
dbodyflag = 1; dbodyflag = 1;
body(flag,1,line); body(flag,1,line);
} else {
} else error->one(FLERR,fmt::format("Unknown section '{}' in molecule " // Error: Either a too long/short section or a typo in the keyword
"file", keyword));
parse_keyword(1,line,keyword); if (utils::strmatch(keyword,"^[A-Za-z ]+$"))
error->one(FLERR,fmt::format("Unknown section '{}' in molecule "
"file\n",keyword));
else error->one(FLERR,fmt::format("Unexpected line in molecule file "
"while looking for the next "
"section:\n{}",line));
}
keyword = parse_keyword(1,line);
} }
// error check // error check
@ -678,11 +689,14 @@ void Molecule::coords(char *line)
for (int i = 0; i < natoms; i++) { for (int i = 0; i < natoms; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 4) error->one(FLERR,"Invalid Coords section in molecule file"); if (values.count() != 4)
error->all(FLERR,fmt::format("Invalid line in Coords section of "
"molecule file: {}",line));
int iatom = values.next_int() - 1; int iatom = values.next_int() - 1;
if (iatom < 0 || iatom >= natoms) error->one(FLERR,"Invalid Coords section in molecule file"); if (iatom < 0 || iatom >= natoms)
error->all(FLERR,"Invalid atom index in Coords section of molecule file");
count[iatom]++; count[iatom]++;
x[iatom][0] = values.next_double(); x[iatom][0] = values.next_double();
x[iatom][1] = values.next_double(); x[iatom][1] = values.next_double();
@ -693,17 +707,19 @@ void Molecule::coords(char *line)
x[iatom][2] *= sizescale; x[iatom][2] *= sizescale;
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Coords section in molecule file\n" error->all(FLERR,fmt::format("Invalid line in Coords section of "
"{}", e.what())); "molecule file: {}\n{}",e.what(),line));
} }
for (int i = 0; i < natoms; i++) for (int i = 0; i < natoms; i++)
if (count[i] == 0) error->all(FLERR,"Invalid Coords section in molecule file"); if (count[i] == 0) error->all(FLERR,fmt::format("Atom {} missing in Coords "
"section of molecule file",i+1));
if (domain->dimension == 2) { if (domain->dimension == 2) {
for (int i = 0; i < natoms; i++) for (int i = 0; i < natoms; i++)
if (x[i][2] != 0.0) if (x[i][2] != 0.0)
error->all(FLERR,"Molecule file z coord must be 0.0 for 2d"); error->all(FLERR,fmt::format("Z coord in molecule file for atom {} "
"must be 0.0 for 2d-simulation.",i+1));
} }
} }
@ -719,29 +735,33 @@ void Molecule::types(char *line)
for (int i = 0; i < natoms; i++) { for (int i = 0; i < natoms; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 2) error->one(FLERR,"Invalid Types section in molecule file"); if (values.count() != 2)
error->all(FLERR,fmt::format("Invalid line in Types section of "
"molecule file: {}",line));
int iatom = values.next_int() - 1; int iatom = values.next_int() - 1;
if (iatom < 0 || iatom >= natoms) error->one(FLERR,"Invalid Types section in molecule file"); if (iatom < 0 || iatom >= natoms)
error->all(FLERR,"Invalid atom index in Types section of molecule file");
count[iatom]++; count[iatom]++;
type[iatom] = values.next_int(); type[iatom] = values.next_int();
type[iatom] += toffset; type[iatom] += toffset;
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Types section in molecule file\n" error->all(FLERR, fmt::format("Invalid line in Types section of "
"{}", e.what())); "molecule file: {}\n{}", e.what(),line));
} }
for (int i = 0; i < natoms; i++) for (int i = 0; i < natoms; i++) {
if (count[i] == 0) error->all(FLERR,"Invalid Types section in molecule file"); if (count[i] == 0) error->all(FLERR,fmt::format("Atom {} missing in Types "
"section of molecule file",i+1));
for (int i = 0; i < natoms; i++)
if ((type[i] <= 0) || (domain->box_exist && (type[i] > atom->ntypes))) if ((type[i] <= 0) || (domain->box_exist && (type[i] > atom->ntypes)))
error->all(FLERR,"Invalid atom type in molecule file"); error->all(FLERR,fmt::format("Invalid atom type {} for atom {} "
"in molecule file",type[i],i+1));
for (int i = 0; i < natoms; i++)
ntypes = MAX(ntypes,type[i]); ntypes = MAX(ntypes,type[i]);
}
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
@ -755,26 +775,31 @@ void Molecule::molecules(char *line)
try { try {
for (int i = 0; i < natoms; i++) { for (int i = 0; i < natoms; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 2) error->one(FLERR,"Invalid Molecules section in molecule file"); if (values.count() != 2)
error->all(FLERR,fmt::format("Invalid line in Molecules section of "
"molecule file: {}",line));
int iatom = values.next_int() - 1; int iatom = values.next_int() - 1;
if (iatom < 0 || iatom >= natoms) error->one(FLERR,"Invalid Molecules section in molecule file"); if (iatom < 0 || iatom >= natoms)
error->all(FLERR,"Invalid atom index in Molecules section of molecule file");
count[iatom]++; count[iatom]++;
molecule[iatom] = values.next_tagint(); molecule[iatom] = values.next_tagint();
// molecule[iatom] += moffset; // placeholder for possible molecule offset // molecule[iatom] += moffset; // placeholder for possible molecule offset
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Molecules section in molecule file\n" error->all(FLERR, fmt::format("Invalid line in Molecules section of "
"{}", e.what())); "molecule file: {}\n{}",e.what(),line));
} }
for (int i = 0; i < natoms; i++) for (int i = 0; i < natoms; i++)
if (count[i] == 0) error->all(FLERR,"Invalid Molecules section in molecule file"); if (count[i] == 0) error->all(FLERR,fmt::format("Atom {} missing in Molecules "
"section of molecule file",i+1));
for (int i = 0; i < natoms; i++) for (int i = 0; i < natoms; i++)
if (molecule[i] <= 0) if (molecule[i] < 0)
error->all(FLERR,"Invalid molecule ID in molecule file"); error->all(FLERR,fmt::format("Invalid molecule ID {} for atom {} "
"in molecule file",molecule[i],i+1));
for (int i = 0; i < natoms; i++) for (int i = 0; i < natoms; i++)
nmolecules = MAX(nmolecules,molecule[i]); nmolecules = MAX(nmolecules,molecule[i]);
@ -790,23 +815,26 @@ void Molecule::fragments(char *line)
for (int i = 0; i < nfragments; i++) { for (int i = 0; i < nfragments; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if ((int)values.count() > natoms+1) if ((int)values.count() > natoms+1)
error->one(FLERR,"Invalid atom ID in Fragments section of molecule file"); error->all(FLERR,"Too many atoms per fragment in Fragments "
"section of molecule file");
fragmentnames[i] = values.next_string(); fragmentnames[i] = values.next_string();
while (values.has_next()) { while (values.has_next()) {
int atomID = values.next_int(); int iatom = values.next_int()-1;
if (atomID <= 0 || atomID > natoms) if (iatom < 0 || iatom >= natoms)
error->one(FLERR,"Invalid atom ID in Fragments section of molecule file"); error->all(FLERR,fmt::format("Invalid atom ID {} for fragment {} in "
fragmentmask[i][atomID-1] = 1; "Fragments section of molecule file",
iatom+1, fragmentnames[i]));
fragmentmask[i][iatom] = 1;
} }
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid atom ID in Fragments section of molecule file\n" error->all(FLERR, fmt::format("Invalid atom ID in Fragments section of "
"{}", e.what())); "molecule file: {}\n{}", e.what(),line));
} }
} }
@ -821,21 +849,26 @@ void Molecule::charges(char *line)
for (int i = 0; i < natoms; i++) { for (int i = 0; i < natoms; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if ((int)values.count() != 2) error->one(FLERR,"Invalid Charges section in molecule file"); if ((int)values.count() != 2)
error->all(FLERR,fmt::format("Invalid line in Charges section of "
"molecule file: {}",line));
int iatom = values.next_int() - 1; int iatom = values.next_int() - 1;
if (iatom < 0 || iatom >= natoms) error->one(FLERR,"Invalid Charges section in molecule file"); if (iatom < 0 || iatom >= natoms)
error->all(FLERR,"Invalid atom index in Charges section of molecule file");
count[iatom]++; count[iatom]++;
q[iatom] = values.next_double(); q[iatom] = values.next_double();
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Charges section in molecule file\n" error->all(FLERR, fmt::format("Invalid line in Charges section of "
"{}", e.what())); "molecule file: {}.\n{}",e.what(),line));
} }
for (int i = 0; i < natoms; i++) for (int i = 0; i < natoms; i++)
if (count[i] == 0) error->all(FLERR,"Invalid Charges section in molecule file"); if (count[i] == 0) error->all(FLERR,fmt::format("Atom {} missing in Charges "
"section of molecule file",i+1));
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
@ -850,11 +883,13 @@ void Molecule::diameters(char *line)
for (int i = 0; i < natoms; i++) { for (int i = 0; i < natoms; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 2) error->one(FLERR,"Invalid Diameters section in molecule file"); if (values.count() != 2)
error->all(FLERR,fmt::format("Invalid line in Diameters section of "
"molecule file: {}",line));
int iatom = values.next_int() - 1; int iatom = values.next_int() - 1;
if (iatom < 0 || iatom >= natoms) error->one(FLERR,"Invalid Diameters section in molecule file"); if (iatom < 0 || iatom >= natoms)
error->all(FLERR,"Invalid atom index in Diameters section of molecule file");
count[iatom]++; count[iatom]++;
radius[iatom] = values.next_double(); radius[iatom] = values.next_double();
radius[iatom] *= sizescale; radius[iatom] *= sizescale;
@ -862,16 +897,17 @@ void Molecule::diameters(char *line)
maxradius = MAX(maxradius,radius[iatom]); maxradius = MAX(maxradius,radius[iatom]);
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Diameters section in molecule file\n" error->all(FLERR, fmt::format("Invalid line in Diameters section of "
"{}", e.what())); "molecule file: {}\n{}",e.what(),line));
} }
for (int i = 0; i < natoms; i++) for (int i = 0; i < natoms; i++) {
if (count[i] == 0) error->all(FLERR,"Invalid Diameters section in molecule file"); if (count[i] == 0) error->all(FLERR,fmt::format("Atom {} missing in Diameters "
"section of molecule file",i+1));
for (int i = 0; i < natoms; i++)
if (radius[i] < 0.0) if (radius[i] < 0.0)
error->all(FLERR,"Invalid atom diameter in molecule file"); error->all(FLERR,fmt::format("Invalid atom diameter {} for atom {} "
"in molecule file", radius[i], i+1));
}
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
@ -885,25 +921,30 @@ void Molecule::masses(char *line)
for (int i = 0; i < natoms; i++) { for (int i = 0; i < natoms; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 2) error->one(FLERR,"Invalid Masses section in molecule file"); if (values.count() != 2)
error->all(FLERR,fmt::format("Invalid line in Masses section of "
"molecule file: {}",line));
int iatom = values.next_int() - 1; int iatom = values.next_int() - 1;
if (iatom < 0 || iatom >= natoms) error->one(FLERR,"Invalid Masses section in molecule file"); if (iatom < 0 || iatom >= natoms)
error->all(FLERR,"Invalid atom index in Masses section of molecule file");
count[iatom]++; count[iatom]++;
rmass[iatom] = values.next_double(); rmass[iatom] = values.next_double();
rmass[iatom] *= sizescale*sizescale*sizescale; rmass[iatom] *= sizescale*sizescale*sizescale;
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Masses section in molecule file\n" error->all(FLERR, fmt::format("Invalid line in Masses section of "
"{}", e.what())); "molecule file: {}\n{}",e.what(),line));
} }
for (int i = 0; i < natoms; i++) for (int i = 0; i < natoms; i++) {
if (count[i] == 0) error->all(FLERR,"Invalid Masses section in molecule file"); if (count[i] == 0) error->all(FLERR,fmt::format("Atom {} missing in Masses "
"section of molecule file",i+1));
for (int i = 0; i < natoms; i++) if (rmass[i] <= 0.0)
if (rmass[i] <= 0.0) error->all(FLERR,"Invalid atom mass in molecule file"); error->all(FLERR,fmt::format("Invalid atom mass {} for atom {} "
"in molecule file", radius[i], i+1));
}
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
@ -929,24 +970,26 @@ void Molecule::bonds(int flag, char *line)
readline(line); readline(line);
try { try {
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 4) error->one(FLERR,"Invalid Bonds section in molecule file"); if (values.count() != 4)
error->all(FLERR,fmt::format("Invalid line in Bonds section of "
"molecule file: {}",line));
values.next_int(); values.next_int();
itype = values.next_int(); itype = values.next_int();
atom1 = values.next_tagint(); atom1 = values.next_tagint();
atom2 = values.next_tagint(); atom2 = values.next_tagint();
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Bonds section in molecule file\n" error->all(FLERR, fmt::format("Invalid line in Bonds section of "
"{}", e.what())); "molecule file: {}\n{}",e.what(),line));
} }
itype += boffset; itype += boffset;
if ((atom1 <= 0) || (atom1 > natoms) || if ((atom1 <= 0) || (atom1 > natoms) ||
(atom2 <= 0) || (atom2 > natoms) || (atom1 == atom2)) (atom2 <= 0) || (atom2 > natoms) || (atom1 == atom2))
error->one(FLERR,"Invalid atom ID in Bonds section of molecule file"); error->all(FLERR,"Invalid atom ID in Bonds section of molecule file");
if ((itype <= 0) || (domain->box_exist && (itype > atom->nbondtypes))) if ((itype <= 0) || (domain->box_exist && (itype > atom->nbondtypes)))
error->one(FLERR,"Invalid bond type in Bonds section of molecule file"); error->all(FLERR,"Invalid bond type in Bonds section of molecule file");
if (flag) { if (flag) {
m = atom1-1; m = atom1-1;
@ -997,16 +1040,18 @@ void Molecule::angles(int flag, char *line)
readline(line); readline(line);
try { try {
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 5) error->one(FLERR,"Invalid Angles section in molecule file"); if (values.count() != 5)
error->all(FLERR,fmt::format("Invalid line in Angles section of "
"molecule file: {}",line));
values.next_int(); values.next_int();
itype = values.next_int(); itype = values.next_int();
atom1 = values.next_tagint(); atom1 = values.next_tagint();
atom2 = values.next_tagint(); atom2 = values.next_tagint();
atom3 = values.next_tagint(); atom3 = values.next_tagint();
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Angles section in molecule file\n" error->all(FLERR, fmt::format("Invalid line in Angles section of "
"{}", e.what())); "molecule file: {}\n{}",e.what(),line));
} }
itype += aoffset; itype += aoffset;
@ -1015,9 +1060,9 @@ void Molecule::angles(int flag, char *line)
(atom2 <= 0) || (atom2 > natoms) || (atom2 <= 0) || (atom2 > natoms) ||
(atom3 <= 0) || (atom3 > natoms) || (atom3 <= 0) || (atom3 > natoms) ||
(atom1 == atom2) || (atom1 == atom3) || (atom2 == atom3)) (atom1 == atom2) || (atom1 == atom3) || (atom2 == atom3))
error->one(FLERR,"Invalid atom ID in Angles section of molecule file"); error->all(FLERR,"Invalid atom ID in Angles section of molecule file");
if ((itype <= 0) || (domain->box_exist && (itype > atom->nangletypes))) if ((itype <= 0) || (domain->box_exist && (itype > atom->nangletypes)))
error->one(FLERR,"Invalid angle type in Angles section of molecule file"); error->all(FLERR,"Invalid angle type in Angles section of molecule file");
if (flag) { if (flag) {
m = atom2-1; m = atom2-1;
@ -1081,8 +1126,11 @@ void Molecule::dihedrals(int flag, char *line)
readline(line); readline(line);
try { try {
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 6) error->one(FLERR,"Invalid Dihedrals section in molecule file"); if (values.count() != 6)
error->all(FLERR,fmt::format("Invalid line in Dihedrals section of "
"molecule file: {}",line));
values.next_int(); values.next_int();
itype = values.next_int(); itype = values.next_int();
atom1 = values.next_tagint(); atom1 = values.next_tagint();
@ -1090,8 +1138,8 @@ void Molecule::dihedrals(int flag, char *line)
atom3 = values.next_tagint(); atom3 = values.next_tagint();
atom4 = values.next_tagint(); atom4 = values.next_tagint();
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Dihedrals section in molecule file\n" error->all(FLERR, fmt::format("Invalid line in Dihedrals section of "
"{}", e.what())); "molecule file: {}\n{}",e.what(),line));
} }
itype += doffset; itype += doffset;
@ -1102,10 +1150,10 @@ void Molecule::dihedrals(int flag, char *line)
(atom4 <= 0) || (atom4 > natoms) || (atom4 <= 0) || (atom4 > natoms) ||
(atom1 == atom2) || (atom1 == atom3) || (atom1 == atom4) || (atom1 == atom2) || (atom1 == atom3) || (atom1 == atom4) ||
(atom2 == atom3) || (atom2 == atom4) || (atom3 == atom4)) (atom2 == atom3) || (atom2 == atom4) || (atom3 == atom4))
error->one(FLERR, error->all(FLERR,
"Invalid atom ID in dihedrals section of molecule file"); "Invalid atom ID in dihedrals section of molecule file");
if ((itype <= 0) || (domain->box_exist && (itype > atom->ndihedraltypes))) if ((itype <= 0) || (domain->box_exist && (itype > atom->ndihedraltypes)))
error->one(FLERR,"Invalid dihedral type in Dihedrals section of molecule file"); error->all(FLERR,"Invalid dihedral type in Dihedrals section of molecule file");
if (flag) { if (flag) {
m = atom2-1; m = atom2-1;
@ -1180,8 +1228,10 @@ void Molecule::impropers(int flag, char *line)
readline(line); readline(line);
try { try {
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 6) error->one(FLERR,"Invalid Impropers section in molecule file"); if (values.count() != 6)
error->all(FLERR,fmt::format("Invalid line in Impropers section of "
"molecule file: {}",line));
values.next_int(); values.next_int();
itype = values.next_int(); itype = values.next_int();
atom1 = values.next_tagint(); atom1 = values.next_tagint();
@ -1189,8 +1239,8 @@ void Molecule::impropers(int flag, char *line)
atom3 = values.next_tagint(); atom3 = values.next_tagint();
atom4 = values.next_tagint(); atom4 = values.next_tagint();
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Impropers section in molecule file\n" error->all(FLERR, fmt::format("Invalid line in Impropers section of "
"{}", e.what())); "molecule file: {}\n{}",e.what(),line));
} }
itype += ioffset; itype += ioffset;
@ -1201,10 +1251,10 @@ void Molecule::impropers(int flag, char *line)
(atom4 <= 0) || (atom4 > natoms) || (atom4 <= 0) || (atom4 > natoms) ||
(atom1 == atom2) || (atom1 == atom3) || (atom1 == atom4) || (atom1 == atom2) || (atom1 == atom3) || (atom1 == atom4) ||
(atom2 == atom3) || (atom2 == atom4) || (atom3 == atom4)) (atom2 == atom3) || (atom2 == atom4) || (atom3 == atom4))
error->one(FLERR, error->all(FLERR,
"Invalid atom ID in impropers section of molecule file"); "Invalid atom ID in impropers section of molecule file");
if ((itype <= 0) || (domain->box_exist && (itype > atom->nimpropertypes))) if ((itype <= 0) || (domain->box_exist && (itype > atom->nimpropertypes)))
error->one(FLERR,"Invalid improper type in Impropers section of molecule file"); error->all(FLERR,"Invalid improper type in Impropers section of molecule file");
if (flag) { if (flag) {
m = atom2-1; m = atom2-1;
@ -1273,15 +1323,17 @@ void Molecule::nspecial_read(int flag, char *line)
int c1, c2, c3; int c1, c2, c3;
try { try {
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 4) error->one(FLERR,"Invalid Special Bond Counts section in molecule file"); if (values.count() != 4)
error->all(FLERR,fmt::format("Invalid line in Special Bond Counts section of "
"molecule file: {}",line));
values.next_int(); values.next_int();
c1 = values.next_tagint(); c1 = values.next_tagint();
c2 = values.next_tagint(); c2 = values.next_tagint();
c3 = values.next_tagint(); c3 = values.next_tagint();
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Special Bond Counts section in molecule file\n" error->all(FLERR, fmt::format("Invalid line in Special Bond Counts section of "
"{}", e.what())); "molecule file: {}\n{}",e.what(),line));
} }
if (flag) { if (flag) {
@ -1302,11 +1354,11 @@ void Molecule::special_read(char *line)
for (int i = 0; i < natoms; i++) { for (int i = 0; i < natoms; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
int nwords = values.count(); int nwords = values.count();
if (nwords != nspecial[i][2]+1) if (nwords != nspecial[i][2]+1)
error->one(FLERR,"Molecule file special list " error->all(FLERR,"Molecule file special list "
"does not match special count"); "does not match special count");
values.next_int(); // ignore values.next_int(); // ignore
@ -1315,12 +1367,13 @@ void Molecule::special_read(char *line)
special[i][m-1] = values.next_tagint(); special[i][m-1] = values.next_tagint();
if (special[i][m-1] <= 0 || special[i][m-1] > natoms || if (special[i][m-1] <= 0 || special[i][m-1] > natoms ||
special[i][m-1] == i+1) special[i][m-1] == i+1)
error->one(FLERR,"Invalid special atom index in molecule file"); error->all(FLERR,"Invalid atom index in Special Bonds "
"section of molecule file");
} }
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Molecule file special list\n" error->all(FLERR, fmt::format("Invalid line in Special Bonds section of "
"{}", e.what())); "molecule file: {}\n{}",e.what(),line));
} }
} }
@ -1351,7 +1404,7 @@ void Molecule::special_generate()
nspecial[i][0]++; nspecial[i][0]++;
nspecial[atom2][0]++; nspecial[atom2][0]++;
if (count[i] >= atom->maxspecial || count[atom2] >= atom->maxspecial) if (count[i] >= atom->maxspecial || count[atom2] >= atom->maxspecial)
error->one(FLERR,"Molecule auto special bond generation overflow"); error->all(FLERR,"Molecule auto special bond generation overflow");
tmpspecial[i][count[i]++] = atom2 + 1; tmpspecial[i][count[i]++] = atom2 + 1;
tmpspecial[atom2][count[atom2]++] = i + 1; tmpspecial[atom2][count[atom2]++] = i + 1;
} }
@ -1363,7 +1416,7 @@ void Molecule::special_generate()
atom1 = i; atom1 = i;
atom2 = bond_atom[i][j]; atom2 = bond_atom[i][j];
if (count[atom1] >= atom->maxspecial) if (count[atom1] >= atom->maxspecial)
error->one(FLERR,"Molecule auto special bond generation overflow"); error->all(FLERR,"Molecule auto special bond generation overflow");
tmpspecial[i][count[atom1]++] = atom2; tmpspecial[i][count[atom1]++] = atom2;
} }
} }
@ -1386,7 +1439,7 @@ void Molecule::special_generate()
} }
if (!dedup) { if (!dedup) {
if (count[i] >= atom->maxspecial) if (count[i] >= atom->maxspecial)
error->one(FLERR,"Molecule auto special bond generation overflow"); error->all(FLERR,"Molecule auto special bond generation overflow");
tmpspecial[i][count[i]++] = tmpspecial[tmpspecial[i][m]-1][j]; tmpspecial[i][count[i]++] = tmpspecial[tmpspecial[i][m]-1][j];
nspecial[i][1]++; nspecial[i][1]++;
} }
@ -1410,7 +1463,7 @@ void Molecule::special_generate()
} }
if (!dedup) { if (!dedup) {
if (count[i] >= atom->maxspecial) if (count[i] >= atom->maxspecial)
error->one(FLERR,"Molecule auto special bond generation overflow"); error->all(FLERR,"Molecule auto special bond generation overflow");
tmpspecial[i][count[i]++] = tmpspecial[tmpspecial[i][m]-1][j]; tmpspecial[i][count[i]++] = tmpspecial[tmpspecial[i][m]-1][j];
nspecial[i][2]++; nspecial[i][2]++;
} }
@ -1440,22 +1493,22 @@ void Molecule::shakeflag_read(char *line)
for (int i = 0; i < natoms; i++) { for (int i = 0; i < natoms; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
if (values.count() != 2) if (values.count() != 2)
error->one(FLERR,"Invalid Shake Flags section in molecule file"); error->all(FLERR,"Invalid Shake Flags section in molecule file");
values.next_int(); values.next_int();
shake_flag[i] = values.next_int(); shake_flag[i] = values.next_int();
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid Shake Flags section in molecule file\n" error->all(FLERR, fmt::format("Invalid Shake Flags section in molecule file\n"
"{}", e.what())); "{}", e.what()));
} }
for (int i = 0; i < natoms; i++) for (int i = 0; i < natoms; i++)
if (shake_flag[i] < 0 || shake_flag[i] > 4) if (shake_flag[i] < 0 || shake_flag[i] > 4)
error->one(FLERR,"Invalid shake flag in molecule file"); error->all(FLERR,"Invalid shake flag in molecule file");
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
@ -1469,7 +1522,7 @@ void Molecule::shakeatom_read(char *line)
for (int i = 0; i < natoms; i++) { for (int i = 0; i < natoms; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
nmatch = values.count(); nmatch = values.count();
switch (shake_flag[i]) { switch (shake_flag[i]) {
@ -1511,15 +1564,15 @@ void Molecule::shakeatom_read(char *line)
break; break;
default: default:
error->one(FLERR,"Invalid shake atom in molecule file"); error->all(FLERR,"Invalid shake atom in molecule file");
} }
if (nmatch != nwant) if (nmatch != nwant)
error->one(FLERR,"Invalid shake atom in molecule file"); error->all(FLERR,"Invalid shake atom in molecule file");
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR,fmt::format("Invalid shake atom in molecule file\n" error->all(FLERR,fmt::format("Invalid shake atom in molecule file\n"
"{}", e.what())); "{}", e.what()));
} }
@ -1528,7 +1581,7 @@ void Molecule::shakeatom_read(char *line)
if (m == 1) m = 3; if (m == 1) m = 3;
for (int j = 0; j < m; j++) for (int j = 0; j < m; j++)
if (shake_atom[i][j] <= 0 || shake_atom[i][j] > natoms) if (shake_atom[i][j] <= 0 || shake_atom[i][j] > natoms)
error->one(FLERR,"Invalid shake atom in molecule file"); error->all(FLERR,"Invalid shake atom in molecule file");
} }
} }
@ -1543,7 +1596,7 @@ void Molecule::shaketype_read(char *line)
for (int i = 0; i < natoms; i++) { for (int i = 0; i < natoms; i++) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
nmatch = values.count(); nmatch = values.count();
switch (shake_flag[i]) { switch (shake_flag[i]) {
@ -1582,14 +1635,14 @@ void Molecule::shaketype_read(char *line)
break; break;
default: default:
error->one(FLERR,"Invalid shake type data in molecule file"); error->all(FLERR,"Invalid shake type data in molecule file");
} }
if (nmatch != nwant) if (nmatch != nwant)
error->one(FLERR,"Invalid shake type data in molecule file"); error->all(FLERR,"Invalid shake type data in molecule file");
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid shake type data in molecule file\n", error->all(FLERR, fmt::format("Invalid shake type data in molecule file\n",
"{}", e.what())); "{}", e.what()));
} }
@ -1598,10 +1651,10 @@ void Molecule::shaketype_read(char *line)
if (m == 1) m = 3; if (m == 1) m = 3;
for (int j = 0; j < m-1; j++) for (int j = 0; j < m-1; j++)
if (shake_type[i][j] <= 0) if (shake_type[i][j] <= 0)
error->one(FLERR,"Invalid shake bond type in molecule file"); error->all(FLERR,"Invalid shake bond type in molecule file");
if (shake_flag[i] == 1) if (shake_flag[i] == 1)
if (shake_type[i][2] <= 0) if (shake_type[i][2] <= 0)
error->one(FLERR,"Invalid shake angle type in molecule file"); error->all(FLERR,"Invalid shake angle type in molecule file");
} }
} }
@ -1621,13 +1674,13 @@ void Molecule::body(int flag, int pflag, char *line)
while (nword < nparam) { while (nword < nparam) {
readline(line); readline(line);
ValueTokenizer values(line); ValueTokenizer values(utils::trim_comment(line));
int ncount = values.count(); int ncount = values.count();
if (ncount == 0) if (ncount == 0)
error->one(FLERR,"Too few values in body section of molecule file"); error->all(FLERR,"Too few values in body section of molecule file");
if (nword+ncount > nparam) if (nword+ncount > nparam)
error->one(FLERR,"Too many values in body section of molecule file"); error->all(FLERR,"Too many values in body section of molecule file");
if (flag) { if (flag) {
if (pflag == 0) { if (pflag == 0) {
@ -1642,7 +1695,7 @@ void Molecule::body(int flag, int pflag, char *line)
} else nword += ncount; } else nword += ncount;
} }
} catch (TokenizerException &e) { } catch (TokenizerException &e) {
error->one(FLERR, fmt::format("Invalid body params in molecule file\n", error->all(FLERR, fmt::format("Invalid body params in molecule file\n",
"{}", e.what())); "{}", e.what()));
} }
} }
@ -1939,18 +1992,6 @@ void Molecule::deallocate()
memory->destroy(dbodyparams); memory->destroy(dbodyparams);
} }
/* ----------------------------------------------------------------------
open molecule file
------------------------------------------------------------------------- */
void Molecule::open(char *file)
{
fp = fopen(file,"r");
if (fp == nullptr)
error->one(FLERR,fmt::format("Cannot open molecule file {}: {}",
file, utils::getsyserror()));
}
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
read and bcast a line read and bcast a line
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
@ -1973,8 +2014,9 @@ void Molecule::readline(char *line)
flag = 1, line has already been read flag = 1, line has already been read
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Molecule::parse_keyword(int flag, char *line, char *keyword) std::string Molecule::parse_keyword(int flag, char *line)
{ {
char line2[MAXLINE];
if (flag) { if (flag) {
// read upto non-blank line plus 1 following line // read upto non-blank line plus 1 following line
@ -1986,42 +2028,38 @@ void Molecule::parse_keyword(int flag, char *line, char *keyword)
while (eof == 0 && strspn(line," \t\n\r") == strlen(line)) { while (eof == 0 && strspn(line," \t\n\r") == strlen(line)) {
if (fgets(line,MAXLINE,fp) == nullptr) eof = 1; if (fgets(line,MAXLINE,fp) == nullptr) eof = 1;
} }
if (fgets(keyword,MAXLINE,fp) == nullptr) eof = 1; if (fgets(line2,MAXLINE,fp) == nullptr) eof = 1;
} }
// if eof, set keyword empty and return // if eof, set keyword empty and return
MPI_Bcast(&eof,1,MPI_INT,0,world); MPI_Bcast(&eof,1,MPI_INT,0,world);
if (eof) { if (eof) {
keyword[0] = '\0'; return std::string("");
return;
} }
// bcast keyword line to all procs // bcast keyword line to all procs
int n; MPI_Bcast(line,MAXLINE,MPI_CHAR,0,world);
if (me == 0) n = strlen(line) + 1;
MPI_Bcast(&n,1,MPI_INT,0,world);
MPI_Bcast(line,n,MPI_CHAR,0,world);
} }
// copy non-whitespace and non-comment portion of line into keyword // return non-whitespace and non-comment portion of line
int start = strspn(line," \t\n\r"); return utils::trim(utils::trim_comment(line));
int stop = strcspn(line,"#") - 1;
while (line[stop] == ' ' || line[stop] == '\t'
|| line[stop] == '\n' || line[stop] == '\r') stop--;
line[stop+1] = '\0';
strcpy(keyword,&line[start]);
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
skip N lines of file skip N lines of file. Check if non-numeric content (e.g. keyword).
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Molecule::skip_lines(int n, char *line) void Molecule::skip_lines(int n, char *line, const std::string &section)
{ {
for (int i = 0; i < n; i++) readline(line); for (int i = 0; i < n; i++) {
readline(line);
if (utils::strmatch(utils::trim(utils::trim_comment(line)),"^[A-Za-z ]+$"))
error->one(FLERR,fmt::format("Unexpected line in molecule file while "
"skipping {} section:\n{}",section,line));
}
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------

5
src/molecule.h
View File

@ -161,10 +161,9 @@ class Molecule : protected Pointers {
void allocate(); void allocate();
void deallocate(); void deallocate();
void open(char *);
void readline(char *); void readline(char *);
void parse_keyword(int, char *, char *); std::string parse_keyword(int, char *);
void skip_lines(int, char *); void skip_lines(int, char *, const std::string &);
// void print(); // void print();
}; };

20
src/utils.cpp
View File

@ -854,7 +854,7 @@ std::vector<std::string> utils::split_words(const std::string &text)
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
bool utils::is_integer(const std::string &str) { bool utils::is_integer(const std::string &str) {
if (str.size() == 0) { if (str.empty()) {
return false; return false;
} }
@ -870,7 +870,7 @@ bool utils::is_integer(const std::string &str) {
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
bool utils::is_double(const std::string &str) { bool utils::is_double(const std::string &str) {
if (str.size() == 0) { if (str.empty()) {
return false; return false;
} }
@ -883,6 +883,22 @@ bool utils::is_double(const std::string &str) {
return true; return true;
} }
/* ----------------------------------------------------------------------
Return whether string is a valid ID string
------------------------------------------------------------------------- */
bool utils::is_id(const std::string &str) {
if (str.empty()) {
return false;
}
for (auto c : str) {
if (isalnum(c) || (c == '_')) continue;
return false;
}
return true;
}
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
strip off leading part of path, return just the filename strip off leading part of path, return just the filename
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */

8
src/utils.h
View File

@ -335,6 +335,14 @@ namespace LAMMPS_NS {
bool is_double(const std::string &str); bool is_double(const std::string &str);
/** Check if string is a valid ID
* ID strings may contain only letters, numbers, and underscores.
*
* \param str string that should be checked
* \return true, if string contains valid id, false otherwise */
bool is_id(const std::string &str);
/** Try to detect pathname from FILE pointer. /** Try to detect pathname from FILE pointer.
* *
* Currently only supported on Linux, otherwise will report "(unknown)". * Currently only supported on Linux, otherwise will report "(unknown)".

4
unittest/formats/CMakeLists.txt
View File

@ -7,6 +7,10 @@ add_executable(test_image_flags test_image_flags.cpp)
target_link_libraries(test_image_flags PRIVATE lammps GTest::GMock GTest::GTest) target_link_libraries(test_image_flags PRIVATE lammps GTest::GMock GTest::GTest)
add_test(NAME ImageFlags COMMAND test_image_flags WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}) add_test(NAME ImageFlags COMMAND test_image_flags WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR})
add_executable(test_molecule_file test_molecule_file.cpp)
target_link_libraries(test_molecule_file PRIVATE lammps GTest::GMock GTest::GTest)
add_test(NAME MoleculeFile COMMAND test_molecule_file WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR})
add_executable(test_pair_unit_convert test_pair_unit_convert.cpp) add_executable(test_pair_unit_convert test_pair_unit_convert.cpp)
target_link_libraries(test_pair_unit_convert PRIVATE lammps GTest::GMock GTest::GTest) target_link_libraries(test_pair_unit_convert PRIVATE lammps GTest::GMock GTest::GTest)
add_test(NAME PairUnitConvert COMMAND test_pair_unit_convert WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}) add_test(NAME PairUnitConvert COMMAND test_pair_unit_convert WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR})

313
unittest/formats/test_molecule_file.cpp Normal file
View File

@ -0,0 +1,313 @@
/* ----------------------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
https://lammps.sandia.gov/, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#include "atom.h"
#include "info.h"
#include "input.h"
#include "lammps.h"
#include "molecule.h"
#include "utils.h"
#include "gmock/gmock.h"
#include "gtest/gtest.h"
#include <cstdio>
#include <mpi.h>
#include <string>
using namespace LAMMPS_NS;
using testing::MatchesRegex;
using testing::StrEq;
using utils::split_words;
#define test_name test_info_->name()
#if defined(OMPI_MAJOR_VERSION)
const bool have_openmpi = true;
#else
const bool have_openmpi = false;
#endif
#define TEST_FAILURE(errmsg, ...) \
if (Info::has_exceptions()) { \
::testing::internal::CaptureStdout(); \
ASSERT_ANY_THROW({__VA_ARGS__}); \
auto mesg = ::testing::internal::GetCapturedStdout(); \
ASSERT_THAT(mesg, MatchesRegex(errmsg)); \
} else { \
if (!have_openmpi) { \
::testing::internal::CaptureStdout(); \
ASSERT_DEATH({__VA_ARGS__}, ""); \
auto mesg = ::testing::internal::GetCapturedStdout(); \
ASSERT_THAT(mesg, MatchesRegex(errmsg)); \
} \
}
static void create_molecule_files()
{
// create molecule files
const char h2o_file[] = "# Water molecule. SPC/E model.\n\n3 atoms\n2 bonds\n1 angles\n\n"
"Coords\n\n1 1.12456 0.09298 1.27452\n"
"2 1.53683 0.75606 1.89928\n3 0.49482 0.56390 0.65678\n\n"
"Types\n\n1 1\n2 2\n3 2\n\n"
"Charges\n\n1 -0.8472\n2 0.4236\n3 0.4236\n\n"
"Bonds\n\n1 1 1 2\n2 1 1 3\n\n"
"Angles\n\n1 1 2 1 3\n\n"
"Shake Flags\n\n1 1\n2 1\n3 1\n\n"
"Shake Atoms\n\n1 1 2 3\n2 1 2 3\n3 1 2 3\n\n"
"Shake Bond Types\n\n1 1 1 1\n2 1 1 1\n3 1 1 1\n\n"
"Special Bond Counts\n\n1 2 0 0\n2 1 1 0\n3 1 1 0\n\n"
"Special Bonds\n\n1 2 3\n2 1 3\n3 1 2\n\n";
const char co2_file[] = "# CO2 molecule file. TraPPE model.\n\n"
"3 atoms\n2 bonds\n1 angles\n\n"
"Coords\n\n1 0.0 0.0 0.0\n2 -1.16 0.0 0.0\n3 1.16 0.0 0.0\n\n"
"Types\n\n1 1\n2 2\n3 2\n\n"
"Charges\n\n1 0.7\n2 -0.35\n3 -0.35\n\n"
"Bonds\n\n1 1 1 2\n2 1 1 3\n\n"
"Angles\n\n1 1 2 1 3\n\n"
"Special Bond Counts\n\n1 2 0 0\n2 1 1 0\n3 1 1 0\n\n"
"Special Bonds\n\n1 2 3\n2 1 3\n3 1 2\n\n";
FILE *fp = fopen("tmp.h2o.mol", "w");
if (fp) {
fputs(h2o_file, fp);
fclose(fp);
}
rename("tmp.h2o.mol", "h2o.mol");
fp = fopen("tmp.co2.mol", "w");
if (fp) {
fputs(co2_file, fp);
fclose(fp);
}
rename("tmp.co2.mol", "co2.mol");
}
// whether to print verbose output (i.e. not capturing LAMMPS screen output).
bool verbose = false;
class MoleculeFileTest : public ::testing::Test {
protected:
LAMMPS *lmp;
void SetUp() override
{
const char *args[] = {"MoleculeFileTest", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
int argc = sizeof(args) / sizeof(char *);
if (!verbose) ::testing::internal::CaptureStdout();
lmp = new LAMMPS(argc, argv, MPI_COMM_WORLD);
create_molecule_files();
if (!verbose) ::testing::internal::GetCapturedStdout();
ASSERT_NE(lmp, nullptr);
}
void TearDown() override
{
if (!verbose) ::testing::internal::CaptureStdout();
delete lmp;
if (!verbose) ::testing::internal::GetCapturedStdout();
remove("h2o.mol");
remove("co2.mol");
}
void run_mol_cmd(const std::string &name, const std::string &args, const std::string &content)
{
std::string file = name + ".mol";
FILE *fp = fopen(file.c_str(), "w");
fputs(content.c_str(), fp);
fclose(fp);
lmp->input->one(fmt::format("molecule {} {} {}", name, file, args));
remove(file.c_str());
}
};
TEST_F(MoleculeFileTest, nofile)
{
TEST_FAILURE(".*Cannot open molecule file nofile.mol.*",
lmp->input->one("molecule 1 nofile.mol"););
}
TEST_F(MoleculeFileTest, badid)
{
TEST_FAILURE(".*Molecule template ID must have only "
"alphanumeric or underscore characters.*",
lmp->input->one("molecule @mol nofile.mol"););
}
TEST_F(MoleculeFileTest, badargs)
{
TEST_FAILURE(".*Illegal molecule command.*",
run_mol_cmd(test_name, "offset 1 2 3 4",
"Comment\n1 atoms\n\n Coords\n\n 1 0.0 0.0 0.0\n"););
TEST_FAILURE(
".*Illegal molecule command.*",
run_mol_cmd(test_name, "toff", "Comment\n1 atoms\n\n Coords\n\n 1 0.0 0.0 0.0\n"););
TEST_FAILURE(
".*Illegal molecule command.*",
run_mol_cmd(test_name, "boff", "Comment\n1 atoms\n\n Coords\n\n 1 0.0 0.0 0.0\n"););
TEST_FAILURE(
".*Illegal molecule command.*",
run_mol_cmd(test_name, "aoff", "Comment\n1 atoms\n\n Coords\n\n 1 0.0 0.0 0.0\n"););
TEST_FAILURE(
".*Illegal molecule command.*",
run_mol_cmd(test_name, "doff", "Comment\n1 atoms\n\n Coords\n\n 1 0.0 0.0 0.0\n"););
TEST_FAILURE(
".*Illegal molecule command.*",
run_mol_cmd(test_name, "ioff", "Comment\n1 atoms\n\n Coords\n\n 1 0.0 0.0 0.0\n"););
TEST_FAILURE(
".*Illegal molecule command.*",
run_mol_cmd(test_name, "scale", "Comment\n1 atoms\n\n Coords\n\n 1 0.0 0.0 0.0\n"););
remove("badargs.mol");
}
TEST_F(MoleculeFileTest, noatom)
{
TEST_FAILURE(".*ERROR: No atoms or invalid atom count in molecule file.*",
run_mol_cmd(test_name, "",
"Comment\n0 atoms\n1 bonds\n\n"
" Coords\n\nBonds\n\n 1 1 2\n"););
remove("noatom.mol");
}
TEST_F(MoleculeFileTest, empty)
{
TEST_FAILURE(".*ERROR: Unexpected end of molecule file.*",
run_mol_cmd(test_name, "", "Comment\n\n"););
remove("empty.mol");
}
TEST_F(MoleculeFileTest, nospecial)
{
TEST_FAILURE(".*ERROR: Cannot auto-generate special bonds before simulation box is defined.*",
run_mol_cmd(test_name, "",
"Comment\n3 atoms\n\n2 bonds\n\n"
" Coords\n\n 1 1.0 1.0 1.0\n 2 1.0 1.0 0.0\n 3 1.0 0.0 1.0\n"
" Bonds\n\n 1 1 1 2\n 2 1 1 3\n"););
remove("nospecial.mol");
}
TEST_F(MoleculeFileTest, minimal)
{
::testing::internal::CaptureStdout();
run_mol_cmd(test_name, "", "Comment\n1 atoms\n\n Coords\n\n 1 0.0 0.0 0.0\n");
auto output = ::testing::internal::GetCapturedStdout();
if (verbose) std::cout << output;
ASSERT_THAT(output, MatchesRegex(".*Read molecule template.*1 molecules.*1 atoms.*0 bonds.*"));
}
TEST_F(MoleculeFileTest, twomols)
{
::testing::internal::CaptureStdout();
run_mol_cmd(test_name, "",
"Comment\n2 atoms\n\n"
" Coords\n\n 1 0.0 0.0 0.0\n 2 0.0 0.0 1.0\n"
" Molecules\n\n 1 1\n 2 2\n\n Types\n\n 1 1\n 2 2\n\n");
auto output = ::testing::internal::GetCapturedStdout();
if (verbose) std::cout << output;
ASSERT_THAT(output, MatchesRegex(".*Read molecule template.*2 molecules.*2 atoms "
"with max type 2.*0 bonds.*"));
}
TEST_F(MoleculeFileTest, twofiles)
{
::testing::internal::CaptureStdout();
lmp->input->one("molecule twomols h2o.mol co2.mol offset 2 1 1 0 0");
auto output = ::testing::internal::GetCapturedStdout();
if (verbose) std::cout << output;
ASSERT_THAT(output, MatchesRegex(".*Read molecule template twomols:.*1 molecules.*3 atoms "
"with max type 2.*2 bonds with max type 1.*"
"1 angles with max type 1.*0 dihedrals.*"
".*Read molecule template twomols:.*1 molecules.*3 atoms "
"with max type 4.*2 bonds with max type 2.*"
"1 angles with max type 2.*0 dihedrals.*"));
}
TEST_F(MoleculeFileTest, bonds)
{
::testing::internal::CaptureStdout();
lmp->input->one("atom_style bond");
lmp->input->one("region box block 0 1 0 1 0 1");
lmp->input->one("create_box 2 box bond/types 2 extra/bond/per/atom 2 "
"extra/special/per/atom 4");
run_mol_cmd(test_name, "",
"Comment\n"
"4 atoms\n"
"2 bonds\n\n"
" Coords\n\n"
" 1 1.0 1.0 1.0\n"
" 2 1.0 1.0 0.0\n"
" 3 1.0 0.0 1.0\n"
" 4 1.0 0.0 0.0\n"
" Types\n\n"
" 1 1\n"
" 2 1\n"
" 3 2\n"
" 4 2\n\n"
" Bonds\n\n"
" 1 1 1 2\n"
" 2 2 1 3\n\n");
auto output = ::testing::internal::GetCapturedStdout();
if (verbose) std::cout << output;
ASSERT_THAT(output, MatchesRegex(".*Read molecule template.*1 molecules.*4 atoms.*type.*2.*"
"2 bonds.*type.*2.*0 angles.*"));
::testing::internal::CaptureStdout();
lmp->input->one("mass * 2.0");
lmp->input->one("create_atoms 0 single 0.5 0.5 0.5 mol bonds 67235");
output = ::testing::internal::GetCapturedStdout();
if (verbose) std::cout << output;
ASSERT_THAT(output, MatchesRegex(".*Created 4 atoms.*"));
::testing::internal::CaptureStdout();
Molecule *mol = lmp->atom->molecules[0];
ASSERT_EQ(mol->natoms, 4);
ASSERT_EQ(lmp->atom->natoms, 4);
mol->compute_mass();
mol->compute_com();
ASSERT_DOUBLE_EQ(mol->masstotal, 8.0);
EXPECT_DOUBLE_EQ(mol->com[0], 1.0);
EXPECT_DOUBLE_EQ(mol->com[1], 0.5);
EXPECT_DOUBLE_EQ(mol->com[2], 0.5);
EXPECT_DOUBLE_EQ(mol->maxextent, sqrt(2));
EXPECT_EQ(mol->comatom, 1);
output = ::testing::internal::GetCapturedStdout();
if (verbose) std::cout << output;
}
int main(int argc, char **argv)
{
MPI_Init(&argc, &argv);
::testing::InitGoogleMock(&argc, argv);
if (have_openmpi && !LAMMPS_NS::Info::has_exceptions())
std::cout << "Warning: using OpenMPI without exceptions. "
"Death tests will be skipped\n";
// handle arguments passed via environment variable
if (const char *var = getenv("TEST_ARGS")) {
std::vector<std::string> env = split_words(var);
for (auto arg : env) {
if (arg == "-v") {
verbose = true;
}
}
}
if ((argc > 1) && (strcmp(argv[1], "-v") == 0)) verbose = true;
int rv = RUN_ALL_TESTS();
MPI_Finalize();
return rv;
}

55
unittest/utils/test_utils.cpp
View File

@ -287,6 +287,61 @@ TEST(Utils, signed_double_and_d_exponential)
ASSERT_FALSE(utils::is_double("-10D-22")); ASSERT_FALSE(utils::is_double("-10D-22"));
} }
TEST(Utils, valid_id1)
{
ASSERT_TRUE(utils::is_id("abc"));
}
TEST(Utils, valid_id2)
{
ASSERT_TRUE(utils::is_id("123"));
}
TEST(Utils, valid_id3)
{
ASSERT_TRUE(utils::is_id("abc123"));
}
TEST(Utils, valid_id4)
{
ASSERT_TRUE(utils::is_id("abc_123"));
}
TEST(Utils, valid_id5)
{
ASSERT_TRUE(utils::is_id("123_abc"));
}
TEST(Utils, valid_id6)
{
ASSERT_TRUE(utils::is_id("_123"));
}
TEST(Utils, valid_id7)
{
ASSERT_TRUE(utils::is_id("___"));
}
TEST(Utils, invalid_id1)
{
ASSERT_FALSE(utils::is_id("+abc"));
}
TEST(Utils, invalid_id2)
{
ASSERT_FALSE(utils::is_id("a[1]"));
}
TEST(Utils, invalid_id3)
{
ASSERT_FALSE(utils::is_id("b(c)"));
}
TEST(Utils, invalid_id4)
{
ASSERT_FALSE(utils::is_id("a$12"));
}
TEST(Utils, strmatch_beg) TEST(Utils, strmatch_beg)
{ {
ASSERT_TRUE(utils::strmatch("rigid/small/omp", "^rigid")); ASSERT_TRUE(utils::strmatch("rigid/small/omp", "^rigid"));