From 5f580b66d71f91d2b4a398bae75eb1d34e46c8ac Mon Sep 17 00:00:00 2001
From: sjplimp
Date: Mon, 6 Feb 2012 15:18:58 +0000
Subject: [PATCH 01/12] git-svn-id:
svn://svn.icms.temple.edu/lammps-ro/trunk@7701
f3b2605a-c512-4ea7-a41b-209d697bcdaa
---
src/finish.cpp | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/src/finish.cpp b/src/finish.cpp
index 482aa95abc..61a7e13096 100644
--- a/src/finish.cpp
+++ b/src/finish.cpp
@@ -53,7 +53,7 @@ void Finish::end(int flag)
// recompute natoms in case atoms have been lost
bigint nblocal = atom->nlocal;
- MPI_Allreduce(&nblocal,&atom->nlocal,1,MPI_LMP_BIGINT,MPI_SUM,world);
+ MPI_Allreduce(&nblocal,&atom->natoms,1,MPI_LMP_BIGINT,MPI_SUM,world);
// choose flavors of statistical output
// flag determines caller
From c81b8c5d042d4556096da5c76c6596c2698402eb Mon Sep 17 00:00:00 2001
From: sjplimp
Date: Mon, 6 Feb 2012 16:58:46 +0000
Subject: [PATCH 02/12] git-svn-id:
svn://svn.icms.temple.edu/lammps-ro/trunk@7702
f3b2605a-c512-4ea7-a41b-209d697bcdaa
---
src/change_box.cpp | 441 ++++++++++++++++++++++++++++++++++++++++---
src/change_box.h | 61 +++++-
src/displace_box.cpp | 418 ----------------------------------------
src/displace_box.h | 105 -----------
src/domain.cpp | 31 ++-
src/domain.h | 3 +-
src/input.cpp | 2 +-
src/pair_hybrid.cpp | 2 +-
src/read_restart.cpp | 6 +-
9 files changed, 503 insertions(+), 566 deletions(-)
delete mode 100644 src/displace_box.cpp
delete mode 100644 src/displace_box.h
diff --git a/src/change_box.cpp b/src/change_box.cpp
index 34f2f67b0b..be83182742 100644
--- a/src/change_box.cpp
+++ b/src/change_box.cpp
@@ -11,17 +11,28 @@
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
+#include "lmptype.h"
+#include "mpi.h"
+#include "math.h"
+#include "stdlib.h"
#include "string.h"
#include "change_box.h"
-#include "domain.h"
+#include "atom.h"
#include "modify.h"
#include "fix.h"
+#include "domain.h"
+#include "lattice.h"
+#include "comm.h"
+#include "irregular.h"
#include "output.h"
+#include "group.h"
#include "error.h"
using namespace LAMMPS_NS;
-enum{ORTHO,TRICLINIC};
+enum{XYZ,TILT,BOUNDARY,ORTHO,TRICLINIC,SET,REMAP};
+enum{FINAL,DELTA,SCALE};
+enum{X,Y,Z,YZ,XZ,XY};
/* ---------------------------------------------------------------------- */
@@ -31,34 +42,414 @@ ChangeBox::ChangeBox(LAMMPS *lmp) : Pointers(lmp) {}
void ChangeBox::command(int narg, char **arg)
{
+ int i;
+
if (domain->box_exist == 0)
error->all(FLERR,"Change_box command before simulation box is defined");
- if (narg != 1) error->all(FLERR,"Illegal change_box command");
+ if (narg < 2) error->all(FLERR,"Illegal change_box command");
+ if (modify->nfix_restart_peratom)
+ error->all(FLERR,"Cannot change_box after "
+ "reading restart file with per-atom info");
- int style;
- if (strcmp(arg[0],"ortho") == 0) style = ORTHO;
- else if (strcmp(arg[0],"triclinic") == 0) style = TRICLINIC;
- else error->all(FLERR,"Illegal change_box command");
+ if (comm->me == 0 && screen) fprintf(screen,"Changing box ...\n");
- if (style == ORTHO && domain->triclinic == 0)
- error->all(FLERR,"Change_box operation is invalid");
- if (style == TRICLINIC && domain->triclinic == 1)
- error->all(FLERR,"Change_box operation is invalid");
- if (style == ORTHO &&
- (domain->xy != 0.0 || domain->yz != 0.0 || domain->xz != 0.0))
- error->all(FLERR,"Cannot change box to orthogonal when tilt is non-zero");
+ // group
- if (output->ndump)
- error->all(FLERR,"Cannot change box with dumps defined");
- for (int i = 0; i < modify->nfix; i++)
- if (modify->fix[i]->no_change_box)
- error->all(FLERR,"Cannot change box with certain fixes defined");
+ int igroup = group->find(arg[0]);
+ if (igroup == -1) error->all(FLERR,"Could not find change_box group ID");
+ int groupbit = group->bitmask[igroup];
- if (style == ORTHO) domain->triclinic = 0;
- else domain->triclinic = 1;
- domain->xy = domain->yz = domain->xz = 0.0;
+ // parse operation arguments
+ // allocate ops to max possible length
+ // volume option does not increment nops
- domain->set_global_box();
- if (style == TRICLINIC) domain->set_lamda_box();
- domain->set_local_box();
+ int dimension = domain->dimension;
+
+ ops = new Operation[narg-1];
+ nops = 0;
+
+ int index;
+ int iarg = 1;
+ while (iarg < narg) {
+ if (strcmp(arg[iarg],"x") == 0 || strcmp(arg[iarg],"y") == 0 ||
+ strcmp(arg[iarg],"z") == 0) {
+ if (iarg+2 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].style = XYZ;
+ if (strcmp(arg[iarg],"x") == 0) ops[nops].dim = X;
+ else if (strcmp(arg[iarg],"y") == 0) ops[nops].dim = Y;
+ else if (strcmp(arg[iarg],"z") == 0) ops[nops].dim = Z;
+
+ if (dimension == 2 && ops[nops].dim == Z)
+ error->all(FLERR,"Cannot change_box in z dimension for 2d simulation");
+
+ if (strcmp(arg[iarg+1],"final") == 0) {
+ if (iarg+4 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].flavor = FINAL;
+ ops[nops].flo = atof(arg[iarg+2]);
+ ops[nops].fhi = atof(arg[iarg+3]);
+ ops[nops].vdim1 = ops[nops].vdim2 = -1;
+ nops++;
+ iarg += 4;
+ } else if (strcmp(arg[iarg+1],"delta") == 0) {
+ if (iarg+4 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].flavor = DELTA;
+ ops[nops].dlo = atof(arg[iarg+2]);
+ ops[nops].dhi = atof(arg[iarg+3]);
+ ops[nops].vdim1 = ops[nops].vdim2 = -1;
+ nops++;
+ iarg += 4;
+ } else if (strcmp(arg[iarg+1],"scale") == 0) {
+ if (iarg+3 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].flavor = SCALE;
+ ops[nops].scale = atof(arg[iarg+2]);
+ ops[nops].vdim1 = ops[nops].vdim2 = -1;
+ nops++;
+ iarg += 3;
+ } else if (strcmp(arg[iarg+1],"volume") == 0) {
+ if (nops == 0 || ops[nops-1].style != XYZ ||
+ ops[nops].dim == ops[nops-1].dim)
+ error->all(FLERR,"Change_box volume used incorrectly");
+ if (ops[nops-1].vdim2 >= 0)
+ error->all(FLERR,"Change_box volume used incorrectly");
+ else if (ops[nops-1].vdim1 >= 0) ops[nops-1].vdim2 = ops[nops].dim;
+ else ops[nops-1].vdim1 = ops[nops].dim;
+ iarg += 2;
+
+ } else error->all(FLERR,"Illegal change_box command");
+
+ } else if (strcmp(arg[iarg],"xy") == 0 || strcmp(arg[iarg],"xz") == 0 ||
+ strcmp(arg[iarg],"yz") == 0) {
+ if (iarg+2 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].style = TILT;
+ if (strcmp(arg[iarg],"xy") == 0) ops[nops].dim = XY;
+ else if (strcmp(arg[iarg],"xz") == 0) ops[nops].dim = XZ;
+ else if (strcmp(arg[iarg],"yz") == 0) ops[nops].dim = YZ;
+
+ if (dimension == 2 && (ops[nops].dim == XZ || ops[nops].dim == YZ))
+ error->all(FLERR,"Cannot change_box in xz or yz for 2d simulation");
+
+ if (strcmp(arg[iarg+1],"final") == 0) {
+ if (iarg+3 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].flavor = FINAL;
+ ops[nops].ftilt = atof(arg[iarg+2]);
+ nops++;
+ iarg += 3;
+ } else if (strcmp(arg[iarg+1],"delta") == 0) {
+ if (iarg+3 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].flavor = DELTA;
+ ops[nops].dtilt = atof(arg[iarg+2]);
+ nops++;
+ iarg += 3;
+ } else error->all(FLERR,"Illegal change_box command");
+
+ } else if (strcmp(arg[iarg],"boundary") == 0) {
+ if (iarg+4 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].style = BOUNDARY;
+ ops[nops].boundindex = iarg+1;
+ nops++;
+ iarg += 4;
+
+ } else if (strcmp(arg[iarg],"ortho") == 0) {
+ if (iarg+2 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].style = ORTHO;
+ nops++;
+ iarg += 1;
+
+ } else if (strcmp(arg[iarg],"triclinic") == 0) {
+ if (iarg+2 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].style = TRICLINIC;
+ nops++;
+ iarg += 1;
+
+ } else if (strcmp(arg[iarg],"set") == 0) {
+ if (iarg+1 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].style = SET;
+ nops++;
+ iarg += 1;
+
+ } else if (strcmp(arg[iarg],"remap") == 0) {
+ if (iarg+1 > narg) error->all(FLERR,"Illegal change_box command");
+ ops[nops].style = REMAP;
+ nops++;
+ iarg += 1;
+
+ } else break;
+ }
+
+ if (nops == 0) error->all(FLERR,"Illegal change_box command");
+
+ // read options from end of input line
+
+ options(narg-iarg,&arg[iarg]);
+
+ // compute scale factors if FINAL,DELTA used since they have distance units
+
+ int flag = 0;
+ for (int i = 0; i < nops; i++)
+ if (ops[i].style == FINAL || ops[i].style == DELTA) flag = 1;
+
+ if (flag && scaleflag && domain->lattice == NULL)
+ error->all(FLERR,"Use of change_box with undefined lattice");
+
+ if (flag && scaleflag) {
+ scale[0] = domain->lattice->xlattice;
+ scale[1] = domain->lattice->ylattice;
+ scale[2] = domain->lattice->zlattice;
+ }
+ else scale[0] = scale[1] = scale[2] = 1.0;
+
+ // perform sequence of operations
+ // first insure atoms are in current box & update box via shrink-wrap
+ // then save current box state so can remap atoms from it
+
+ domain->pbc();
+ domain->reset_box();
+ save_box_state();
+
+ for (int i = 0; i < nops; i++) {
+ if (ops[i].style == XYZ) {
+ double volume;
+ if (domain->dimension == 2) volume = domain->xprd * domain->yprd;
+ else volume = domain->xprd * domain->yprd * domain->zprd;
+
+ if (ops[i].flavor == FINAL) {
+ domain->boxlo[ops[i].dim] = scale[ops[i].dim]*ops[i].flo;
+ domain->boxhi[ops[i].dim] = scale[ops[i].dim]*ops[i].fhi;
+ if (ops[i].vdim1)
+ volume_preserve(ops[i].vdim1,ops[i].vdim2,volume);
+ domain->set_initial_box();
+ domain->set_global_box();
+ domain->set_local_box();
+ domain->print_box(" ");
+
+ } else if (ops[i].flavor == DELTA) {
+ domain->boxlo[ops[i].dim] += scale[ops[i].dim]*ops[i].dlo;
+ domain->boxhi[ops[i].dim] += scale[ops[i].dim]*ops[i].dhi;
+ if (ops[i].vdim1)
+ volume_preserve(ops[i].vdim1,ops[i].vdim2,volume);
+ domain->set_initial_box();
+ domain->set_global_box();
+ domain->set_local_box();
+ domain->print_box(" ");
+
+ } else if (ops[i].flavor == SCALE) {
+ double mid = 0.5 *
+ (domain->boxlo[ops[i].dim] + domain->boxhi[ops[i].dim]);
+ double delta = domain->boxlo[ops[i].dim] - mid;
+ domain->boxlo[ops[i].dim] = mid + ops[i].scale*delta;
+ delta = domain->boxhi[ops[i].dim] - mid;
+ domain->boxhi[ops[i].dim] = mid + ops[i].scale*delta;
+ if (ops[i].vdim1)
+ volume_preserve(ops[i].vdim1,ops[i].vdim2,volume);
+ domain->set_initial_box();
+ domain->set_global_box();
+ domain->set_local_box();
+ domain->print_box(" ");
+ }
+
+ } else if (ops[i].style == TILT) {
+ if (domain->triclinic == 0)
+ error->all(FLERR,"Cannot change box tilt factors for orthogonal box");
+
+ if (ops[i].flavor == FINAL) {
+ if (ops[i].dim == XY) domain->xy = scale[X]*ops[i].ftilt;
+ else if (ops[i].dim == XZ) domain->xz = scale[X]*ops[i].ftilt;
+ else if (ops[i].dim == YZ) domain->yz = scale[Y]*ops[i].ftilt;
+ domain->set_initial_box();
+ domain->set_global_box();
+ domain->set_local_box();
+ domain->print_box(" ");
+
+ } else if (ops[i].flavor == DELTA) {
+ if (ops[i].dim == XY) domain->xy += scale[X]*ops[i].dtilt;
+ else if (ops[i].dim == XZ) domain->xz += scale[X]*ops[i].dtilt;
+ else if (ops[i].dim == YZ) domain->yz += scale[Y]*ops[i].dtilt;
+ domain->set_initial_box();
+ domain->set_global_box();
+ domain->set_local_box();
+ domain->print_box(" ");
+ }
+
+ } else if (ops[i].style == BOUNDARY) {
+ domain->set_boundary(3,&arg[ops[i].boundindex],1);
+ if (domain->dimension == 2 && domain->zperiodic == 0)
+ error->all(FLERR,
+ "Cannot run 2d simulation with nonperiodic Z dimension");
+ domain->set_initial_box();
+ domain->set_global_box();
+ domain->set_local_box();
+
+ } else if (ops[i].style == ORTHO) {
+ if (domain->xy != 0.0 || domain->yz != 0.0 || domain->xz != 0.0)
+ error->all(FLERR,
+ "Cannot change box to orthogonal when tilt is non-zero");
+ if (output->ndump)
+ error->all(FLERR,
+ "Cannot change box ortho/triclinic with dumps defined");
+ for (int i = 0; i < modify->nfix; i++)
+ if (modify->fix[i]->no_change_box)
+ error->all(FLERR,
+ "Cannot change box ortho/triclinic with "
+ "certain fixes defined");
+ domain->triclinic = 0;
+ domain->set_initial_box();
+ domain->set_global_box();
+ domain->set_local_box();
+ domain->print_box(" ");
+
+ } else if (ops[i].style == TRICLINIC) {
+ if (output->ndump)
+ error->all(FLERR,
+ "Cannot change box ortho/triclinic with dumps defined");
+ for (int i = 0; i < modify->nfix; i++)
+ if (modify->fix[i]->no_change_box)
+ error->all(FLERR,
+ "Cannot change box ortho/triclinic with "
+ "certain fixes defined");
+ domain->triclinic = 1;
+ domain->set_lamda_box();
+ domain->set_initial_box();
+ domain->set_global_box();
+ domain->set_local_box();
+ domain->print_box(" ");
+
+ } else if (ops[i].style == SET) {
+ save_box_state();
+
+ } else if (ops[i].style == REMAP) {
+
+ // convert atoms to lamda coords, using last box state
+ // convert atoms back to box coords, using current box state
+ // save current box state
+
+ double **x = atom->x;
+ int *mask = atom->mask;
+ int nlocal = atom->nlocal;
+
+ for (i = 0; i < nlocal; i++)
+ if (mask[i] & groupbit)
+ domain->x2lamda(x[i],x[i],boxlo,h_inv);
+
+ for (i = 0; i < nlocal; i++)
+ if (mask[i] & groupbit)
+ domain->lamda2x(x[i],x[i]);
+
+ save_box_state();
+ }
+ }
+
+ // clean up
+
+ delete [] ops;
+
+ // apply shrink-wrap boundary conditions
+
+ if (domain->nonperiodic == 2) domain->reset_box();
+
+ // move atoms back inside simulation box and to new processors
+ // use remap() instead of pbc()
+ // in case box moved a long distance relative to atoms
+ // use irregular() in case box moved a long distance relative to atoms
+
+ double **x = atom->x;
+ int *image = atom->image;
+ int nlocal = atom->nlocal;
+ for (i = 0; i < nlocal; i++) domain->remap(x[i],image[i]);
+
+ if (domain->triclinic) domain->x2lamda(atom->nlocal);
+ domain->reset_box();
+ Irregular *irregular = new Irregular(lmp);
+ irregular->migrate_atoms();
+ delete irregular;
+ if (domain->triclinic) domain->lamda2x(atom->nlocal);
+
+ // check if any atoms were lost
+
+ domain->print_box(" ");
+
+ bigint natoms;
+ bigint nblocal = atom->nlocal;
+ MPI_Allreduce(&nblocal,&natoms,1,MPI_LMP_BIGINT,MPI_SUM,world);
+ if (natoms != atom->natoms) {
+ char str[128];
+ sprintf(str,"Lost atoms via change_box: original " BIGINT_FORMAT
+ " current " BIGINT_FORMAT,atom->natoms,natoms);
+ error->warning(FLERR,str);
+ }
+}
+
+/* ----------------------------------------------------------------------
+ parse optional parameters
+------------------------------------------------------------------------- */
+
+void ChangeBox::options(int narg, char **arg)
+{
+ if (narg < 0) error->all(FLERR,"Illegal change_box command");
+
+ scaleflag = 1;
+
+ int iarg = 0;
+ while (iarg < narg) {
+ if (strcmp(arg[iarg],"units") == 0) {
+ if (iarg+2 > narg) error->all(FLERR,"Illegal change_box command");
+ if (strcmp(arg[iarg+1],"box") == 0) scaleflag = 0;
+ else if (strcmp(arg[iarg+1],"lattice") == 0) scaleflag = 1;
+ else error->all(FLERR,"Illegal change_box command");
+ iarg += 2;
+ } else error->all(FLERR,"Illegal change_box command");
+ }
+}
+
+/* ----------------------------------------------------------------------
+ save current box state for converting atoms to lamda coords
+------------------------------------------------------------------------- */
+
+void ChangeBox::save_box_state()
+{
+ boxlo[0] = domain->boxlo[0];
+ boxlo[1] = domain->boxlo[1];
+ boxlo[2] = domain->boxlo[2];
+
+ for (int i = 0; i < 6; i++)
+ h_inv[i] = domain->h_inv[i];
+}
+
+/* ----------------------------------------------------------------------
+ reset box lengths of dim1/2 to preserve old volume
+ which changed due to change in 3rd dim
+------------------------------------------------------------------------- */
+
+void ChangeBox::volume_preserve(int dim1, int dim2, double oldvol)
+{
+ double newvol;
+ if (domain->dimension == 2) newvol = domain->xprd * domain->yprd;
+ else newvol = domain->xprd * domain->yprd * domain->zprd;
+
+ double scale = oldvol/newvol;
+ double mid,delta;
+
+ // change dim1 only
+
+ if (dim2 < 0) {
+ mid = 0.5 * (domain->boxlo[dim1] + domain->boxhi[dim1]);
+ delta = domain->boxlo[dim1] - mid;
+ domain->boxlo[dim1] = mid + scale*delta;
+ delta = domain->boxhi[dim1] - mid;
+ domain->boxhi[dim1] = mid + scale*delta;
+
+ // change dim1 and dim2, keeping their relative aspect constant
+ // means both are scaled by sqrt(scale)
+
+ } else {
+ mid = 0.5 * (domain->boxlo[dim1] + domain->boxhi[dim1]);
+ delta = domain->boxlo[dim1] - mid;
+ domain->boxlo[dim1] = mid + sqrt(scale)*delta;
+ delta = domain->boxhi[dim1] - mid;
+ domain->boxhi[dim1] = mid + sqrt(scale)*delta;
+ mid = 0.5 * (domain->boxlo[dim2] + domain->boxhi[dim2]);
+ delta = domain->boxlo[dim2] - mid;
+ domain->boxlo[dim2] = mid + sqrt(scale)*delta;
+ delta = domain->boxhi[dim2] - mid;
+ domain->boxhi[dim2] = mid + sqrt(scale)*delta;
+ }
}
diff --git a/src/change_box.h b/src/change_box.h
index 9399dcad8e..9b2eac97ae 100644
--- a/src/change_box.h
+++ b/src/change_box.h
@@ -28,6 +28,28 @@ class ChangeBox : protected Pointers {
public:
ChangeBox(class LAMMPS *);
void command(int, char **);
+
+ private:
+ int scaleflag;
+ double scale[3];
+
+ struct Operation {
+ int style,flavor;
+ int dim,boundindex;
+ int vdim1,vdim2;
+ double flo,fhi,ftilt;
+ double dlo,dhi,dtilt;
+ double scale;
+ };
+
+ Operation *ops;
+ int nops;
+
+ double boxlo[3],h_inv[6];
+
+ void options(int, char **);
+ void save_box_state();
+ void volume_preserve(int, int, double);
};
}
@@ -37,7 +59,7 @@ class ChangeBox : protected Pointers {
/* ERROR/WARNING messages:
-E: Change_box command before simulation box is defined
+E: Displace_box command before simulation box is defined
Self-explanatory.
@@ -47,22 +69,41 @@ Self-explanatory. Check the input script syntax and compare to the
documentation for the command. You can use -echo screen as a
command-line option when running LAMMPS to see the offending line.
-E: Change_box operation is invalid
+E: Cannot displace_box after reading restart file with per-atom info
-Cannot change orthogonal box to orthogonal or a triclinic box to
-triclinic.
+This is because the restart file info cannot be migrated with the
+atoms. You can get around this by performing a 0-timestep run which
+will assign the restart file info to actual atoms.
-E: Cannot change box to orthogonal when tilt is non-zero
+E: Could not find displace_box group ID
-Self-explanatory
+Group ID used in the displace_box command does not exist.
-E: Cannot change box with dumps defined
+E: Displace_box tilt factors require triclinic box
+
+Cannot use tilt factors unless the simulation box is
+non-orthogonal.
+
+E: Cannot displace_box on a non-periodic boundary
Self-explanatory.
-E: Cannot change box with certain fixes defined
+E: Use of displace_box with undefined lattice
-The change_box command cannot be used when fix ave/spatial or
-fix/deform are defined .
+Must use lattice command with displace_box command if units option is
+set to lattice.
+
+E: Fix deform volume setting is invalid
+
+Cannot use volume style unless other dimensions are being controlled.
+
+E: Induced tilt by displace_box is too large
+
+The final tilt value must be between -1/2 and 1/2 of the perpendicular
+box length.
+
+E: Lost atoms via displace_box: original %ld current %ld
+
+UNDOCUMENTED
*/
diff --git a/src/displace_box.cpp b/src/displace_box.cpp
deleted file mode 100644
index 79f8865788..0000000000
--- a/src/displace_box.cpp
+++ /dev/null
@@ -1,418 +0,0 @@
-/* ----------------------------------------------------------------------
- LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
- http://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 "lmptype.h"
-#include "mpi.h"
-#include "math.h"
-#include "stdlib.h"
-#include "string.h"
-#include "displace_box.h"
-#include "atom.h"
-#include "modify.h"
-#include "domain.h"
-#include "lattice.h"
-#include "comm.h"
-#include "irregular.h"
-#include "group.h"
-#include "error.h"
-
-using namespace LAMMPS_NS;
-
-enum{NONE,FINAL,DELTA,SCALE,VOLUME};
-enum{ONE_FROM_ONE,ONE_FROM_TWO,TWO_FROM_ONE};
-enum{NO_REMAP,X_REMAP};
-
-/* ---------------------------------------------------------------------- */
-
-DisplaceBox::DisplaceBox(LAMMPS *lmp) : Pointers(lmp) {}
-
-/* ---------------------------------------------------------------------- */
-
-void DisplaceBox::command(int narg, char **arg)
-{
- int i;
-
- if (domain->box_exist == 0)
- error->all(FLERR,"Displace_box command before simulation box is defined");
- if (narg < 2) error->all(FLERR,"Illegal displace_box command");
- if (modify->nfix_restart_peratom)
- error->all(FLERR,"Cannot displace_box after "
- "reading restart file with per-atom info");
-
- if (comm->me == 0 && screen) fprintf(screen,"Displacing box ...\n");
-
- // group
-
- int igroup = group->find(arg[0]);
- if (igroup == -1) error->all(FLERR,"Could not find displace_box group ID");
- int groupbit = group->bitmask[igroup];
-
- // set defaults
-
- set = new Set[6];
- set[0].style = set[1].style = set[2].style =
- set[3].style = set[4].style = set[5].style = NONE;
-
- // parse arguments
-
- int triclinic = domain->triclinic;
-
- int index;
- int iarg = 1;
- while (iarg < narg) {
- if (strcmp(arg[iarg],"x") == 0 || strcmp(arg[iarg],"y") == 0 ||
- strcmp(arg[iarg],"z") == 0) {
- if (strcmp(arg[iarg],"x") == 0) index = 0;
- else if (strcmp(arg[iarg],"y") == 0) index = 1;
- else if (strcmp(arg[iarg],"z") == 0) index = 2;
-
- if (iarg+2 > narg) error->all(FLERR,"Illegal displace_box command");
- if (strcmp(arg[iarg+1],"final") == 0) {
- if (iarg+4 > narg) error->all(FLERR,"Illegal displace_box command");
- set[index].style = FINAL;
- set[index].flo = atof(arg[iarg+2]);
- set[index].fhi = atof(arg[iarg+3]);
- iarg += 4;
- } else if (strcmp(arg[iarg+1],"delta") == 0) {
- if (iarg+4 > narg) error->all(FLERR,"Illegal displace_box command");
- set[index].style = DELTA;
- set[index].dlo = atof(arg[iarg+2]);
- set[index].dhi = atof(arg[iarg+3]);
- iarg += 4;
- } else if (strcmp(arg[iarg+1],"scale") == 0) {
- if (iarg+3 > narg) error->all(FLERR,"Illegal displace_box command");
- set[index].style = SCALE;
- set[index].scale = atof(arg[iarg+2]);
- iarg += 3;
- } else if (strcmp(arg[iarg+1],"volume") == 0) {
- set[index].style = VOLUME;
- iarg += 2;
- } else error->all(FLERR,"Illegal displace_box command");
-
- } else if (strcmp(arg[iarg],"xy") == 0 || strcmp(arg[iarg],"xz") == 0 ||
- strcmp(arg[iarg],"yz") == 0) {
- if (triclinic == 0)
- error->all(FLERR,"Displace_box tilt factors require triclinic box");
- if (strcmp(arg[iarg],"xy") == 0) index = 5;
- else if (strcmp(arg[iarg],"xz") == 0) index = 4;
- else if (strcmp(arg[iarg],"yz") == 0) index = 3;
- if (iarg+2 > narg) error->all(FLERR,"Illegal displace_box command");
- if (strcmp(arg[iarg+1],"final") == 0) {
- if (iarg+3 > narg) error->all(FLERR,"Illegal displace_box command");
- set[index].style = FINAL;
- set[index].ftilt = atof(arg[iarg+2]);
- iarg += 3;
- } else if (strcmp(arg[iarg+1],"delta") == 0) {
- if (iarg+3 > narg) error->all(FLERR,"Illegal displace_box command");
- set[index].style = DELTA;
- set[index].dtilt = atof(arg[iarg+2]);
- iarg += 3;
- } else error->all(FLERR,"Illegal displace_box command");
-
- } else break;
- }
-
- // read options from end of input line
-
- options(narg-iarg,&arg[iarg]);
-
- // check periodicity
-
- if ((set[0].style && domain->xperiodic == 0) ||
- (set[1].style && domain->yperiodic == 0) ||
- (set[2].style && domain->zperiodic == 0))
- error->all(FLERR,"Cannot displace_box on a non-periodic boundary");
-
- if (set[3].style && (domain->yperiodic == 0 || domain->zperiodic == 0))
- error->all(FLERR,"Cannot displace_box on a non-periodic boundary");
- if (set[4].style && (domain->xperiodic == 0 || domain->zperiodic == 0))
- error->all(FLERR,"Cannot displace_box on a non-periodic boundary");
- if (set[5].style && (domain->xperiodic == 0 || domain->yperiodic == 0))
- error->all(FLERR,"Cannot displace_box on a non-periodic boundary");
-
- // apply scaling to FINAL,DELTA since they have distance units
-
- int flag = 0;
- for (int i = 0; i < 6; i++)
- if (set[i].style == FINAL || set[i].style == DELTA) flag = 1;
-
- if (flag && scaleflag && domain->lattice == NULL)
- error->all(FLERR,"Use of displace_box with undefined lattice");
-
- double xscale,yscale,zscale;
- if (flag && scaleflag) {
- xscale = domain->lattice->xlattice;
- yscale = domain->lattice->ylattice;
- zscale = domain->lattice->zlattice;
- }
- else xscale = yscale = zscale = 1.0;
-
- // for 3,4,5 scaling is in 1st dimension, e.g. x for xz
-
- double map[6];
- map[0] = xscale; map[1] = yscale; map[2] = zscale;
- map[3] = yscale; map[4] = xscale; map[5] = xscale;
-
- for (int i = 0; i < 3; i++) {
- if (set[i].style == FINAL) {
- set[i].flo *= map[i];
- set[i].fhi *= map[i];
- } else if (set[i].style == DELTA) {
- set[i].dlo *= map[i];
- set[i].dhi *= map[i];
- }
- }
-
- for (int i = 3; i < 6; i++) {
- if (set[i].style == FINAL) set[i].ftilt *= map[i];
- else if (set[i].style == DELTA) set[i].dtilt *= map[i];
- }
-
- // set initial/final values for box size and shape
- // final = initial if no setting
-
- for (int i = 0; i < 3; i++) {
- set[i].lo_stop = set[i].lo_start = domain->boxlo[i];
- set[i].hi_stop = set[i].hi_start = domain->boxhi[i];
-
- if (set[i].style == FINAL) {
- set[i].lo_stop = set[i].flo;
- set[i].hi_stop = set[i].fhi;
- } else if (set[i].style == DELTA) {
- set[i].lo_stop = set[i].lo_start + set[i].dlo;
- set[i].hi_stop = set[i].hi_start + set[i].dhi;
- } else if (set[i].style == SCALE) {
- set[i].lo_stop = 0.5*(set[i].lo_start+set[i].hi_start) -
- 0.5*set[i].scale*(set[i].hi_start-set[i].lo_start);
- set[i].hi_stop = 0.5*(set[i].lo_start+set[i].hi_start) +
- 0.5*set[i].scale*(set[i].hi_start-set[i].lo_start);
- }
- }
-
- for (int i = 3; i < 6; i++) {
- if (i == 5) set[i].tilt_start = domain->xy;
- else if (i == 4) set[i].tilt_start = domain->xz;
- else if (i == 3) set[i].tilt_start = domain->yz;
- set[i].tilt_stop = set[i].tilt_start;
-
- if (set[i].style == FINAL) {
- set[i].tilt_stop = set[i].ftilt;
- } else if (set[i].style == DELTA) {
- set[i].tilt_stop = set[i].tilt_start + set[i].dtilt;
- }
- }
-
- // for VOLUME, setup links to other dims
- // fixed, dynamic1,2, vol_start
-
- for (int i = 0; i < 3; i++) {
- set[i].vol_start = domain->xprd * domain->yprd * domain->zprd;
-
- if (set[i].style != VOLUME) continue;
- int other1 = (i+1) % 3;
- int other2 = (i+2) % 3;
-
- if (set[other1].style == NONE) {
- if (set[other2].style == NONE || set[other2].style == VOLUME)
- error->all(FLERR,"Fix deform volume setting is invalid");
- set[i].substyle = ONE_FROM_ONE;
- set[i].fixed = other1;
- set[i].dynamic1 = other2;
- } else if (set[other2].style == NONE) {
- if (set[other1].style == NONE || set[other1].style == VOLUME)
- error->all(FLERR,"Fix deform volume setting is invalid");
- set[i].substyle = ONE_FROM_ONE;
- set[i].fixed = other2;
- set[i].dynamic1 = other1;
- } else if (set[other1].style == VOLUME) {
- if (set[other2].style == NONE || set[other2].style == VOLUME)
- error->all(FLERR,"Fix deform volume setting is invalid");
- set[i].substyle = TWO_FROM_ONE;
- set[i].fixed = other1;
- set[i].dynamic1 = other2;
- } else if (set[other2].style == VOLUME) {
- if (set[other1].style == NONE || set[other1].style == VOLUME)
- error->all(FLERR,"Fix deform volume setting is invalid");
- set[i].substyle = TWO_FROM_ONE;
- set[i].fixed = other2;
- set[i].dynamic1 = other1;
- } else {
- set[i].substyle = ONE_FROM_TWO;
- set[i].dynamic2 = other1;
- set[i].dynamic2 = other2;
- }
- }
-
- // set new box size for VOLUME dims that are linked to other dims
-
- for (int i = 0; i < 3; i++) {
- if (set[i].style != VOLUME) continue;
-
- if (set[i].substyle == ONE_FROM_ONE) {
- set[i].lo_stop = 0.5*(set[i].lo_start+set[i].hi_start) -
- 0.5*(set[i].vol_start /
- (set[set[i].dynamic1].hi_stop -
- set[set[i].dynamic1].lo_stop) /
- (set[set[i].fixed].hi_start-set[set[i].fixed].lo_start));
- set[i].hi_stop = 0.5*(set[i].lo_start+set[i].hi_start) +
- 0.5*(set[i].vol_start /
- (set[set[i].dynamic1].hi_stop -
- set[set[i].dynamic1].lo_stop) /
- (set[set[i].fixed].hi_start-set[set[i].fixed].lo_start));
-
- } else if (set[i].substyle == ONE_FROM_TWO) {
- set[i].lo_stop = 0.5*(set[i].lo_start+set[i].hi_start) -
- 0.5*(set[i].vol_start /
- (set[set[i].dynamic1].hi_stop -
- set[set[i].dynamic1].lo_stop) /
- (set[set[i].dynamic2].hi_stop -
- set[set[i].dynamic2].lo_stop));
- set[i].hi_stop = 0.5*(set[i].lo_start+set[i].hi_start) +
- 0.5*(set[i].vol_start /
- (set[set[i].dynamic1].hi_stop -
- set[set[i].dynamic1].lo_stop) /
- (set[set[i].dynamic2].hi_stop -
- set[set[i].dynamic2].lo_stop));
-
- } else if (set[i].substyle == TWO_FROM_ONE) {
- set[i].lo_stop = 0.5*(set[i].lo_start+set[i].hi_start) -
- 0.5*sqrt(set[i].vol_start /
- (set[set[i].dynamic1].hi_stop -
- set[set[i].dynamic1].lo_stop) /
- (set[set[i].fixed].hi_start -
- set[set[i].fixed].lo_start) *
- (set[i].hi_start - set[i].lo_start));
- set[i].hi_stop = 0.5*(set[i].lo_start+set[i].hi_start) +
- 0.5*sqrt(set[i].vol_start /
- (set[set[i].dynamic1].hi_stop -
- set[set[i].dynamic1].lo_stop) /
- (set[set[i].fixed].hi_start -
- set[set[i].fixed].lo_start) *
- (set[i].hi_start - set[i].lo_start));
- }
- }
-
- // check that final tilt is not illegal value
-
- double xprd_stop = set[0].hi_stop - set[0].lo_stop;
- double yprd_stop = set[0].hi_stop - set[0].lo_stop;
-
- if (set[3].tilt_stop < -0.5*yprd_stop || set[3].tilt_stop > 0.5*yprd_stop ||
- set[4].tilt_stop < -0.5*xprd_stop || set[4].tilt_stop > 0.5*xprd_stop ||
- set[5].tilt_stop < -0.5*xprd_stop || set[5].tilt_stop > 0.5*xprd_stop)
- error->all(FLERR,"Induced tilt by displace_box is too large");
-
- // convert atoms to lamda coords
-
- if (remapflag == X_REMAP) {
- double **x = atom->x;
- int *mask = atom->mask;
- int nlocal = atom->nlocal;
-
- for (i = 0; i < nlocal; i++)
- if (mask[i] & groupbit)
- domain->x2lamda(x[i],x[i]);
- }
-
- // reset global and local box to new size/shape
-
- domain->boxlo[0] = set[0].lo_stop;
- domain->boxlo[1] = set[1].lo_stop;
- domain->boxlo[2] = set[2].lo_stop;
- domain->boxhi[0] = set[0].hi_stop;
- domain->boxhi[1] = set[1].hi_stop;
- domain->boxhi[2] = set[2].hi_stop;
-
- if (triclinic) {
- domain->yz = set[3].tilt_stop;
- domain->xz = set[4].tilt_stop;
- domain->xy = set[5].tilt_stop;
- }
-
- domain->set_global_box();
- domain->set_local_box();
-
- // convert atoms back to box coords
-
- if (remapflag == X_REMAP) {
- double **x = atom->x;
- int *mask = atom->mask;
- int nlocal = atom->nlocal;
-
- for (i = 0; i < nlocal; i++)
- if (mask[i] & groupbit)
- domain->lamda2x(x[i],x[i]);
- }
-
- // move atoms back inside simulation box and to new processors
- // use remap() instead of pbc()
- // in case box moved a long distance relative to atoms
- // use irregular() in case box moved a long distance relative to atoms
-
- double **x = atom->x;
- int *image = atom->image;
- int nlocal = atom->nlocal;
- for (i = 0; i < nlocal; i++) domain->remap(x[i],image[i]);
-
- if (domain->triclinic) domain->x2lamda(atom->nlocal);
- domain->reset_box();
- Irregular *irregular = new Irregular(lmp);
- irregular->migrate_atoms();
- delete irregular;
- if (domain->triclinic) domain->lamda2x(atom->nlocal);
-
- // clean up
-
- delete [] set;
-
- // check if any atoms were lost
-
- bigint natoms;
- bigint nblocal = atom->nlocal;
- MPI_Allreduce(&nblocal,&natoms,1,MPI_LMP_BIGINT,MPI_SUM,world);
- if (natoms != atom->natoms) {
- char str[128];
- sprintf(str,"Lost atoms via displace_box: original " BIGINT_FORMAT
- " current " BIGINT_FORMAT,atom->natoms,natoms);
- error->all(FLERR,str);
- }
-}
-
-/* ----------------------------------------------------------------------
- parse optional parameters at end of displace_box input line
-------------------------------------------------------------------------- */
-
-void DisplaceBox::options(int narg, char **arg)
-{
- if (narg < 0) error->all(FLERR,"Illegal displace_box command");
-
- remapflag = X_REMAP;
- scaleflag = 1;
-
- int iarg = 0;
- while (iarg < narg) {
- if (strcmp(arg[iarg],"remap") == 0) {
- if (iarg+2 > narg) error->all(FLERR,"Illegal displace_box command");
- if (strcmp(arg[iarg+1],"x") == 0) remapflag = X_REMAP;
- else if (strcmp(arg[iarg+1],"none") == 0) remapflag = NO_REMAP;
- else error->all(FLERR,"Illegal displace_box command");
- iarg += 2;
- } else if (strcmp(arg[iarg],"units") == 0) {
- if (iarg+2 > narg) error->all(FLERR,"Illegal displace_box command");
- if (strcmp(arg[iarg+1],"box") == 0) scaleflag = 0;
- else if (strcmp(arg[iarg+1],"lattice") == 0) scaleflag = 1;
- else error->all(FLERR,"Illegal displace_box command");
- iarg += 2;
- } else error->all(FLERR,"Illegal displace_box command");
- }
-}
diff --git a/src/displace_box.h b/src/displace_box.h
deleted file mode 100644
index 44df201c36..0000000000
--- a/src/displace_box.h
+++ /dev/null
@@ -1,105 +0,0 @@
-/* ----------------------------------------------------------------------
- LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
- http://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.
-------------------------------------------------------------------------- */
-
-#ifdef COMMAND_CLASS
-
-CommandStyle(displace_box,DisplaceBox)
-
-#else
-
-#ifndef LMP_DISPLACE_BOX_H
-#define LMP_DISPLACE_BOX_H
-
-#include "pointers.h"
-
-namespace LAMMPS_NS {
-
-class DisplaceBox : protected Pointers {
- public:
- DisplaceBox(class LAMMPS *);
- void command(int, char **);
-
- private:
- int remapflag,scaleflag;
-
- struct Set {
- int style,substyle;
- double flo,fhi,ftilt;
- double dlo,dhi,dtilt;
- double scale;
- double lo_start,hi_start;
- double lo_stop,hi_stop;
- double tilt_start,tilt_stop;
- double vol_start;
- int fixed,dynamic1,dynamic2;
- };
- Set *set;
-
- void options(int, char **);
-};
-
-}
-
-#endif
-#endif
-
-/* ERROR/WARNING messages:
-
-E: Displace_box command before simulation box is defined
-
-Self-explanatory.
-
-E: Illegal ... command
-
-Self-explanatory. Check the input script syntax and compare to the
-documentation for the command. You can use -echo screen as a
-command-line option when running LAMMPS to see the offending line.
-
-E: Cannot displace_box after reading restart file with per-atom info
-
-This is because the restart file info cannot be migrated with the
-atoms. You can get around this by performing a 0-timestep run which
-will assign the restart file info to actual atoms.
-
-E: Could not find displace_box group ID
-
-Group ID used in the displace_box command does not exist.
-
-E: Displace_box tilt factors require triclinic box
-
-Cannot use tilt factors unless the simulation box is
-non-orthogonal.
-
-E: Cannot displace_box on a non-periodic boundary
-
-Self-explanatory.
-
-E: Use of displace_box with undefined lattice
-
-Must use lattice command with displace_box command if units option is
-set to lattice.
-
-E: Fix deform volume setting is invalid
-
-Cannot use volume style unless other dimensions are being controlled.
-
-E: Induced tilt by displace_box is too large
-
-The final tilt value must be between -1/2 and 1/2 of the perpendicular
-box length.
-
-E: Lost atoms via displace_box: original %ld current %ld
-
-UNDOCUMENTED
-
-*/
diff --git a/src/domain.cpp b/src/domain.cpp
index c729ac16bb..a83f958f6a 100644
--- a/src/domain.cpp
+++ b/src/domain.cpp
@@ -1109,10 +1109,12 @@ int Domain::find_region(char *name)
}
/* ----------------------------------------------------------------------
- boundary settings from the input script
+ (re)set boundary settings
+ flag = 0, called from the input script
+ flag = 1, called from change box command
------------------------------------------------------------------------- */
-void Domain::set_boundary(int narg, char **arg)
+void Domain::set_boundary(int narg, char **arg, int flag)
{
if (narg != 3) error->all(FLERR,"Illegal boundary command");
@@ -1127,7 +1129,10 @@ void Domain::set_boundary(int narg, char **arg)
else if (c == 'f') boundary[idim][iside] = 1;
else if (c == 's') boundary[idim][iside] = 2;
else if (c == 'm') boundary[idim][iside] = 3;
- else error->all(FLERR,"Illegal boundary command");
+ else {
+ if (flag == 0) error->all(FLERR,"Illegal boundary command");
+ if (flag == 1) error->all(FLERR,"Illegal change_box command");
+ }
}
for (int idim = 0; idim < 3; idim++)
@@ -1257,6 +1262,26 @@ void Domain::x2lamda(double *x, double *lamda)
lamda[2] = h_inv[2]*delta[2];
}
+/* ----------------------------------------------------------------------
+ convert box coords to triclinic 0-1 lamda coords for one atom
+ use my_boxlo & my_h_inv stored by caller for previous state of box
+ lamda = H^-1 (x - x0)
+ x and lamda can point to same 3-vector
+------------------------------------------------------------------------- */
+
+void Domain::x2lamda(double *x, double *lamda,
+ double *my_boxlo, double *my_h_inv)
+{
+ double delta[3];
+ delta[0] = x[0] - my_boxlo[0];
+ delta[1] = x[1] - my_boxlo[1];
+ delta[2] = x[2] - my_boxlo[2];
+
+ lamda[0] = my_h_inv[0]*delta[0] + my_h_inv[5]*delta[1] + my_h_inv[4]*delta[2];
+ lamda[1] = my_h_inv[1]*delta[1] + my_h_inv[3]*delta[2];
+ lamda[2] = my_h_inv[2]*delta[2];
+}
+
/* ----------------------------------------------------------------------
convert 8 lamda corner pts of lo/hi box to box coords
return bboxlo/hi = bounding box around 8 corner pts in box coords
diff --git a/src/domain.h b/src/domain.h
index 0496783c1f..d21119cffb 100644
--- a/src/domain.h
+++ b/src/domain.h
@@ -107,13 +107,14 @@ class Domain : protected Pointers {
void add_region(int, char **);
void delete_region(int, char **);
int find_region(char *);
- void set_boundary(int, char **);
+ void set_boundary(int, char **, int);
void print_box(const char *);
virtual void lamda2x(int);
virtual void x2lamda(int);
void lamda2x(double *, double *);
void x2lamda(double *, double *);
+ void x2lamda(double *, double *, double *, double *);
void bbox(double *, double *, double *, double *);
void box_corners();
diff --git a/src/input.cpp b/src/input.cpp
index 865975c68c..9c1905168e 100644
--- a/src/input.cpp
+++ b/src/input.cpp
@@ -930,7 +930,7 @@ void Input::boundary()
{
if (domain->box_exist)
error->all(FLERR,"Boundary command after simulation box is defined");
- domain->set_boundary(narg,arg);
+ domain->set_boundary(narg,arg,0);
}
/* ---------------------------------------------------------------------- */
diff --git a/src/pair_hybrid.cpp b/src/pair_hybrid.cpp
index 2ffd9c096e..73dd0709d5 100644
--- a/src/pair_hybrid.cpp
+++ b/src/pair_hybrid.cpp
@@ -16,6 +16,7 @@
#include "string.h"
#include "ctype.h"
#include "pair_hybrid.h"
+#include "style_pair.h"
#include "atom.h"
#include "force.h"
#include "pair.h"
@@ -763,4 +764,3 @@ double PairHybrid::memory_usage()
for (int m = 0; m < nstyles; m++) bytes += styles[m]->memory_usage();
return bytes;
}
-
diff --git a/src/read_restart.cpp b/src/read_restart.cpp
index 02e313a8d0..c015449ff4 100644
--- a/src/read_restart.cpp
+++ b/src/read_restart.cpp
@@ -464,7 +464,8 @@ void ReadRestart::header()
if (flag == VERSION) {
char *version = read_char();
if (strcmp(version,universe->version) != 0 && me == 0) {
- error->warning(FLERR,"Restart file version does not match LAMMPS version");
+ error->warning(FLERR,
+ "Restart file version does not match LAMMPS version");
if (screen) fprintf(screen," restart file = %s, LAMMPS = %s\n",
version,universe->version);
}
@@ -581,7 +582,8 @@ void ReadRestart::header()
boundary[2][0] != domain->boundary[2][0] ||
boundary[2][1] != domain->boundary[2][1]) {
if (me == 0)
- error->warning(FLERR,"Restart file used different boundary settings, "
+ error->warning(FLERR,
+ "Restart file used different boundary settings, "
"using restart file values");
}
}
From a3c28800f4657db0c78232cd1766e7d4187d8744 Mon Sep 17 00:00:00 2001
From: sjplimp
Date: Mon, 6 Feb 2012 16:58:56 +0000
Subject: [PATCH 03/12] git-svn-id:
svn://svn.icms.temple.edu/lammps-ro/trunk@7703
f3b2605a-c512-4ea7-a41b-209d697bcdaa
---
doc/Section_commands.html | 24 +--
doc/Section_commands.txt | 3 +-
doc/Section_howto.html | 4 +-
doc/Section_howto.txt | 4 +-
doc/boundary.html | 5 +-
doc/boundary.txt | 5 +-
doc/change_box.html | 297 ++++++++++++++++++++++++++++++++++----
doc/change_box.txt | 293 +++++++++++++++++++++++++++++++++----
doc/displace_atoms.html | 36 ++---
doc/displace_atoms.txt | 36 ++---
doc/displace_box.html | 191 ------------------------
doc/displace_box.txt | 181 -----------------------
doc/fix_deform.html | 2 +-
doc/fix_deform.txt | 2 +-
14 files changed, 597 insertions(+), 486 deletions(-)
delete mode 100644 doc/displace_box.html
delete mode 100644 doc/displace_box.txt
diff --git a/doc/Section_commands.html b/doc/Section_commands.html
index 2fa418e8ca..7c1a07a6e7 100644
--- a/doc/Section_commands.html
+++ b/doc/Section_commands.html
@@ -290,7 +290,7 @@ in the command's documentation.
delete_atoms, delete_bonds,
displace_atoms,
-displace_box, minimize,
+change_box, minimize,
neb prd, run, temper
Miscellaneous:
@@ -317,17 +317,17 @@ in the command's documentation.
| angle_coeff | angle_style | atom_modify | atom_style | bond_coeff | bond_style |
| boundary | change_box | clear | communicate | compute | compute_modify |
| create_atoms | create_box | delete_atoms | delete_bonds | dielectric | dihedral_coeff |
-| dihedral_style | dimension | displace_atoms | displace_box | dump | dump image |
-| dump_modify | echo | fix | fix_modify | group | if |
-| improper_coeff | improper_style | include | jump | kspace_modify | kspace_style |
-| label | lattice | log | mass | minimize | min_modify |
-| min_style | neb | neigh_modify | neighbor | newton | next |
-| package | pair_coeff | pair_modify | pair_style | pair_write | partition |
-| prd | print | processors | quit | read_data | read_restart |
-| region | replicate | reset_timestep | restart | run | run_style |
-| set | shell | special_bonds | suffix | tad | temper |
-| thermo | thermo_modify | thermo_style | timestep | uncompute | undump |
-| unfix | units | variable | velocity | write_restart
+ |
| dihedral_style | dimension | displace_atoms | dump | dump image | dump_modify |
+| echo | fix | fix_modify | group | if | improper_coeff |
+| improper_style | include | jump | kspace_modify | kspace_style | label |
+| lattice | log | mass | minimize | min_modify | min_style |
+| neb | neigh_modify | neighbor | newton | next | package |
+| pair_coeff | pair_modify | pair_style | pair_write | partition | prd |
+| print | processors | quit | read_data | read_restart | region |
+| replicate | reset_timestep | restart | run | run_style | set |
+| shell | special_bonds | suffix | tad | temper | thermo |
+| thermo_modify | thermo_style | timestep | uncompute | undump | unfix |
+| units | variable | velocity | write_restart
|
diff --git a/doc/Section_commands.txt b/doc/Section_commands.txt
index 773a5ed031..08d67735b0 100644
--- a/doc/Section_commands.txt
+++ b/doc/Section_commands.txt
@@ -286,7 +286,7 @@ Actions:
"delete_atoms"_delete_atoms.html, "delete_bonds"_delete_bonds.html,
"displace_atoms"_displace_atoms.html,
-"displace_box"_displace_box.html, "minimize"_minimize.html,
+"change_box"_change_box.html, "minimize"_minimize.html,
"neb"_neb.html "prd"_prd.html, "run"_run.html, "temper"_temper.html
Miscellaneous:
@@ -330,7 +330,6 @@ in the command's documentation.
"dihedral_style"_dihedral_style.html,
"dimension"_dimension.html,
"displace_atoms"_displace_atoms.html,
-"displace_box"_displace_box.html,
"dump"_dump.html,
"dump image"_dump_image.html,
"dump_modify"_dump_modify.html,
diff --git a/doc/Section_howto.html b/doc/Section_howto.html
index 221685fea3..041124e229 100644
--- a/doc/Section_howto.html
+++ b/doc/Section_howto.html
@@ -1660,8 +1660,8 @@ derivative w.r.t. e_i, it follows that C_ij is also symmetric, with at
most 7*6/2 = 21 distinct elements.
At zero temperature, it is easy to estimate these derivatives by
-deforming the cell in one of the six directions using the command
-displace_box and measuring the change in the
+deforming the simulation box in one of the six directions using the
+change_box command and measuring the change in the
stress tensor. A general-purpose script that does this is given in the
examples/elastic directory described in this
section.
diff --git a/doc/Section_howto.txt b/doc/Section_howto.txt
index b19befe7f2..1f14a0a51b 100644
--- a/doc/Section_howto.txt
+++ b/doc/Section_howto.txt
@@ -1647,8 +1647,8 @@ derivative w.r.t. e_i, it follows that C_ij is also symmetric, with at
most 7*6/2 = 21 distinct elements.
At zero temperature, it is easy to estimate these derivatives by
-deforming the cell in one of the six directions using the command
-"displace_box"_displace_box.html and measuring the change in the
+deforming the simulation box in one of the six directions using the
+"change_box"_change_box.html command and measuring the change in the
stress tensor. A general-purpose script that does this is given in the
examples/elastic directory described in "this
section"_Section_example.html.
diff --git a/doc/boundary.html b/doc/boundary.html
index 38544ed2d8..602a3d4f0d 100644
--- a/doc/boundary.html
+++ b/doc/boundary.html
@@ -77,7 +77,10 @@ triclinic representations.
Restrictions:
This command cannot be used after the simulation box is defined by a
-read_data or create_box command.
+read_data or create_box command or
+read_restart command. See the
+change_box command for how to change the simulation
+box boundaries after it has been defined.
For 2d simulations, the z dimension must be periodic.
diff --git a/doc/boundary.txt b/doc/boundary.txt
index cfdb43ce21..ec559d7f2f 100644
--- a/doc/boundary.txt
+++ b/doc/boundary.txt
@@ -72,7 +72,10 @@ triclinic representations.
[Restrictions:]
This command cannot be used after the simulation box is defined by a
-"read_data"_read_data.html or "create_box"_create_box.html command.
+"read_data"_read_data.html or "create_box"_create_box.html command or
+"read_restart"_read_restart.html command. See the
+"change_box"_change_box.html command for how to change the simulation
+box boundaries after it has been defined.
For 2d simulations, the z dimension must be periodic.
diff --git a/doc/change_box.html b/doc/change_box.html
index e258ec2337..b504671ee9 100644
--- a/doc/change_box.html
+++ b/doc/change_box.html
@@ -13,51 +13,292 @@
Syntax:
-change_box style
+change_box group-ID parameter args ... keyword args ...
-
style = ortho or triclinic
+- group-ID = ID of group of atoms to (optionally) displace
-
ortho = convert simulation box from non-orthogonal (triclinic) to orthogonal
- triclinic = convert simulation box from orthogonal to non-orthogonal (triclinic)
+
- one or more parameter/arg pairs may be appended
+
+
parameter = x or y or z or xy or xz or yz or boundary or ortho or triclinic or set or remap
+ x, y, z args = style value(s)
+ style = final or delta or scale or volume
+ final values = lo hi
+ lo hi = box boundaries after displacement (distance units)
+ delta values = dlo dhi
+ dlo dhi = change in box boundaries after displacement (distance units)
+ scale values = factor
+ factor = multiplicative factor for change in box length after displacement
+ volume value = none = adjust this dim to preserve volume of system
+ xy, xz, yz args = style value
+ style = final or delta
+ final value = tilt
+ tilt = tilt factor after displacement (distance units)
+ delta value = dtilt
+ dtilt = change in tilt factor after displacement (distance units)
+ boundary args = x y z
+ x,y,z = p or s or f or m, one or two letters
+ p is periodic
+ f is non-periodic and fixed
+ s is non-periodic and shrink-wrapped
+ m is non-periodic and shrink-wrapped with a minimum value
+ ortho args = none = change box to orthogonal
+ (triclinic args = none = change box to triclinic
+ set args = none = store state of current box
+ remap args = none = remap atom coords from last saved state to current box
+ - zero or more keyword/value pairs may be appended
+
+
- keyword = units
+
+
units value = lattice or box
+ lattice = distances are defined in lattice units
+ box = distances are defined in simulation box units
+
+
+
Examples:
-change_box ortho
-change_box triclinic
+change_box all xy final -2.0 z final 0.0 5.0 boundary p p f remap units box
+change_box all x scale 1.1 y volume z volume remap
Description:
-By default LAMMPS runs a simulation in an orthogonal, axis-aligned
-simulation box. LAMMPS can also run simulations in non-orthogonal
-(triclinic) simulation boxes. A box is
-defined as either orthogonal or non-orthogonal when it is created via
-the create_box, read_data, or
-read_restart commands.
+
Change the volume and/or shape and/or boundary conditions for the
+simulation box. Orthogonal simulation boxes have 3 adjustable size
+parameters (x,y,z). Triclinic (non-orthogonal) simulation boxes have
+6 adjustable size/shape parameters (x,y,z,xy,xz,yz). Any or all of
+them can be adjusted independently by this command. Thus it can be
+used to expand or contract a box, or to apply a shear strain to a
+non-orthogonal box. It can also be used to change the boundary
+conditions for the simulation box, similar to the
+boundary command.
-This command allows you to toggle the existing simulation box from
-orthogonal to non-orthogonal and vice versa. For example, an initial
+
The size and shape of the initial simulation box are specified by the
+create_box or read_data or
+read_restart command used to setup the simulation.
+The size and shape may be altered by subsequent runs, e.g. by use of
+the fix npt or fix deform commands.
+The create_box, read data, and
+read_restart commands also determine whether the
+simulation box is orthogonal or triclinic and their doc pages explain
+the meaning of the xy,xz,yz tilt factors.
+
+See Section_howto 12 of the doc pages
+for a geometric description of triclinic boxes, as defined by LAMMPS,
+and how to transform these parameters to and from other commonly used
+triclinic representations.
+
+The keywords used in this command are applied sequentially to the
+simulation box and the atoms in it, in the order specified.
+
+Before the sequence of keywords are invoked, the current box
+size/shape are stored, in case a remap keyword is used to map the
+atom coordinates from the previous box size/shape to the current one.
+
+After all the keywords have been processed, any shrink-wrap boundary
+conditions are invoked (see the boundary command)
+which may change simulation box boundaries, and atoms are migrated to
+new owning processors.
+
+IMPORTANT NOTE: It is possible to lose atoms with this command.
+E.g. by changing the box without remapping the atoms, and having atoms
+end up outside of non-periodic boundaries. It is also possible when
+remapping atoms to put them (nearly) on top of each other which will
+lead to bad dynamics. E.g. by converting a boundary from non-periodic
+to periodic.
+
+IMPORTANT NOTE: The simulation box size/shape can be changed by
+arbitrarily large amounts by this command. This is not a problem,
+except that the mapping of processors to the simulation box is not
+changed from its initial 3d configuration; see the
+processors command. Thus, if the box size/shape
+changes dramatically, the mapping of processors to the simulation box
+may not end up as optimal as the initial mapping attempted to be.
+
+
+
+For the x, y, and z parameters, this is the meaning of their
+styles and values.
+
+For style final, the final lo and hi box boundaries of a dimension
+are specified. The values can be in lattice or box distance units.
+See the discussion of the units keyword below.
+
+For style delta, plus or minus changes in the lo/hi box boundaries
+of a dimension are specified. The values can be in lattice or box
+distance units. See the discussion of the units keyword below.
+
+For style scale, a multiplicative factor to apply to the box length
+of a dimension is specified. For example, if the initial box length
+is 10, and the factor is 1.1, then the final box length will be 11. A
+factor less than 1.0 means compression.
+
+The volume style changes the specified dimension in such a way that
+the overall box volume remains constant with respect to the operation
+performed by the preceding keyword. The volume style can only be
+used following a keyword that changed the volume, which is any of the
+x, y, z keywords. If the preceding keyword "key" had a volume
+style, then both it and the current keyword apply to the keyword
+preceding "key". I.e. this sequence of keywords is allowed:
+
+change_box all x scale 1.1 y volume z volume
+
+The volume style changes the associated dimension so that the
+overall box volume is unchanged relative to its value before the
+preceding keyword was invoked.
+
+If "x scale 1.1 z volume" is used, then the z box length will
+shrink by the same 1.1 factor the x box length was increased by.
+
+If "x scale 1.1 y volume z volume" is used, then the y,z box lengths
+will each shrink by sqrt(1.1) to keep the volume constant. In this
+case, the y,z box lengths shrink so as to keep their relative aspect
+ratio constant.
+
+If "x scale 1.1 z volume y scale 1.1 z volume" is used, then
+the final box will be a factor of 10% larger in x and y, and a
+factor of 21% smaller in z, so as to keep the volume constant.
+
+IMPORTANT NOTE: For solids or liquids, when one dimension of the box
+is expanded, it may be physically undesirable to hold the other 2 box
+lengths constant since that implies a density change. For solids,
+adjusting the other dimensions via the volume style may make
+physical sense (just as for a liquid), but may not be correct for
+materials and potentials whose Poisson ratio is not 0.5.
+
+For the scale and volume styles, the box length is expanded or
+compressed around its mid point.
+
+
+
+For the xy, xz, and yz parameters, this is the meaning of their
+styles and values. Note that changing the tilt factors of a triclinic
+box does not change its volume.
+
+For style final, the final tilt factor is specified. The value
+can be in lattice or box distance units. See the discussion of the
+units keyword below.
+
+For style delta, a plus or minus change in the tilt factor is
+specified. The value can be in lattice or box distance units. See
+the discussion of the units keyword below.
+
+All of these styles change the xy, xz, yz tilt factors. In LAMMPS,
+tilt factors (xy,xz,yz) for triclinic boxes are required to be no more
+than half the distance of the parallel box length. For example, if
+xlo = 2 and xhi = 12, then the x box length is 10 and the xy tilt
+factor must be between -5 and 5. Similarly, both xz and yz must be
+between -(xhi-xlo)/2 and +(yhi-ylo)/2. Note that this is not a
+limitation, since if the maximum tilt factor is 5 (as in this
+example), then configurations with tilt = ..., -15, -5, 5, 15, 25,
+... are all equivalent. Any tilt factor specified by this command
+must be within these limits.
+
+
+
+The boundary keyword takes arguments that have exactly the same
+meaning as they do for the boundary command. In each
+dimension, a single letter assigns the same style to both the lower
+and upper face of the box. Two letters assigns the first style to the
+lower face and the second style to the upper face.
+
+The style p means the box is periodic; the other styles mean
+non-periodic. For style f, the position of the face is fixed. For
+style s, the position of the face is set so as to encompass the
+atoms in that dimension (shrink-wrapping), no matter how far they
+move. For style m, shrink-wrapping occurs, but is bounded by the
+current box edge in that dimension, so that the box will become no
+smaller. See the boundary command for more
+explanation of these style options.
+
+Note that the "boundary" command itself can only be used before the
+simulation box is defined via a read_data or
+create_box or read_restart
+command. This command allows the boundary conditions to be changed
+later in your input script. Also note that the
+read_restart will change boundary conditions to
+match what is stored in the restart file. So if you wish to change
+them, you should use the change_box command after the read_restart
+command.
+
+
+
+The ortho and triclinic keywords convert the simulation box to be
+orthogonal or triclinic (non-orthongonal). See this
+section for a discussion of how non-orthongal
+boxes are represented in LAMMPS.
+
+The simulation box is defined as either orthogonal or triclinic when
+it is created via the create_box,
+read_data, or read_restart
+commands.
+
+These keywords allow you to toggle the existing simulation box from
+orthogonal to triclinic and vice versa. For example, an initial
equilibration simulation can be run in an orthogonal box, the box can
-be toggled to non-orthogonal, and then a non-equilibrium MD (NEMD)
+be toggled to triclinic, and then a non-equilibrium MD (NEMD)
simulation can be run with deformation
via the fix deform command.
-Note that if the simulation box is currently non-orthogonal and has
-non-zero tilt in xy, yz, or xz, then it cannot be converted to an
-orthogonal box.
+
If the simulation box is currently triclinic and has non-zero tilt in
+xy, yz, or xz, then it cannot be converted to an orthogonal box.
+
+
+The set keyword saves the current box size/shape. This can be
+useful if you wish to use the remap keyword more than once or if you
+wish it to be applied to an intermediate box size/shape in a sequence
+of keyword operations. Note that the box size/shape is saved before
+any of the keywords are processed, i.e. the box size/shape at the time
+the create_box command is encountered in the input script.
+
+The remap keyword remaps atom coordinates from the last saved box
+size/shape to the current box state. For example, if you stretch the
+box in the x dimension or tilt it in the xy plane via the x and xy
+keywords, then the remap commmand will dilate or tilt the atoms to
+conform to the new box size/shape, as if the atoms moved with the box
+as it deformed.
+
+Note that this operation is performed without regard to periodic
+boundaries. Any shrink-wrapping of non-periodic boundaries (see the
+boundary command occurs after all keywords, including
+this one, have been processed.
+
+Only atoms in the specified group are remapped.
+
+IMPORTANT NOTE: If you do not explicitly specify the remap keyword,
+atom coordinates will not be changed even though the box size/shape
+changes. This may be the behavior you desire, but can also cause
+atoms to be lost.
+
+
+
+The units keyword determines the meaning of the distance units used
+to define various arguments. A box value selects standard distance
+units as defined by the units command, e.g. Angstroms for
+units = real or metal. A lattice value means the distance units are
+in lattice spacings. The lattice command must have
+been previously used to define the lattice spacing.
+
+
+
Restrictions:
-At the point in the input script when this command is issued, no
-dumps can be active, nor can a fix
-ave/spatial or fix deform be
-active. This is because these commands test whether the simulation
-box is orthogonal when they are first issued. Note that these
-commands can appear in your script before a change_box command is
-issued, so long as an undump or unfix
-command is also used to turn them off.
+
If you use the ortho or triclinic keywords, then at the point in
+the input script when this command is issued, no dumps can
+be active, nor can a fix ave/spatial or fix
+deform be active. This is because these commands
+test whether the simulation box is orthogonal when they are first
+issued. Note that these commands can be used in your script before a
+change_box command is issued, so long as an undump or
+unfix command is also used to turn them off.
-Related commands: none
+
Related commands:
-Default: none
+
fix deform, boundary
+
+Default:
+
+The option default is units = lattice.