sync with current develop

This commit is contained in:
Steve Plimpton
2022-04-18 17:29:23 -06:00
4880 changed files with 381040 additions and 501360 deletions

View File

@ -29,21 +29,18 @@
#include "modify.h"
#include "molecule.h"
#include "neighbor.h"
#include "tokenizer.h"
#include "update.h"
#include "variable.h"
#include "library.h"
#include <algorithm>
#include <cmath>
#include <cstring>
#ifdef LMP_INTEL
#include "neigh_request.h"
#endif
#ifdef LMP_GPU
#include "fix_gpu.h"
#include <cmath>
#endif
using namespace LAMMPS_NS;
@ -53,6 +50,15 @@ using namespace MathConst;
#define DELTA_PERATOM 64
#define EPSILON 1.0e-6
/* ----------------------------------------------------------------------
one instance per AtomVec style in style_atom.h
------------------------------------------------------------------------- */
template <typename T> static AtomVec *avec_creator(LAMMPS *lmp)
{
return new T(lmp);
}
/* ---------------------------------------------------------------------- */
/** \class LAMMPS_NS::Atom
@ -122,6 +128,7 @@ Atom::Atom(LAMMPS *lmp) : Pointers(lmp)
omega = angmom = torque = nullptr;
radius = rmass = nullptr;
ellipsoid = line = tri = body = nullptr;
quat = nullptr;
// molecular systems
@ -266,6 +273,7 @@ Atom::Atom(LAMMPS *lmp) : Pointers(lmp)
map_hash = nullptr;
unique_tags = nullptr;
reset_image_flag[0] = reset_image_flag[1] = reset_image_flag[2] = false;
atom_style = nullptr;
avec = nullptr;
@ -275,7 +283,7 @@ Atom::Atom(LAMMPS *lmp) : Pointers(lmp)
#define ATOM_CLASS
#define AtomStyle(key,Class) \
(*avec_map)[#key] = &avec_creator<Class>;
#include "style_atom.h"
#include "style_atom.h" // IWYU pragma: keep
#undef AtomStyle
#undef ATOM_CLASS
}
@ -357,7 +365,7 @@ Atom::~Atom()
// delete mapping data structures
map_delete();
Atom::map_delete();
delete unique_tags;
}
@ -425,6 +433,10 @@ void Atom::peratom_create()
add_peratom("tri",&tri,INT,0);
add_peratom("body",&body,INT,0);
// BPM package
add_peratom("quat",&quat,DOUBLE,4);
// MOLECULE package
add_peratom("molecule",&molecule,tagintsize,0);
@ -650,6 +662,7 @@ void Atom::set_atomflag_defaults()
// identical list as 2nd customization in atom.h
sphere_flag = ellipsoid_flag = line_flag = tri_flag = body_flag = 0;
quat_flag = 0;
peri_flag = electron_flag = 0;
wavepacket_flag = sph_flag = 0;
molecule_flag = molindex_flag = molatom_flag = 0;
@ -677,7 +690,7 @@ void Atom::set_atomflag_defaults()
void Atom::create_avec(const std::string &style, int narg, char **arg, int trysuffix)
{
delete [] atom_style;
delete[] atom_style;
if (avec) delete avec;
atom_style = nullptr;
avec = nullptr;
@ -701,11 +714,9 @@ void Atom::create_avec(const std::string &style, int narg, char **arg, int trysu
std::string estyle = style + "/";
if (sflag == 1) estyle += lmp->suffix;
else estyle += lmp->suffix2;
atom_style = new char[estyle.size()+1];
strcpy(atom_style,estyle.c_str());
atom_style = utils::strdup(estyle);
} else {
atom_style = new char[style.size()+1];
strcpy(atom_style,style.c_str());
atom_style = utils::strdup(style);
}
// if molecular system:
@ -755,16 +766,6 @@ AtomVec *Atom::new_avec(const std::string &style, int trysuffix, int &sflag)
return nullptr;
}
/* ----------------------------------------------------------------------
one instance per AtomVec style in style_atom.h
------------------------------------------------------------------------- */
template <typename T>
AtomVec *Atom::avec_creator(LAMMPS *lmp)
{
return new T(lmp);
}
/* ---------------------------------------------------------------------- */
void Atom::init()
@ -813,7 +814,7 @@ AtomVec *Atom::style_match(const char *style)
{
if (strcmp(atom_style,style) == 0) return avec;
else if (strcmp(atom_style,"hybrid") == 0) {
auto avec_hybrid = (AtomVecHybrid *) avec;
auto avec_hybrid = dynamic_cast<AtomVecHybrid *>( avec);
for (int i = 0; i < avec_hybrid->nstyles; i++)
if (strcmp(avec_hybrid->keywords[i],style) == 0)
return avec_hybrid->styles[i];
@ -836,17 +837,13 @@ void Atom::modify_params(int narg, char **arg)
if (strcmp(arg[iarg],"id") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal atom_modify command");
if (domain->box_exist)
error->all(FLERR,
"Atom_modify id command after simulation box is defined");
if (strcmp(arg[iarg+1],"yes") == 0) tag_enable = 1;
else if (strcmp(arg[iarg+1],"no") == 0) tag_enable = 0;
else error->all(FLERR,"Illegal atom_modify command");
error->all(FLERR,"Atom_modify id command after simulation box is defined");
tag_enable = utils::logical(FLERR,arg[iarg+1],false,lmp);
iarg += 2;
} else if (strcmp(arg[iarg],"map") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal atom_modify command");
if (domain->box_exist)
error->all(FLERR,
"Atom_modify map command after simulation box is defined");
error->all(FLERR,"Atom_modify map command after simulation box is defined");
if (strcmp(arg[iarg+1],"array") == 0) map_user = 1;
else if (strcmp(arg[iarg+1],"hash") == 0) map_user = 2;
else if (strcmp(arg[iarg+1],"yes") == 0) map_user = 3;
@ -902,10 +899,10 @@ void Atom::tag_check()
MPI_Allreduce(&min,&minall,1,MPI_LMP_TAGINT,MPI_MIN,world);
MPI_Allreduce(&max,&maxall,1,MPI_LMP_TAGINT,MPI_MAX,world);
if (minall < 0) error->all(FLERR,"One or more Atom IDs is negative");
if (maxall >= MAXTAGINT) error->all(FLERR,"One or more atom IDs is too big");
if (minall < 0) error->all(FLERR,"One or more Atom IDs are negative");
if (maxall >= MAXTAGINT) error->all(FLERR,"One or more atom IDs are too big");
if (maxall > 0 && minall == 0)
error->all(FLERR,"One or more atom IDs is zero");
error->all(FLERR,"One or more atom IDs are zero");
if (maxall > 0 && tag_enable == 0)
error->all(FLERR,"Non-zero atom IDs with atom_modify id = no");
if (maxall == 0 && natoms && tag_enable)
@ -1068,12 +1065,13 @@ void Atom::deallocate_topology()
void Atom::data_atoms(int n, char *buf, tagint id_offset, tagint mol_offset,
int type_offset, int shiftflag, double *shift)
{
int m,xptr,iptr;
int xptr,iptr;
imageint imagedata;
double xdata[3],lamda[3];
double *coord;
char *next;
// use the first line to detect and validate the number of words/tokens per line
next = strchr(buf,'\n');
*next = '\0';
int nwords = utils::trim_and_count_words(buf);
@ -1082,8 +1080,6 @@ void Atom::data_atoms(int n, char *buf, tagint id_offset, tagint mol_offset,
if (nwords != avec->size_data_atom && nwords != avec->size_data_atom + 3)
error->all(FLERR,"Incorrect atom format in data file");
char **values = new char*[nwords];
// set bounds for my proc
// if periodic and I am lo/hi proc, adjust bounds by EPSILON
// insures all data atoms will be owned even with round-off
@ -1152,18 +1148,12 @@ void Atom::data_atoms(int n, char *buf, tagint id_offset, tagint mol_offset,
// remap atom into simulation box
// if atom is in my sub-domain, unpack its values
int flagx = 0, flagy = 0, flagz = 0;
for (int i = 0; i < n; i++) {
next = strchr(buf,'\n');
for (m = 0; m < nwords; m++) {
buf += strspn(buf," \t\n\r\f");
buf[strcspn(buf," \t\n\r\f")] = '\0';
if (strlen(buf) == 0)
error->all(FLERR,"Incorrect atom format in data file");
values[m] = buf;
buf += strlen(buf)+1;
}
*next = '\0';
auto values = Tokenizer(utils::trim_comment(buf)).as_vector();
if ((int)values.size() != nwords)
error->all(FLERR, "Incorrect atom format in data file: {}", utils::trim(buf));
int imx = 0, imy = 0, imz = 0;
if (imageflag) {
@ -1172,9 +1162,9 @@ void Atom::data_atoms(int n, char *buf, tagint id_offset, tagint mol_offset,
imz = utils::inumeric(FLERR,values[iptr+2],false,lmp);
if ((domain->dimension == 2) && (imz != 0))
error->all(FLERR,"Z-direction image flag must be 0 for 2d-systems");
if ((!domain->xperiodic) && (imx != 0)) { flagx = 1; imx = 0; }
if ((!domain->yperiodic) && (imy != 0)) { flagy = 1; imy = 0; }
if ((!domain->zperiodic) && (imz != 0)) { flagz = 1; imz = 0; }
if ((!domain->xperiodic) && (imx != 0)) { reset_image_flag[0] = true; imx = 0; }
if ((!domain->yperiodic) && (imy != 0)) { reset_image_flag[1] = true; imy = 0; }
if ((!domain->zperiodic) && (imz != 0)) { reset_image_flag[2] = true; imz = 0; }
}
imagedata = ((imageint) (imx + IMGMAX) & IMGMASK) |
(((imageint) (imy + IMGMAX) & IMGMASK) << IMGBITS) |
@ -1210,24 +1200,6 @@ void Atom::data_atoms(int n, char *buf, tagint id_offset, tagint mol_offset,
buf = next + 1;
}
// warn if reading data with non-zero image flags for non-periodic boundaries.
// we may want to turn this into an error at some point, since this essentially
// creates invalid position information that works by accident most of the time.
if (comm->me == 0) {
if (flagx)
error->warning(FLERR,"Non-zero imageflag(s) in x direction for "
"non-periodic boundary reset to zero");
if (flagy)
error->warning(FLERR,"Non-zero imageflag(s) in y direction for "
"non-periodic boundary reset to zero");
if (flagz)
error->warning(FLERR,"Non-zero imageflag(s) in z direction for "
"non-periodic boundary reset to zero");
}
delete [] values;
}
/* ----------------------------------------------------------------------
@ -1238,8 +1210,7 @@ void Atom::data_atoms(int n, char *buf, tagint id_offset, tagint mol_offset,
void Atom::data_vels(int n, char *buf, tagint id_offset)
{
int j,m;
tagint tagdata;
int m;
char *next;
next = strchr(buf,'\n');
@ -1250,31 +1221,24 @@ void Atom::data_vels(int n, char *buf, tagint id_offset)
if (nwords != avec->size_data_vel)
error->all(FLERR,"Incorrect velocity format in data file");
char **values = new char*[nwords];
// loop over lines of atom velocities
// tokenize the line into values
// if I own atom tag, unpack its values
for (int i = 0; i < n; i++) {
next = strchr(buf,'\n');
*next = '\0';
auto values = Tokenizer(utils::trim_comment(buf)).as_vector();
if ((int)values.size() != nwords)
error->all(FLERR, "Incorrect atom format in data file: {}", utils::trim(buf));
for (j = 0; j < nwords; j++) {
buf += strspn(buf," \t\n\r\f");
buf[strcspn(buf," \t\n\r\f")] = '\0';
values[j] = buf;
buf += strlen(buf)+1;
}
tagdata = ATOTAGINT(values[0]) + id_offset;
tagint tagdata = utils::tnumeric(FLERR,values[0],false,lmp) + id_offset;
if (tagdata <= 0 || tagdata > map_tag_max)
error->one(FLERR,"Invalid atom ID in Velocities section of data file");
if ((m = map(tagdata)) >= 0) avec->data_vel(m,&values[1]);
if ((m = map(tagdata)) >= 0) avec->data_vel(m,values);
buf = next + 1;
}
delete [] values;
}
/* ----------------------------------------------------------------------
@ -1287,18 +1251,25 @@ void Atom::data_vels(int n, char *buf, tagint id_offset)
void Atom::data_bonds(int n, char *buf, int *count, tagint id_offset,
int type_offset)
{
int m,tmp,itype,rv;
int m,itype;
tagint atom1,atom2;
char *next;
int newton_bond = force->newton_bond;
auto location = "Bonds section of data file";
for (int i = 0; i < n; i++) {
next = strchr(buf,'\n');
*next = '\0';
rv = sscanf(buf,"%d %d " TAGINT_FORMAT " " TAGINT_FORMAT,
&tmp,&itype,&atom1,&atom2);
if (rv != 4)
error->one(FLERR,"Incorrect format of Bonds section in data file");
try {
ValueTokenizer values(utils::trim_comment(buf));
values.next_int();
itype = values.next_int();
atom1 = values.next_tagint();
atom2 = values.next_tagint();
if (values.has_next()) throw TokenizerException("Too many tokens","");
} catch (TokenizerException &e) {
error->one(FLERR,"{} in {}: {}", e.what(), location, utils::trim(buf));
}
if (id_offset) {
atom1 += id_offset;
atom2 += id_offset;
@ -1307,9 +1278,9 @@ void Atom::data_bonds(int n, char *buf, int *count, tagint id_offset,
if ((atom1 <= 0) || (atom1 > map_tag_max) ||
(atom2 <= 0) || (atom2 > map_tag_max) || (atom1 == atom2))
error->one(FLERR,"Invalid atom ID in Bonds section of data file");
error->one(FLERR,"Invalid atom ID in {}: {}", location, utils::trim(buf));
if (itype <= 0 || itype > nbondtypes)
error->one(FLERR,"Invalid bond type in Bonds section of data file");
error->one(FLERR,"Invalid bond type in {}: {}", location, utils::trim(buf));
if ((m = map(atom1)) >= 0) {
if (count) count[m]++;
else {
@ -1344,18 +1315,26 @@ void Atom::data_bonds(int n, char *buf, int *count, tagint id_offset,
void Atom::data_angles(int n, char *buf, int *count, tagint id_offset,
int type_offset)
{
int m,tmp,itype,rv;
int m,itype;
tagint atom1,atom2,atom3;
char *next;
int newton_bond = force->newton_bond;
auto location = "Angles section of data file";
for (int i = 0; i < n; i++) {
next = strchr(buf,'\n');
*next = '\0';
rv = sscanf(buf,"%d %d " TAGINT_FORMAT " " TAGINT_FORMAT " " TAGINT_FORMAT,
&tmp,&itype,&atom1,&atom2,&atom3);
if (rv != 5)
error->one(FLERR,"Incorrect format of Angles section in data file");
try {
ValueTokenizer values(utils::trim_comment(buf));
values.next_int();
itype = values.next_int();
atom1 = values.next_tagint();
atom2 = values.next_tagint();
atom3 = values.next_tagint();
if (values.has_next()) throw TokenizerException("Too many tokens","");
} catch (TokenizerException &e) {
error->one(FLERR,"{} in {}: {}", e.what(), location, utils::trim(buf));
}
if (id_offset) {
atom1 += id_offset;
atom2 += id_offset;
@ -1367,9 +1346,9 @@ void Atom::data_angles(int n, char *buf, int *count, tagint id_offset,
(atom2 <= 0) || (atom2 > map_tag_max) ||
(atom3 <= 0) || (atom3 > map_tag_max) ||
(atom1 == atom2) || (atom1 == atom3) || (atom2 == atom3))
error->one(FLERR,"Invalid atom ID in Angles section of data file");
error->one(FLERR,"Invalid atom ID in {}: {}", location, utils::trim(buf));
if (itype <= 0 || itype > nangletypes)
error->one(FLERR,"Invalid angle type in Angles section of data file");
error->one(FLERR,"Invalid angle type in {}: {}", location, utils::trim(buf));
if ((m = map(atom2)) >= 0) {
if (count) count[m]++;
else {
@ -1416,19 +1395,27 @@ void Atom::data_angles(int n, char *buf, int *count, tagint id_offset,
void Atom::data_dihedrals(int n, char *buf, int *count, tagint id_offset,
int type_offset)
{
int m,tmp,itype,rv;
int m,itype;
tagint atom1,atom2,atom3,atom4;
char *next;
int newton_bond = force->newton_bond;
auto location = "Dihedrals section of data file";
for (int i = 0; i < n; i++) {
next = strchr(buf,'\n');
*next = '\0';
rv = sscanf(buf,"%d %d " TAGINT_FORMAT " " TAGINT_FORMAT
" " TAGINT_FORMAT " " TAGINT_FORMAT,
&tmp,&itype,&atom1,&atom2,&atom3,&atom4);
if (rv != 6)
error->one(FLERR,"Incorrect format of Dihedrals section in data file");
try {
ValueTokenizer values(utils::trim_comment(buf));
values.next_int();
itype = values.next_int();
atom1 = values.next_tagint();
atom2 = values.next_tagint();
atom3 = values.next_tagint();
atom4 = values.next_tagint();
if (values.has_next()) throw TokenizerException("Too many tokens","");
} catch (TokenizerException &e) {
error->one(FLERR,"{} in {}: {}", e.what(), location, utils::trim(buf));
}
if (id_offset) {
atom1 += id_offset;
atom2 += id_offset;
@ -1443,10 +1430,9 @@ void Atom::data_dihedrals(int n, char *buf, int *count, tagint id_offset,
(atom4 <= 0) || (atom4 > map_tag_max) ||
(atom1 == atom2) || (atom1 == atom3) || (atom1 == atom4) ||
(atom2 == atom3) || (atom2 == atom4) || (atom3 == atom4))
error->one(FLERR,"Invalid atom ID in Dihedrals section of data file");
error->one(FLERR, "Invalid atom ID in {}: {}", location, utils::trim(buf));
if (itype <= 0 || itype > ndihedraltypes)
error->one(FLERR,
"Invalid dihedral type in Dihedrals section of data file");
error->one(FLERR, "Invalid dihedral type in {}: {}", location, utils::trim(buf));
if ((m = map(atom2)) >= 0) {
if (count) count[m]++;
else {
@ -1507,19 +1493,27 @@ void Atom::data_dihedrals(int n, char *buf, int *count, tagint id_offset,
void Atom::data_impropers(int n, char *buf, int *count, tagint id_offset,
int type_offset)
{
int m,tmp,itype,rv;
int m,itype;
tagint atom1,atom2,atom3,atom4;
char *next;
int newton_bond = force->newton_bond;
auto location = "Impropers section of data file";
for (int i = 0; i < n; i++) {
next = strchr(buf,'\n');
*next = '\0';
rv = sscanf(buf,"%d %d "
TAGINT_FORMAT " " TAGINT_FORMAT " " TAGINT_FORMAT " " TAGINT_FORMAT,
&tmp,&itype,&atom1,&atom2,&atom3,&atom4);
if (rv != 6)
error->one(FLERR,"Incorrect format of Impropers section in data file");
try {
ValueTokenizer values(utils::trim_comment(buf));
values.next_int();
itype = values.next_int();
atom1 = values.next_tagint();
atom2 = values.next_tagint();
atom3 = values.next_tagint();
atom4 = values.next_tagint();
if (values.has_next()) throw TokenizerException("Too many tokens","");
} catch (TokenizerException &e) {
error->one(FLERR,"{} in {}: {}", e.what(), location, utils::trim(buf));
}
if (id_offset) {
atom1 += id_offset;
atom2 += id_offset;
@ -1534,10 +1528,9 @@ void Atom::data_impropers(int n, char *buf, int *count, tagint id_offset,
(atom4 <= 0) || (atom4 > map_tag_max) ||
(atom1 == atom2) || (atom1 == atom3) || (atom1 == atom4) ||
(atom2 == atom3) || (atom2 == atom4) || (atom3 == atom4))
error->one(FLERR,"Invalid atom ID in Impropers section of data file");
error->one(FLERR, "Invalid atom ID in {}: {}", location, utils::trim(buf));
if (itype <= 0 || itype > nimpropertypes)
error->one(FLERR,
"Invalid improper type in Impropers section of data file");
error->one(FLERR, "Invalid improper type in {}: {}", location, utils::trim(buf));
if ((m = map(atom2)) >= 0) {
if (count) count[m]++;
else {
@ -1596,7 +1589,7 @@ void Atom::data_impropers(int n, char *buf, int *count, tagint id_offset,
void Atom::data_bonus(int n, char *buf, AtomVec *avec_bonus, tagint id_offset)
{
int j,m,tagdata;
int m;
char *next;
next = strchr(buf,'\n');
@ -1607,35 +1600,28 @@ void Atom::data_bonus(int n, char *buf, AtomVec *avec_bonus, tagint id_offset)
if (nwords != avec_bonus->size_data_bonus)
error->all(FLERR,"Incorrect bonus data format in data file");
char **values = new char*[nwords];
// loop over lines of bonus atom data
// tokenize the line into values
// if I own atom tag, unpack its values
for (int i = 0; i < n; i++) {
next = strchr(buf,'\n');
*next = '\0';
auto values = Tokenizer(utils::trim_comment(buf)).as_vector();
if ((int)values.size() != nwords)
error->all(FLERR, "Incorrect atom format in data file: {}", utils::trim(buf));
for (j = 0; j < nwords; j++) {
buf += strspn(buf," \t\n\r\f");
buf[strcspn(buf," \t\n\r\f")] = '\0';
values[j] = buf;
buf += strlen(buf)+1;
}
tagdata = ATOTAGINT(values[0]) + id_offset;
tagint tagdata = utils::tnumeric(FLERR,values[0],false,lmp) + id_offset;
if (tagdata <= 0 || tagdata > map_tag_max)
error->one(FLERR,"Invalid atom ID in Bonus section of data file");
// ok to call child's data_atom_bonus() method thru parent avec_bonus,
// since data_bonus() was called with child ptr, and method is virtual
if ((m = map(tagdata)) >= 0) avec_bonus->data_atom_bonus(m,&values[1]);
if ((m = map(tagdata)) >= 0) avec_bonus->data_atom_bonus(m,values);
buf = next + 1;
}
delete [] values;
}
/* ----------------------------------------------------------------------
@ -1647,12 +1633,8 @@ void Atom::data_bonus(int n, char *buf, AtomVec *avec_bonus, tagint id_offset)
void Atom::data_bodies(int n, char *buf, AtomVec *avec_body, tagint id_offset)
{
int j,m,nvalues,tagdata,ninteger,ndouble;
int maxint = 0;
int maxdouble = 0;
int *ivalues = nullptr;
double *dvalues = nullptr;
std::vector<int> ivalues;
std::vector<double> dvalues;
if (!unique_tags) unique_tags = new std::set<tagint>;
@ -1661,69 +1643,51 @@ void Atom::data_bodies(int n, char *buf, AtomVec *avec_body, tagint id_offset)
// else skip values
for (int i = 0; i < n; i++) {
buf += strspn(buf," \t\n\r\f");
buf[strcspn(buf," \t\n\r\f")] = '\0';
tagdata = utils::tnumeric(FLERR,buf,false,lmp) + id_offset;
buf += strlen(buf)+1;
char *next = strchr(buf,'\n');
*next = '\0';
if (tagdata <= 0 || tagdata > map_tag_max)
error->one(FLERR,"Invalid atom ID in Bodies section of data file");
auto values = Tokenizer(utils::trim_comment(buf)).as_vector();
tagint tagdata = utils::tnumeric(FLERR,values[0],false,lmp) + id_offset;
int ninteger = utils::inumeric(FLERR,values[1],false,lmp);
int ndouble = utils::inumeric(FLERR,values[2],false,lmp);
if (unique_tags->find(tagdata) == unique_tags->end())
unique_tags->insert(tagdata);
else
error->one(FLERR,"Duplicate atom ID in Bodies section of data file");
buf += strspn(buf," \t\n\r\f");
buf[strcspn(buf," \t\n\r\f")] = '\0';
ninteger = utils::inumeric(FLERR,buf,false,lmp);
buf += strlen(buf)+1;
buf = next + 1;
int m = map(tagdata);
if (m >= 0) {
ivalues.resize(ninteger);
dvalues.resize(ndouble);
buf += strspn(buf," \t\n\r\f");
buf[strcspn(buf," \t\n\r\f")] = '\0';
ndouble = utils::inumeric(FLERR,buf,false,lmp);
buf += strlen(buf)+1;
if ((m = map(tagdata)) >= 0) {
if (ninteger > maxint) {
delete [] ivalues;
maxint = ninteger;
ivalues = new int[maxint];
}
if (ndouble > maxdouble) {
delete [] dvalues;
maxdouble = ndouble;
dvalues = new double[maxdouble];
}
for (j = 0; j < ninteger; j++) {
for (int j = 0; j < ninteger; j++) {
buf += strspn(buf," \t\n\r\f");
buf[strcspn(buf," \t\n\r\f")] = '\0';
ivalues[j] = utils::inumeric(FLERR,buf,false,lmp);
buf += strlen(buf)+1;
}
for (j = 0; j < ndouble; j++) {
for (int j = 0; j < ndouble; j++) {
buf += strspn(buf," \t\n\r\f");
buf[strcspn(buf," \t\n\r\f")] = '\0';
dvalues[j] = utils::numeric(FLERR,buf,false,lmp);
buf += strlen(buf)+1;
}
avec_body->data_body(m,ninteger,ndouble,ivalues,dvalues);
avec_body->data_body(m,ninteger,ndouble,ivalues.data(),dvalues.data());
} else {
nvalues = ninteger + ndouble; // number of values to skip
for (j = 0; j < nvalues; j++) {
int nvalues = ninteger + ndouble; // number of values to skip
for (int j = 0; j < nvalues; j++) {
buf += strspn(buf," \t\n\r\f");
buf[strcspn(buf," \t\n\r\f")] = '\0';
buf += strlen(buf)+1;
}
}
buf += strspn(buf," \t\n\r\f");
}
delete [] ivalues;
delete [] dvalues;
}
/* ----------------------------------------------------------------------
@ -1735,8 +1699,7 @@ void Atom::data_bodies(int n, char *buf, AtomVec *avec_body, tagint id_offset)
void Atom::data_fix_compute_variable(int nprev, int nnew)
{
for (int m = 0; m < modify->nfix; m++) {
Fix *fix = modify->fix[m];
for (const auto &fix : modify->get_fix_list()) {
if (fix->create_attribute)
for (int i = nprev; i < nnew; i++)
fix->set_arrays(i);
@ -1779,17 +1742,20 @@ void Atom::set_mass(const char *file, int line, const char *str, int type_offset
int itype;
double mass_one;
int n = sscanf(str,"%d %lg",&itype,&mass_one);
if (n != 2) error->all(file,line,"Invalid mass line in data file");
itype += type_offset;
try {
ValueTokenizer values(utils::trim_comment(str));
itype = values.next_int() + type_offset;
mass_one = values.next_double();
if (values.has_next()) throw TokenizerException("Too many tokens", "");
if (itype < 1 || itype > ntypes)
error->all(file,line,"Invalid type for mass set");
if (itype < 1 || itype > ntypes) throw TokenizerException("Invalid atom type", "");
if (mass_one <= 0.0) throw TokenizerException("Invalid mass value", "");
} catch (TokenizerException &e) {
error->all(file,line,"{} in Masses section of data file: {}", e.what(), utils::trim(str));
}
mass[itype] = mass_one;
mass_setflag[itype] = 1;
if (mass[itype] <= 0.0) error->all(file,line,"Invalid mass value");
}
/* ----------------------------------------------------------------------
@ -1893,7 +1859,7 @@ int Atom::shape_consistency(int itype,
double one[3] = {-1.0, -1.0, -1.0};
double *shape;
auto avec_ellipsoid = (AtomVecEllipsoid *) style_match("ellipsoid");
auto avec_ellipsoid = dynamic_cast<AtomVecEllipsoid *>( style_match("ellipsoid"));
auto bonus = avec_ellipsoid->bonus;
int flag = 0;
@ -1939,10 +1905,9 @@ void Atom::add_molecule(int narg, char **arg)
int ifile = 1;
int index = 1;
while (1) {
while (true) {
molecules = (Molecule **)
memory->srealloc(molecules,(nmolecule+1)*sizeof(Molecule *),
"atom::molecules");
memory->srealloc(molecules,(nmolecule+1)*sizeof(Molecule *), "atom::molecules");
molecules[nmolecule] = new Molecule(lmp,narg,arg,index);
molecules[nmolecule]->nset = 0;
molecules[nmolecule-ifile+1]->nset++;
@ -1972,8 +1937,7 @@ int Atom::find_molecule(char *id)
called by fixes and commands that add molecules
------------------------------------------------------------------------- */
void Atom::add_molecule_atom(Molecule *onemol, int iatom,
int ilocal, tagint offset)
void Atom::add_molecule_atom(Molecule *onemol, int iatom, int ilocal, tagint offset)
{
if (onemol->qflag && q_flag) q[ilocal] = onemol->q[iatom];
if (onemol->radiusflag && radius_flag) radius[ilocal] = onemol->radius[iatom];
@ -1988,6 +1952,19 @@ void Atom::add_molecule_atom(Molecule *onemol, int iatom,
onemol->avec_body->set_quat(ilocal,onemol->quat_external);
}
// initialize custom per-atom properties to zero if present
for (int i = 0; i < nivector; ++i)
ivector[i][ilocal] = 0;
for (int i = 0; i < ndvector; ++i)
dvector[i][ilocal] = 0.0;
for (int i = 0; i < niarray; ++i)
for (int j = 0; j < icols[i]; ++j)
iarray[i][ilocal][j] = 0;
for (int i = 0; i < ndarray; ++i)
for (int j = 0; j < dcols[i]; ++j)
darray[i][ilocal][j] = 0.0;
if (molecular != Atom::MOLECULAR) return;
// add bond topology info
@ -2225,12 +2202,7 @@ void Atom::setup_sort_bins()
bininvz = nbinz / (bboxhi[2]-bboxlo[2]);
#ifdef LMP_INTEL
int intel_neigh = 0;
if (neighbor->nrequest) {
if (neighbor->requests[0]->intel) intel_neigh = 1;
} else if (neighbor->old_nrequest)
if (neighbor->old_requests[0]->intel) intel_neigh = 1;
if (intel_neigh && userbinsize == 0.0) {
if (neighbor->has_intel_request() && userbinsize == 0.0) {
if (neighbor->binsizeflag) bininv = 1.0/neighbor->binsize_user;
double nx_low = neighbor->bboxlo[0];
@ -2273,15 +2245,13 @@ void Atom::setup_sort_bins()
#ifdef LMP_GPU
if (userbinsize == 0.0) {
int ifix = modify->find_fix("package_gpu");
if (ifix >= 0) {
FixGPU *fix = dynamic_cast<FixGPU *>(modify->get_fix_by_id("package_gpu"));
if (fix) {
const double subx = domain->subhi[0] - domain->sublo[0];
const double suby = domain->subhi[1] - domain->sublo[1];
const double subz = domain->subhi[2] - domain->sublo[2];
FixGPU *fix = static_cast<FixGPU *>(modify->fix[ifix]);
binsize = fix->binsize(subx, suby, subz, atom->nlocal,
neighbor->cutneighmax);
binsize = fix->binsize(subx, suby, subz, atom->nlocal,neighbor->cutneighmax);
bininv = 1.0 / binsize;
nbinx = static_cast<int> (ceil(subx * bininv));
@ -2483,7 +2453,7 @@ This function is called, e.g. from :doc:`fix property/atom <fix_property_atom>`.
*/
int Atom::add_custom(const char *name, int flag, int cols)
{
int index;
int index = -1;
if ((flag == 0) && (cols == 0)) {
index = nivector;
@ -2523,7 +2493,8 @@ int Atom::add_custom(const char *name, int flag, int cols)
dcols = (int *) memory->srealloc(dcols,ndarray*sizeof(int),"atom:dcols");
dcols[index] = cols;
}
if (index < 0)
error->all(FLERR,"Invalid call to Atom::add_custom()");
return index;
}
@ -2544,27 +2515,27 @@ void Atom::remove_custom(int index, int flag, int cols)
{
if (flag == 0 && cols == 0) {
memory->destroy(ivector[index]);
ivector[index] = NULL;
ivector[index] = nullptr;
delete [] ivname[index];
ivname[index] = NULL;
ivname[index] = nullptr;
} else if (flag == 1 && cols == 0) {
memory->destroy(dvector[index]);
dvector[index] = NULL;
dvector[index] = nullptr;
delete [] dvname[index];
dvname[index] = NULL;
dvname[index] = nullptr;
} else if (flag == 0 && cols) {
memory->destroy(iarray[index]);
iarray[index] = NULL;
iarray[index] = nullptr;
delete [] ianame[index];
ianame[index] = NULL;
ianame[index] = nullptr;
} else if (flag == 1 && cols) {
memory->destroy(darray[index]);
darray[index] = NULL;
darray[index] = nullptr;
delete [] daname[index];
daname[index] = NULL;
daname[index] = nullptr;
}
}
@ -2672,6 +2643,10 @@ length of the data area, and a short description.
- int
- 1
- 1 if the particle is a body particle, 0 if not
* - quat
- double
- 4
- four quaternion components of the particles
* - i_name
- int
- 1
@ -2727,6 +2702,7 @@ void *Atom::extract(const char *name)
if (strcmp(name,"line") == 0) return (void *) line;
if (strcmp(name,"tri") == 0) return (void *) tri;
if (strcmp(name,"body") == 0) return (void *) body;
if (strcmp(name,"quat") == 0) return (void *) quat;
if (strcmp(name,"vfrac") == 0) return (void *) vfrac;
if (strcmp(name,"s0") == 0) return (void *) s0;
@ -2797,9 +2773,9 @@ void *Atom::extract(const char *name)
if (!array) index = find_custom(&name[2],flag,cols);
else index = find_custom(&name[3],flag,cols);
if (index < 0) return NULL;
if (which != flag) return NULL;
if ((!array && cols) || (array && !cols)) return NULL;
if (index < 0) return nullptr;
if (which != flag) return nullptr;
if ((!array && cols) || (array && !cols)) return nullptr;
if (!which && !array) return (void *) ivector[index];
if (which && !array) return (void *) dvector[index];
@ -2849,6 +2825,7 @@ int Atom::extract_datatype(const char *name)
if (strcmp(name,"line") == 0) return LAMMPS_INT;
if (strcmp(name,"tri") == 0) return LAMMPS_INT;
if (strcmp(name,"body") == 0) return LAMMPS_INT;
if (strcmp(name,"quat") == 0) return LAMMPS_DOUBLE_2D;
if (strcmp(name,"vfrac") == 0) return LAMMPS_DOUBLE;
if (strcmp(name,"s0") == 0) return LAMMPS_DOUBLE;