git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@13332 f3b2605a-c512-4ea7-a41b-209d697bcdaa

This commit is contained in:
sjplimp 2015-03-30 23:17:06 +00:00
parent 8d078d0602
commit bd37304ade
27 changed files with 135 additions and 39 deletions

View File

@ -133,7 +133,6 @@ void ComputeTempAsphere::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -141,7 +140,7 @@ void ComputeTempAsphere::setup()
void ComputeTempAsphere::dof_compute()
{
if (fix_dof) adjust_dof_fix();
fix_dof = modify->adjust_dof_fix(igroup);
// 6 dof for 3d, 3 dof for 2d
// which dof are included also depends on mode
@ -268,6 +267,8 @@ double ComputeTempAsphere::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic || tempbias == 2) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -201,7 +201,6 @@ void ComputeTempCS::setup()
// calculate DOF for temperature
fix_dof = -1;
dof_compute();
}
@ -209,7 +208,7 @@ void ComputeTempCS::setup()
void ComputeTempCS::dof_compute()
{
if (fix_dof) adjust_dof_fix();
fix_dof = modify->adjust_dof_fix(igroup);
int nper = domain->dimension;
double natoms = group->count(igroup);
dof = nper * natoms;
@ -257,6 +256,8 @@ double ComputeTempCS::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -249,6 +249,41 @@ void FixDeposit::init()
error->all(FLERR,"Fix deposit and fix shake not using "
"same molecule template ID");
}
// for finite size spherical particles:
// warn if near < 2 * maxrad of existing and inserted particles
// since may lead to overlaps
// if inserted molecule does not define diameters,
// use AtomVecSphere::create_atom() default radius = 0.5
if (atom->radius_flag) {
double *radius = atom->radius;
int nlocal = atom->nlocal;
double maxrad = 0.0;
for (int i = 0; i < nlocal; i++)
maxrad = MAX(maxrad,radius[i]);
double maxradall;
MPI_Allreduce(&maxrad,&maxradall,1,MPI_DOUBLE,MPI_MAX,world);
double maxradinsert = 0.0;
if (mode == MOLECULE) {
for (int i = 0; i < nmol; i++) {
if (onemols[i]->radiusflag)
maxradinsert = MAX(maxradinsert,onemols[i]->maxradius);
else maxradinsert = MAX(maxradinsert,0.5);
}
} else maxradinsert = 0.5;
double separation = MAX(2.0*maxradinsert,maxradall+maxradinsert);
if (sqrt(nearsq) < separation && comm->me == 0) {
char str[128];
sprintf(str,"Fix deposit near setting < possible overlap separation %g",
separation);
error->warning(FLERR,str);
}
}
}
/* ----------------------------------------------------------------------
@ -389,7 +424,8 @@ void FixDeposit::pre_exchange()
}
}
// if distance to any inserted atom is less than near, try again
// check distance between any existing atom and any inserted atom
// if less than near, try again
// use minimum_image() to account for PBC
double **x = atom->x;

View File

@ -71,6 +71,7 @@ FixPOEMS::FixPOEMS(LAMMPS *lmp, int narg, char **arg) :
time_integrate = 1;
rigid_flag = 1;
virial_flag = 1;
dof_flag = 1;
MPI_Comm_rank(world,&me);

View File

@ -66,6 +66,7 @@ FixRigid::FixRigid(LAMMPS *lmp, int narg, char **arg) :
rigid_flag = 1;
virial_flag = 1;
create_attribute = 1;
dof_flag = 1;
MPI_Comm_rank(world,&me);
MPI_Comm_size(world,&nprocs);

View File

@ -78,6 +78,7 @@ FixRigidSmall::FixRigidSmall(LAMMPS *lmp, int narg, char **arg) :
rigid_flag = 1;
virial_flag = 1;
create_attribute = 1;
dof_flag = 1;
MPI_Comm_rank(world,&me);
MPI_Comm_size(world,&nprocs);
@ -2191,7 +2192,7 @@ void FixRigidSmall::setup_bodies_static()
/* ----------------------------------------------------------------------
one-time initialization of dynamic rigid body attributes
Vcm and angmom, computed explicitly from constituent particles
vcm and angmom, computed explicitly from constituent particles
even if wrong for overlapping particles, is OK,
since is just setting initial time=0 Vcm and angmom of the body
which can be estimated value

View File

@ -55,6 +55,7 @@ FixShake::FixShake(LAMMPS *lmp, int narg, char **arg) :
virial_flag = 1;
create_attribute = 1;
dof_flag = 1;
// error check

View File

@ -63,6 +63,7 @@ FixShakeCuda::FixShakeCuda(LAMMPS* lmp, int narg, char** arg) :
virial_flag = 1;
create_attribute = 1;
dof_flag = 1;
// error check

View File

@ -93,7 +93,6 @@ void ComputeTempDeformEff::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -101,6 +100,7 @@ void ComputeTempDeformEff::setup()
void ComputeTempDeformEff::dof_compute()
{
fix_dof = modify->adjust_dof_fix(igroup);
double natoms = group->count(igroup);
dof = domain->dimension * natoms;
dof -= extra_dof + fix_dof;
@ -120,6 +120,7 @@ void ComputeTempDeformEff::dof_compute()
MPI_Allreduce(&one,&nelectrons,1,MPI_INT,MPI_SUM,world);
// Assume 3/2 k T per nucleus
dof -= domain->dimension * nelectrons;
if (dof > 0) tfactor = force->mvv2e / (dof * force->boltz);
@ -173,6 +174,8 @@ double ComputeTempDeformEff::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -59,7 +59,6 @@ void ComputeTempEff::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -67,7 +66,7 @@ void ComputeTempEff::setup()
void ComputeTempEff::dof_compute()
{
if (fix_dof) adjust_dof_fix();
fix_dof = adjust_dof_fix(igroup);
double natoms = group->count(igroup);
dof = domain->dimension * natoms;
dof -= extra_dof + fix_dof;
@ -121,6 +120,8 @@ double ComputeTempEff::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -51,6 +51,7 @@ FixLbRigidPCSphere::FixLbRigidPCSphere(LAMMPS *lmp, int narg, char **arg) :
rigid_flag = 1;
create_attribute = 1;
virial_flag = 1;
dof_flag = 1;
// perform initial allocation of atom-based arrays
// register with Atom class

View File

@ -70,7 +70,6 @@ void ComputeTempRotate::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -78,7 +77,7 @@ void ComputeTempRotate::setup()
void ComputeTempRotate::dof_compute()
{
if (fix_dof) adjust_dof_fix();
fix_dof = adjust_dof_fix(igroup);
double natoms = group->count(igroup);
int nper = domain->dimension;
dof = nper * natoms;
@ -144,6 +143,8 @@ double ComputeTempRotate::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -143,17 +143,6 @@ void Compute::modify_params(int narg, char **arg)
}
}
/* ----------------------------------------------------------------------
calculate adjustment in DOF due to fixes
------------------------------------------------------------------------- */
void Compute::adjust_dof_fix()
{
fix_dof = 0;
for (int i = 0; i < modify->nfix; i++)
fix_dof += modify->fix[i]->dof(igroup);
}
/* ----------------------------------------------------------------------
reset extra_dof to its default value
------------------------------------------------------------------------- */

