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

This commit is contained in:
sjplimp
2013-09-04 17:03:05 +00:00
parent 06de308627
commit 2a6b344897
4 changed files with 52 additions and 20 deletions

View File

@ -195,7 +195,8 @@ ComputeReduce::ComputeReduce(LAMMPS *lmp, int narg, char **arg) :
"calculate a per-atom array"); "calculate a per-atom array");
if (argindex[i] && if (argindex[i] &&
argindex[i] > modify->compute[icompute]->size_peratom_cols) argindex[i] > modify->compute[icompute]->size_peratom_cols)
error->all(FLERR,"Compute reduce compute array is accessed out-of-range"); error->all(FLERR,
"Compute reduce compute array is accessed out-of-range");
} else if (modify->compute[icompute]->local_flag) { } else if (modify->compute[icompute]->local_flag) {
flavor[i] = LOCAL; flavor[i] = LOCAL;
if (argindex[i] == 0 && if (argindex[i] == 0 &&
@ -207,8 +208,10 @@ ComputeReduce::ComputeReduce(LAMMPS *lmp, int narg, char **arg) :
"calculate a local array"); "calculate a local array");
if (argindex[i] && if (argindex[i] &&
argindex[i] > modify->compute[icompute]->size_local_cols) argindex[i] > modify->compute[icompute]->size_local_cols)
error->all(FLERR,"Compute reduce compute array is accessed out-of-range"); error->all(FLERR,
} else error->all(FLERR,"Compute reduce compute calculates global values"); "Compute reduce compute array is accessed out-of-range");
} else error->all(FLERR,
"Compute reduce compute calculates global values");
} else if (which[i] == FIX) { } else if (which[i] == FIX) {
int ifix = modify->find_fix(ids[i]); int ifix = modify->find_fix(ids[i]);

View File

@ -250,14 +250,6 @@ void Respa::init()
if (modify->nfix == 0 && comm->me == 0) if (modify->nfix == 0 && comm->me == 0)
error->warning(FLERR,"No fixes defined, atoms won't move"); error->warning(FLERR,"No fixes defined, atoms won't move");
// warn about incorrect pressures when using rRESPA with fix SHAKE
int shakeflag = 0;
for (int i = 0; i < modify->nfix; i++)
if (strcmp(modify->fix[i]->style,"shake") == 0) shakeflag = 1;
if (shakeflag && comm->me == 0)
error->warning(FLERR,"Fix shake with rRESPA computes invalid pressures");
// create fix needed for storing atom-based respa level forces // create fix needed for storing atom-based respa level forces
// will delete it at end of run // will delete it at end of run

View File

@ -26,6 +26,9 @@
#include "variable.h" #include "variable.h"
#include "force.h" #include "force.h"
#include "modify.h" #include "modify.h"
#include "fix.h"
#include "fix_rigid.h"
#include "fix_rigid_small.h"
#include "compute.h" #include "compute.h"
#include "compute_temp.h" #include "compute_temp.h"
#include "random_park.h" #include "random_park.h"
@ -86,6 +89,7 @@ void Velocity::command(int narg, char **arg)
rotation_flag = 0; rotation_flag = 0;
loop_flag = ALL; loop_flag = ALL;
scale_flag = 1; scale_flag = 1;
rfix = -1;
// read options from end of input line // read options from end of input line
// change defaults as options specify // change defaults as options specify
@ -617,13 +621,41 @@ void Velocity::ramp(int narg, char **arg)
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
zero linear or angular momentum of a group zero linear or angular momentum of a group
if using rigid/small requires init of entire system since
its methods perform forward/reverse comm,
comm::init needs neighbor::init needs pair::init needs kspace::init, etc
also requires setup_pre_neighbor call to setup bodies
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Velocity::zero(int narg, char **arg) void Velocity::zero(int narg, char **arg)
{ {
if (strcmp(arg[0],"linear") == 0) zero_momentum(); if (strcmp(arg[0],"linear") == 0) {
else if (strcmp(arg[0],"angular") == 0) zero_rotation(); if (rfix < 0) zero_momentum();
else error->all(FLERR,"Illegal velocity command"); else {
if (strcmp(modify->fix[rfix]->style,"rigid") == 0)
((FixRigid *) modify->fix[rfix])->zero_momentum(igroup);
else if (strcmp(modify->fix[rfix]->style,"rigid/small") == 0) {
lmp->init();
((FixRigidSmall *) modify->fix[rfix])->setup_pre_neighbor();
((FixRigidSmall *) modify->fix[rfix])->zero_momentum(igroup);
}
else error->all(FLERR,"Velocity rigid used with non-rigid fix-ID");
}
} else if (strcmp(arg[0],"angular") == 0) {
if (rfix < 0) zero_rotation();
else {
if (strcmp(modify->fix[rfix]->style,"rigid") == 0)
((FixRigid *) modify->fix[rfix])->zero_rotation(igroup);
else if (strcmp(modify->fix[rfix]->style,"rigid/small") == 0) {
lmp->init();
((FixRigidSmall *) modify->fix[rfix])->setup_pre_neighbor();
((FixRigidSmall *) modify->fix[rfix])->zero_rotation(igroup);
}
else error->all(FLERR,"Velocity rigid used with non-rigid fix-ID");
}
} else error->all(FLERR,"Illegal velocity command");
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
@ -654,10 +686,10 @@ void Velocity::rescale(double t_old, double t_new)
void Velocity::zero_momentum() void Velocity::zero_momentum()
{ {
// cannot have 0 atoms in group // cannot have no atoms in group
if (group->count(igroup) == 0) if (group->count(igroup) == 0)
error->all(FLERR,"Cannot zero momentum of 0 atoms"); error->all(FLERR,"Cannot zero momentum of no atoms");
// compute velocity of center-of-mass of group // compute velocity of center-of-mass of group
@ -687,10 +719,10 @@ void Velocity::zero_rotation()
{ {
int i; int i;
// cannot have 0 atoms in group // cannot have no atoms in group
if (group->count(igroup) == 0) if (group->count(igroup) == 0)
error->all(FLERR,"Cannot zero momentum of 0 atoms"); error->all(FLERR,"Cannot zero momentum of no atoms");
// compute omega (angular velocity) of group around center-of-mass // compute omega (angular velocity) of group around center-of-mass
@ -717,7 +749,7 @@ void Velocity::zero_rotation()
for (i = 0; i < nlocal; i++) for (i = 0; i < nlocal; i++)
if (mask[i] & groupbit) { if (mask[i] & groupbit) {
domain->unmap(x[i],image[i],unwrap); domain->unmap(x[i],image[i],unwrap);
dx = unwrap[0]- xcm[0]; dx = unwrap[0] - xcm[0];
dy = unwrap[1] - xcm[1]; dy = unwrap[1] - xcm[1];
dz = unwrap[2] - xcm[2]; dz = unwrap[2] - xcm[2];
v[i][0] -= omega[1]*dz - omega[2]*dy; v[i][0] -= omega[1]*dz - omega[2]*dy;
@ -779,6 +811,11 @@ void Velocity::options(int narg, char **arg)
else if (strcmp(arg[iarg+1],"geom") == 0) loop_flag = GEOM; else if (strcmp(arg[iarg+1],"geom") == 0) loop_flag = GEOM;
else error->all(FLERR,"Illegal velocity command"); else error->all(FLERR,"Illegal velocity command");
iarg += 2; iarg += 2;
} else if (strcmp(arg[iarg],"rigid") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal velocity command");
rfix = modify->find_fix(arg[iarg+1]);
if (rfix < 0) error->all(FLERR,"Fix ID for velocity does not exist");
iarg += 2;
} else if (strcmp(arg[iarg],"units") == 0) { } else if (strcmp(arg[iarg],"units") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal velocity command"); if (iarg+2 > narg) error->all(FLERR,"Illegal velocity command");
if (strcmp(arg[iarg+1],"box") == 0) scale_flag = 0; if (strcmp(arg[iarg+1],"box") == 0) scale_flag = 0;

View File

@ -35,7 +35,7 @@ class Velocity : protected Pointers {
private: private:
int igroup,groupbit; int igroup,groupbit;
int style; int style;
int dist_flag,sum_flag,momentum_flag,rotation_flag,loop_flag,scale_flag; int dist_flag,sum_flag,momentum_flag,rotation_flag,loop_flag,scale_flag,rfix;
double xscale,yscale,zscale; double xscale,yscale,zscale;
class Compute *temperature; class Compute *temperature;