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

resolved merge conflicts with master

Browse Source
This commit is contained in:
Steve Plimpton
2021-08-17 11:55:33 -06:00
parent b2f96d92df dad9942bb8
commit e6af1eecb8
8612 changed files with 510252 additions and 341930 deletions
Download Patch File Download Diff File Expand all files Collapse all files

127
src/RIGID/fix_rigid.cpp
View File

@ -1,6 +1,7 @@
// clang-format off
/* ----------------------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
https://lammps.sandia.gov/, Sandia National Laboratories
https://www.lammps.org/, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
@ -13,31 +14,29 @@
#include "fix_rigid.h"
#include <cmath>
#include <cstring>
#include "math_extra.h"
#include "math_eigen.h"
#include "atom.h"
#include "atom_vec_ellipsoid.h"
#include "atom_vec_line.h"
#include "atom_vec_tri.h"
#include "domain.h"
#include "update.h"
#include "respa.h"
#include "modify.h"
#include "group.h"
#include "comm.h"
#include "random_mars.h"
#include "force.h"
#include "input.h"
#include "variable.h"
#include "math_const.h"
#include "memory.h"
#include "domain.h"
#include "error.h"
#include "force.h"
#include "group.h"
#include "input.h"
#include "math_const.h"
#include "math_eigen.h"
#include "math_extra.h"
#include "memory.h"
#include "modify.h"
#include "random_mars.h"
#include "respa.h"
#include "rigid_const.h"
#include "update.h"
#include "variable.h"
#include <cmath>
#include <cstring>
using namespace LAMMPS_NS;
using namespace FixConst;
@ -64,7 +63,7 @@ FixRigid::FixRigid(LAMMPS *lmp, int narg, char **arg) :
extscalar = 0;
time_integrate = 1;
rigid_flag = 1;
virial_flag = 1;
virial_global_flag = virial_peratom_flag = 1;
thermo_virial = 1;
create_attribute = 1;
dof_flag = 1;
@ -129,11 +128,15 @@ FixRigid::FixRigid(LAMMPS *lmp, int narg, char **arg) :
// determine whether atom-style variable or atom property is used
if (strstr(arg[4],"i_") == arg[4]) {
int flag,cols;
int custom_index = atom->find_custom(arg[4]+2,flag,cols);
if (custom_index < 0 || !flag || cols)
error->all(FLERR,"Fix rigid custom requires custom integer vector");
if (utils::strmatch(arg[4],"^i_")) {
int is_double,cols;
int custom_index = atom->find_custom(arg[4]+2,is_double,cols);
if (custom_index == -1)
error->all(FLERR,"Fix rigid custom requires "
"previously defined property/atom");
else if (is_double)
error->all(FLERR,"Fix rigid custom requires "
"integer-valued property/atom vector");
int minval = INT_MAX;
int *value = atom->ivector[custom_index];
for (i = 0; i < nlocal; i++)
@ -147,7 +150,7 @@ FixRigid::FixRigid(LAMMPS *lmp, int narg, char **arg) :
else
molecule[i] = 0;
} else if (strstr(arg[4],"v_") == arg[4]) {
} else if (utils::strmatch(arg[4],"^v_")) {
int ivariable = input->variable->find(arg[4]+2);
if (ivariable < 0)
error->all(FLERR,"Variable name for fix rigid custom does not exist");
@ -327,7 +330,7 @@ FixRigid::FixRigid(LAMMPS *lmp, int narg, char **arg) :
pstyle = ANISO;
dimension = domain->dimension;
for (int i = 0; i < 3; i++) {
for (i = 0; i < 3; i++) {
p_start[i] = p_stop[i] = p_period[i] = 0.0;
p_flag[i] = 0;
}
@ -498,9 +501,7 @@ FixRigid::FixRigid(LAMMPS *lmp, int narg, char **arg) :
else {
allremap = 0;
delete [] id_dilate;
int n = strlen(arg[iarg+1]) + 1;
id_dilate = new char[n];
strcpy(id_dilate,arg[iarg+1]);
id_dilate = utils::strdup(arg[iarg+1]);
int idilate = group->find(id_dilate);
if (idilate == -1)
error->all(FLERR,
@ -527,9 +528,7 @@ FixRigid::FixRigid(LAMMPS *lmp, int narg, char **arg) :
} else if (strcmp(arg[iarg],"infile") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal fix rigid command");
delete [] inpfile;
int n = strlen(arg[iarg+1]) + 1;
inpfile = new char[n];
strcpy(inpfile,arg[iarg+1]);
inpfile = utils::strdup(arg[iarg+1]);
restart_file = 1;
reinitflag = 0;
iarg += 2;
@ -544,9 +543,7 @@ FixRigid::FixRigid(LAMMPS *lmp, int narg, char **arg) :
} else if (strcmp(arg[iarg],"gravity") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal fix rigid command");
delete [] id_gravity;
int n = strlen(arg[iarg+1]) + 1;
id_gravity = new char[n];
strcpy(id_gravity,arg[iarg+1]);
id_gravity = utils::strdup(arg[iarg+1]);
iarg += 2;
} else error->all(FLERR,"Illegal fix rigid command");
@ -555,7 +552,7 @@ FixRigid::FixRigid(LAMMPS *lmp, int narg, char **arg) :
// set pstat_flag
pstat_flag = 0;
for (int i = 0; i < 3; i++)
for (i = 0; i < 3; i++)
if (p_flag[i]) pstat_flag = 1;
if (pcouple == XYZ || (dimension == 2 && pcouple == XY)) pstyle = ISO;
@ -606,8 +603,7 @@ FixRigid::FixRigid(LAMMPS *lmp, int narg, char **arg) :
for (ibody = 0; ibody < nbody; ibody++) nsum += nrigid[ibody];
if (me == 0)
utils::logmesg(lmp,fmt::format(" {} rigid bodies with {} atoms\n",
nbody,nsum));
utils::logmesg(lmp," {} rigid bodies with {} atoms\n",nbody,nsum);
}
/* ---------------------------------------------------------------------- */
@ -702,8 +698,8 @@ void FixRigid::init()
if (modify->fix[i]->rigid_flag) rflag = 1;
if (rflag && (modify->fmask[i] & POST_FORCE) &&
!modify->fix[i]->rigid_flag)
error->warning(FLERR,fmt::format("Fix {} alters forces after fix "
"rigid",modify->fix[i]->id));
error->warning(FLERR,"Fix {} alters forces after fix rigid",
modify->fix[i]->id);
}
}
@ -724,17 +720,14 @@ void FixRigid::init()
}
}
// error if npt,nph fix comes before rigid fix
for (i = 0; i < modify->nfix; i++) {
if (strcmp(modify->fix[i]->style,"npt") == 0) break;
if (strcmp(modify->fix[i]->style,"nph") == 0) break;
}
// error if a fix changing the box comes before rigid fix
for (i = 0; i < modify->nfix; i++)
if (modify->fix[i]->box_change) break;
if (i < modify->nfix) {
for (int j = i; j < modify->nfix; j++)
if (strcmp(modify->fix[j]->style,"rigid") == 0)
error->all(FLERR,"Rigid fix must come before NPT/NPH fix");
for (int j = i+1; j < modify->nfix; j++)
if (utils::strmatch(modify->fix[j]->style,"^rigid"))
error->all(FLERR,"Rigid fixes must come before any box changing fix");
}
// add gravity forces based on gravity vector from fix
@ -742,8 +735,8 @@ void FixRigid::init()
if (id_gravity) {
int ifix = modify->find_fix(id_gravity);
if (ifix < 0) error->all(FLERR,"Fix rigid cannot find fix gravity ID");
if (strcmp(modify->fix[ifix]->style,"gravity") != 0)
error->all(FLERR,"Fix rigid gravity fix is invalid");
if (!utils::strmatch(modify->fix[ifix]->style,"^gravity"))
error->all(FLERR,"Fix rigid gravity fix ID is not a gravity fix style");
int tmp;
gvec = (double *) modify->fix[ifix]->extract("gvec",tmp);
}
@ -754,7 +747,7 @@ void FixRigid::init()
dtf = 0.5 * update->dt * force->ftm2v;
dtq = 0.5 * update->dt;
if (strstr(update->integrate_style,"respa"))
if (utils::strmatch(update->integrate_style,"^respa"))
step_respa = ((Respa *) update->integrate)->step;
// setup rigid bodies, using current atom info. if reinitflag is not set,
@ -883,8 +876,7 @@ void FixRigid::setup(int vflag)
// virial setup before call to set_v
if (vflag) v_setup(vflag);
else evflag = 0;
v_init(vflag);
// set velocities from angmom & omega
@ -947,8 +939,7 @@ void FixRigid::initial_integrate(int vflag)
// virial setup before call to set_xv
if (vflag) v_setup(vflag);
else evflag = 0;
v_init(vflag);
// set coords/orient and velocity/rotation of atoms in rigid bodies
// from quarternion and omega
@ -2290,8 +2281,8 @@ void FixRigid::readfile(int which, double *vec,
if (me == 0) {
fp = fopen(inpfile,"r");
if (fp == nullptr)
error->one(FLERR,fmt::format("Cannot open fix rigid file {}: {}",
inpfile,utils::getsyserror()));
error->one(FLERR,"Cannot open fix rigid file {}: {}",
inpfile,utils::getsyserror());
while (1) {
eof = fgets(line,MAXLINE,fp);
if (eof == nullptr) error->one(FLERR,"Unexpected end of fix rigid file");
@ -2311,7 +2302,7 @@ void FixRigid::readfile(int which, double *vec,
int nread = 0;
while (nread < nlines) {
nchunk = MIN(nlines-nread,CHUNK);
eofflag = comm->read_lines_from_file(fp,nchunk,MAXLINE,buffer);
eofflag = utils::read_lines_from_file(fp,nchunk,MAXLINE,buffer,me,world);
if (eofflag) error->all(FLERR,"Unexpected end of fix rigid file");
buf = buffer;
@ -2399,8 +2390,8 @@ void FixRigid::write_restart_file(const char *file)
auto outfile = std::string(file) + ".rigid";
FILE *fp = fopen(outfile.c_str(),"w");
if (fp == nullptr)
error->one(FLERR,fmt::format("Cannot open fix rigid restart file {}: {}",
outfile,utils::getsyserror()));
error->one(FLERR,"Cannot open fix rigid restart file {}: {}",
outfile,utils::getsyserror());
fmt::print(fp,"# fix rigid mass, COM, inertia tensor info for "
"{} bodies on timestep {}\n\n",nbody,update->ntimestep);
@ -2448,14 +2439,14 @@ void FixRigid::write_restart_file(const char *file)
double FixRigid::memory_usage()
{
int nmax = atom->nmax;
double bytes = nmax * sizeof(int);
bytes += nmax * sizeof(imageint);
bytes += nmax*3 * sizeof(double);
bytes += maxvatom*6 * sizeof(double); // vatom
double bytes = (double)nmax * sizeof(int);
bytes += (double)nmax * sizeof(imageint);
bytes += (double)nmax*3 * sizeof(double);
bytes += (double)maxvatom*6 * sizeof(double); // vatom
if (extended) {
bytes += nmax * sizeof(int);
if (orientflag) bytes = nmax*orientflag * sizeof(double);
if (dorientflag) bytes = nmax*3 * sizeof(double);
bytes += (double)nmax * sizeof(int);
if (orientflag) bytes = (double)nmax*orientflag * sizeof(double);
if (dorientflag) bytes = (double)nmax*3 * sizeof(double);
}
return bytes;
}