View File

@ -92,7 +92,6 @@ class Compute : protected Pointers {
Compute(class LAMMPS *, int, char **);
virtual ~Compute();
void modify_params(int, char **);
void adjust_dof_fix();
void reset_extra_dof();
virtual void init() = 0;

View File

@ -18,6 +18,8 @@
#include "update.h"
#include "force.h"
#include "domain.h"
#include "modify.h"
#include "comm.h"
#include "group.h"
#include "error.h"
@ -52,7 +54,6 @@ void ComputeTemp::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -60,7 +61,7 @@ void ComputeTemp::setup()
void ComputeTemp::dof_compute()
{
if (fix_dof) adjust_dof_fix();
fix_dof = modify->adjust_dof_fix(igroup);
double natoms = group->count(igroup);
dof = domain->dimension * natoms;
dof -= extra_dof + fix_dof;
@ -96,6 +97,8 @@ double ComputeTemp::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -20,6 +20,7 @@
#include "force.h"
#include "group.h"
#include "domain.h"
#include "modify.h"
#include "lattice.h"
#include "error.h"
@ -62,7 +63,6 @@ void ComputeTempCOM::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -70,7 +70,7 @@ void ComputeTempCOM::setup()
void ComputeTempCOM::dof_compute()
{
if (fix_dof) adjust_dof_fix();
fix_dof = modify->adjust_dof_fix(igroup);
double natoms = group->count(igroup);
int nper = domain->dimension;
dof = nper * natoms;
@ -113,6 +113,8 @@ double ComputeTempCOM::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -88,7 +88,6 @@ void ComputeTempDeform::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -96,7 +95,7 @@ void ComputeTempDeform::setup()
void ComputeTempDeform::dof_compute()
{
if (fix_dof) adjust_dof_fix();
fix_dof = modify->adjust_dof_fix(igroup);
double natoms = group->count(igroup);
dof = domain->dimension * natoms;
dof -= extra_dof + fix_dof;
@ -149,6 +148,8 @@ double ComputeTempDeform::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -18,6 +18,7 @@
#include "update.h"
#include "force.h"
#include "domain.h"
#include "modify.h"
#include "group.h"
#include "memory.h"
#include "error.h"
@ -63,7 +64,6 @@ void ComputeTempPartial::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -74,7 +74,7 @@ void ComputeTempPartial::setup()
void ComputeTempPartial::dof_compute()
{
if (fix_dof) adjust_dof_fix();
fix_dof = modify->adjust_dof_fix(igroup);
double natoms = group->count(igroup);
int nper = xflag+yflag+zflag;
dof = nper * natoms;
@ -120,6 +120,8 @@ double ComputeTempPartial::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -161,7 +161,6 @@ ComputeTempProfile::~ComputeTempProfile()
void ComputeTempProfile::init()
{
fix_dof = -1;
dof_compute();
// ptrs to domain data
@ -189,7 +188,6 @@ void ComputeTempProfile::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -197,7 +195,7 @@ void ComputeTempProfile::setup()
void ComputeTempProfile::dof_compute()
{
if (fix_dof) adjust_dof_fix();
fix_dof = modify->adjust_dof_fix(igroup);
double natoms = group->count(igroup);
int nper = domain->dimension;
dof = nper * natoms;
@ -245,6 +243,8 @@ double ComputeTempProfile::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -119,7 +119,6 @@ void ComputeTempRamp::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -127,7 +126,7 @@ void ComputeTempRamp::setup()
void ComputeTempRamp::dof_compute()
{
if (fix_dof) adjust_dof_fix();
fix_dof = modify->adjust_dof_fix(igroup);
double natoms = group->count(igroup);
int nper = domain->dimension;
dof = nper * natoms;
@ -173,6 +172,8 @@ double ComputeTempRamp::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -20,6 +20,7 @@
#include "force.h"
#include "domain.h"
#include "modify.h"
#include "comm.h"
#include "group.h"
#include "error.h"
@ -114,7 +115,6 @@ void ComputeTempSphere::setup()
{
dynamic = 0;
if (dynamic_user || group->dynamic[igroup]) dynamic = 1;
fix_dof = -1;
dof_compute();
}
@ -124,7 +124,7 @@ void ComputeTempSphere::dof_compute()
{
int count,count_all;
if (fix_dof) adjust_dof_fix();
fix_dof = modify->adjust_dof_fix(igroup);
// 6 or 3 dof for extended/point particles for 3d
// 3 or 2 dof for extended/point particles for 2d
@ -252,6 +252,8 @@ double ComputeTempSphere::compute_scalar()
MPI_Allreduce(&t,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
if (dynamic || tempbias == 2) dof_compute();
if (tfactor == 0.0 && scalar != 0.0)
error->all(FLERR,"Temperature compute degrees of freedom < 0");
scalar *= tfactor;
return scalar;
}

View File

@ -65,6 +65,7 @@ Fix::Fix(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp)
restart_pbc = 0;
wd_header = wd_section = 0;
dynamic_group_allow = 0;
dof_flag = 0;
cudable_comm = 0;
scalar_flag = vector_flag = array_flag = 0;

View File

@ -49,6 +49,7 @@ class Fix : protected Pointers {
int wd_header; // # of header values fix writes to data file
int wd_section; // # of sections fix writes to data file
int dynamic_group_allow; // 1 if can be used with dynamic group, else 0
int dof_flag; // 1 if has dof() method (not min_dof())
int cudable_comm; // 1 if fix has CUDA-enabled communication
int scalar_flag; // 0/1 if compute_scalar() function exists

View File

@ -64,6 +64,7 @@ Modify::Modify(LAMMPS *lmp) : Pointers(lmp)
end_of_step_every = NULL;
list_dofflag = NULL;
list_timeflag = NULL;
nfix_restart_global = 0;
@ -135,6 +136,7 @@ Modify::~Modify()
delete [] list_min_energy;
delete [] end_of_step_every;
delete [] list_dofflag;
delete [] list_timeflag;
restart_deallocate();
@ -184,6 +186,8 @@ void Modify::init()
list_init(MIN_POST_FORCE,n_min_post_force,list_min_post_force);
list_init(MIN_ENERGY,n_min_energy,list_min_energy);
list_init_dofflag(n_dofflag,list_dofflag);
// init each fix
// not sure if now needs to come before compute init
// used to b/c temperature computes called fix->dof() in their init,
@ -898,6 +902,19 @@ int Modify::check_package(const char *package_fix_name)
return 1;
}
/* ----------------------------------------------------------------------
loop over fixes with dof() method
accumulate # of DOFs removed by fixes and return it
called by temperature computes
------------------------------------------------------------------------- */
int Modify::adjust_dof_fix(int igroup)
{
int n = 0;
for (int ifix = 0; ifix < n_dofflag; ifix++)
n += fix[list_dofflag[ifix]]->dof(igroup);
}
/* ----------------------------------------------------------------------
add a new compute
------------------------------------------------------------------------- */
@ -1282,6 +1299,27 @@ void Modify::list_init_thermo_energy(int mask, int &n, int *&list)
if (fmask[i] & mask && fix[i]->thermo_energy) list[n++] = i;
}
/* ----------------------------------------------------------------------
create list of fix indices for thermo energy fixes
only added to list if fix has THERMO_ENERGY mask
and its thermo_energy flag was set via fix_modify
------------------------------------------------------------------------- */
void Modify::list_init_dofflag(int &n, int *&list)
{
delete [] list;
n = 0;
for (int i = 0; i < nfix; i++)
if (fix[i]->dof_flag) n++;
list = new int[n];
n = 0;
for (int i = 0; i < nfix; i++)
if (fix[i]->dof_flag) list[n++] = i;
}
/* ----------------------------------------------------------------------
create list of compute indices for computes which store invocation times
------------------------------------------------------------------------- */

View File

@ -88,6 +88,7 @@ class Modify : protected Pointers {
void delete_fix(const char *);
int find_fix(const char *);
int check_package(const char *);
int adjust_dof_fix(int);
void add_compute(int, char **, int trysuffix=0);
void modify_compute(int, char **);
@ -121,6 +122,9 @@ class Modify : protected Pointers {
int *end_of_step_every;
int n_dofflag; // list of fixes with dof() method
int *list_dofflag;
int n_timeflag; // list of computes that store time invocation
int *list_timeflag;
@ -137,6 +141,7 @@ class Modify : protected Pointers {
void list_init(int, int &, int *&);
void list_init_end_of_step(int, int &, int *&);
void list_init_thermo_energy(int, int &, int *&);
void list_init_dofflag(int &, int *&);
void list_init_compute();
private:

View File

@ -536,10 +536,12 @@ void Molecule::charges(char *line)
void Molecule::diameters(char *line)
{
int tmp;
maxradius = 0.0;
for (int i = 0; i < natoms; i++) {
readline(line);
sscanf(line,"%d %lg",&tmp,&radius[i]);
radius[i] *= 0.5;
maxradius = MAX(maxradius,radius[i]);
}
for (int i = 0; i < natoms; i++)

View File

@ -90,7 +90,8 @@ class Molecule : protected Pointers {
double ex[3],ey[3],ez[3]; // principal axes of molecule in space coords
double quat[4]; // quaternion for orientation of molecule
double molradius; // radius of molecule from COM,
double maxradius; // max radius of any atom in molecule
double molradius; // radius of molecule from geometric center
// including finite-size particle radii
int comatom; // index (1-Natom) of atom closest to COM
double maxextent; // furthest any atom in molecule is from comatom