diff -Naur lammps-20Jul07/doc/change_box.html lammps-3Aug07/doc/change_box.html --- lammps-20Jul07/doc/change_box.html 2007-06-27 16:39:31.000000000 -0600 +++ lammps-3Aug07/doc/change_box.html 2007-08-03 13:53:18.000000000 -0600 @@ -47,11 +47,14 @@

Restrictions:

-

When this command is used, no dumps can be defined, nor -can fix ave/spatial or fix -deform be defined. This is because these commands -test whether the simulation box is orthogonal or not when they are -first issued. +

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 +commmands 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.

Related commands: none

diff -Naur lammps-20Jul07/doc/change_box.txt lammps-3Aug07/doc/change_box.txt --- lammps-20Jul07/doc/change_box.txt 2007-06-27 16:39:31.000000000 -0600 +++ lammps-3Aug07/doc/change_box.txt 2007-08-03 13:53:18.000000000 -0600 @@ -43,11 +43,14 @@ [Restrictions:] -When this command is used, no "dumps"_dump.html can be defined, nor -can "fix ave/spatial"_fix_ave_spatial.html or "fix -deform"_fix_deform.html be defined. This is because these commands -test whether the simulation box is orthogonal or not when they are -first issued. +At the point in the input script when this command is issued, no +"dumps"_dump.html can be active, nor can a "fix +ave/spatial"_fix_ave_spatial.html or "fix deform"_fix_deform.html be +active. This is because these commands test whether the simulation +box is orthogonal when they are first issued. Note that these +commmands can appear in your script before a change_box command is +issued, so long as an "undump"_undump.html or "unfix"_unfix.html +command is also used to turn them off. [Related commands:] none diff -Naur lammps-20Jul07/doc/compute_temp_deform.html lammps-3Aug07/doc/compute_temp_deform.html --- lammps-20Jul07/doc/compute_temp_deform.html 2007-06-22 17:41:35.000000000 -0600 +++ lammps-3Aug07/doc/compute_temp_deform.html 2007-08-03 13:53:18.000000000 -0600 @@ -46,13 +46,13 @@

IMPORTANT NOTE: Fix deform has an option for remapping either atom coordinates or velocities to the changing -simulation box. To use this compute, the fix should NOT remap atom +simulation box. To use this compute, fix deform should NOT remap atom positions, but rather should let atoms respond to the changing box by adjusting their own velocities (or let fix deform remap the atom -velocities). If the fix does remap atom positions, their velocity is -not changed, and thus they do not have the streaming velocity assumed -by this compute. LAMMPS will warn you if this setting is not -consistent. +velocities). If fix deform does remap atom positions, their velocity +is not changed, and thus they do not have the streaming velocity +assumed by this compute. LAMMPS will warn you if the fix deform +setting is not consistent with this compute.

The temperature is calculated by the formula KE = dim/2 N k T, where KE = total kinetic energy of the group of atoms (sum of 1/2 m v^2), diff -Naur lammps-20Jul07/doc/compute_temp_deform.txt lammps-3Aug07/doc/compute_temp_deform.txt --- lammps-20Jul07/doc/compute_temp_deform.txt 2007-06-22 17:41:35.000000000 -0600 +++ lammps-3Aug07/doc/compute_temp_deform.txt 2007-08-03 13:53:18.000000000 -0600 @@ -43,13 +43,13 @@ IMPORTANT NOTE: "Fix deform"_fix_deform.html has an option for remapping either atom coordinates or velocities to the changing -simulation box. To use this compute, the fix should NOT remap atom +simulation box. To use this compute, fix deform should NOT remap atom positions, but rather should let atoms respond to the changing box by adjusting their own velocities (or let fix deform remap the atom -velocities). If the fix does remap atom positions, their velocity is -not changed, and thus they do not have the streaming velocity assumed -by this compute. LAMMPS will warn you if this setting is not -consistent. +velocities). If fix deform does remap atom positions, their velocity +is not changed, and thus they do not have the streaming velocity +assumed by this compute. LAMMPS will warn you if the fix deform +setting is not consistent with this compute. The temperature is calculated by the formula KE = dim/2 N k T, where KE = total kinetic energy of the group of atoms (sum of 1/2 m v^2), diff -Naur lammps-20Jul07/doc/fix_deform.html lammps-3Aug07/doc/fix_deform.html --- lammps-20Jul07/doc/fix_deform.html 2007-07-09 08:09:39.000000000 -0600 +++ lammps-3Aug07/doc/fix_deform.html 2007-08-03 13:53:18.000000000 -0600 @@ -50,9 +50,8 @@ V = change tilt factor at this velocity (distance/time units), effectively an engineering shear strain rate erate value = R - R = engineering shear strain rate (1/time units) - -

      trate value = R
+        R = engineering shear strain rate (1/time units)
+      trate value = R
         R = true shear strain rate (1/time units) 
 
  • zero or more keyword/value pairs may be appended @@ -337,10 +336,9 @@ atom positions and velocities to provide a velocity profile that matches the changing box size/shape. Thus atom coordinates should NOT be remapped by fix deform, but velocities SHOULD be when atoms cross -periodic boundaries, since when atoms cross periodic boundaries since -that is consistent with maintaining the velocity profile created by -fix nvt/sllod. LAMMPS will warn you if this settings is not -consistent. +periodic boundaries, since that is consistent with maintaining the +velocity profile already created by fix nvt/sllod. LAMMPS will warn +you if the remap setting is not consistent with fix nvt/sllod.

    The units keyword determines the meaning of the distance units used to define various arguments. A box value selects standard distance diff -Naur lammps-20Jul07/doc/fix_deform.txt lammps-3Aug07/doc/fix_deform.txt --- lammps-20Jul07/doc/fix_deform.txt 2007-07-09 08:09:39.000000000 -0600 +++ lammps-3Aug07/doc/fix_deform.txt 2007-08-03 13:53:18.000000000 -0600 @@ -43,7 +43,7 @@ V = change tilt factor at this velocity (distance/time units), effectively an engineering shear strain rate {erate} value = R - R = engineering shear strain rate (1/time units) :pre + R = engineering shear strain rate (1/time units) {trate} value = R R = true shear strain rate (1/time units) :pre @@ -326,10 +326,9 @@ atom positions and velocities to provide a velocity profile that matches the changing box size/shape. Thus atom coordinates should NOT be remapped by fix deform, but velocities SHOULD be when atoms cross -periodic boundaries, since when atoms cross periodic boundaries since -that is consistent with maintaining the velocity profile created by -fix nvt/sllod. LAMMPS will warn you if this settings is not -consistent. +periodic boundaries, since that is consistent with maintaining the +velocity profile already created by fix nvt/sllod. LAMMPS will warn +you if the {remap} setting is not consistent with fix nvt/sllod. The {units} keyword determines the meaning of the distance units used to define various arguments. A {box} value selects standard distance diff -Naur lammps-20Jul07/src/ASPHERE/fix_npt_asphere.cpp lammps-3Aug07/src/ASPHERE/fix_npt_asphere.cpp --- lammps-20Jul07/src/ASPHERE/fix_npt_asphere.cpp 2007-07-03 13:52:58.000000000 -0600 +++ lammps-3Aug07/src/ASPHERE/fix_npt_asphere.cpp 2007-08-03 13:40:00.000000000 -0600 @@ -108,9 +108,9 @@ } } - // rescale simulation box and all owned atoms by 1/2 step + // remap simulation box and all owned atoms by 1/2 step - box_dilate(0); + remap(0); // x update by full step only for atoms in NPT group @@ -137,10 +137,10 @@ } } - // rescale simulation box and all owned atoms by 1/2 step + // remap simulation box and all owned atoms by 1/2 step // redo KSpace coeffs since volume has changed - box_dilate(0); + remap(0); if (kspace_flag) force->kspace->setup(); } diff -Naur lammps-20Jul07/src/MAKE/Makefile.liberty lammps-3Aug07/src/MAKE/Makefile.liberty --- lammps-20Jul07/src/MAKE/Makefile.liberty 2006-09-27 13:51:33.000000000 -0600 +++ lammps-3Aug07/src/MAKE/Makefile.liberty 2007-07-20 09:21:54.000000000 -0600 @@ -5,13 +5,24 @@ # System-specific settings -FFTW = /apps/libraries/fftw/icc +# Note: this Makefile builds LAMMPS according to what modules you've loaded +# by default this is openmpi MPI +# if you want to build with mpich MPI instead, do the following: +# (1) put these lines in your .cshrc +# module switch mpi mpi/mpich-mx-1.2.7..4_intel-10.0-f025-c025 +# module load libraries/fftw-2.1.5_mpich-mx-1.2.7..4_intel-9.1-f040-c046 +# (2) your qsub environment also needs these same modules +# this will occur if you use qsub -V to inherit from your login shell + +#FFTW = /home/sjplimp/fftw/fftw_liberty CC = mpiCC -CCFLAGS = -O -DFFT_FFTW -I${FFTW}/fftw +#CCFLAGS = -O -DFFT_FFTW -I${FFTW}/fftw +CCFLAGS = -O -DFFT_FFTW -I${FFTW_INCLUDE} DEPFLAGS = -M LINK = mpiCC -LINKFLAGS = -O -L${FFTW}/fftw/.libs +#LINKFLAGS = -O -L${FFTW}/fftw/.libs +LINKFLAGS = -O -L${FFTW_LIB} USRLIB = -lfftw SYSLIB = -lstdc++ -lm SIZE = size diff -Naur lammps-20Jul07/src/POEMS/fix_poems.cpp lammps-3Aug07/src/POEMS/fix_poems.cpp --- lammps-20Jul07/src/POEMS/fix_poems.cpp 2007-06-21 12:58:43.000000000 -0600 +++ lammps-3Aug07/src/POEMS/fix_poems.cpp 2007-08-03 13:41:05.000000000 -0600 @@ -856,24 +856,17 @@ } /* ---------------------------------------------------------------------- - adjust xcm of each cluster due to box dilation in idim - called by various fixes that change box + adjust xcm of each cluster due to box deformation + called by various fixes that change box size/shape + flag = 0/1 means map from box to lamda coords or vice versa NOTE: cannot do this by changing xcm of each body in cluster - or even 1st body in cluster since POEMS does not observe xcm - but only sets xcm + or even 1st body in cluster + b/c POEMS does not see xcm but only sets xcm so dilation needs to be coordinated with POEMS this routine does nothing for now ------------------------------------------------------------------------- */ -void FixPOEMS::dilate(int idim, - double oldlo, double oldhi, double newlo, double newhi) -{ - // double ratio; - // for (int ibody = 0; ibody < nbody; ibody++) { - // ratio = (xcm[ibody][idim] - oldlo) / (oldhi - oldlo); - // xcm[ibody][idim] = newlo + ratio*(newhi - newlo); - // } -} +void FixPOEMS::deform(int flag) {} /* ---------------------------------------------------------------------- */ diff -Naur lammps-20Jul07/src/POEMS/fix_poems.h lammps-3Aug07/src/POEMS/fix_poems.h --- lammps-20Jul07/src/POEMS/fix_poems.h 2007-01-29 17:22:05.000000000 -0700 +++ lammps-3Aug07/src/POEMS/fix_poems.h 2007-08-03 13:41:05.000000000 -0600 @@ -39,7 +39,7 @@ int unpack_exchange(int, double *); int dof(int); - void dilate(int, double, double, double, double); + void deform(int); private: int me; diff -Naur lammps-20Jul07/src/dump_atom.cpp lammps-3Aug07/src/dump_atom.cpp --- lammps-20Jul07/src/dump_atom.cpp 2007-03-07 17:54:02.000000000 -0700 +++ lammps-3Aug07/src/dump_atom.cpp 2007-08-03 13:39:46.000000000 -0600 @@ -40,6 +40,9 @@ void DumpAtom::init() { + // this error check cannot occur in constructor since scale_flag = 1 + // is the default and it can be changed by dump_modify + if (scale_flag && domain->triclinic) error->all("Cannot dump scaled coords with triclinic box"); diff -Naur lammps-20Jul07/src/fix.h lammps-3Aug07/src/fix.h --- lammps-20Jul07/src/fix.h 2007-06-27 14:57:28.000000000 -0600 +++ lammps-3Aug07/src/fix.h 2007-08-03 13:36:35.000000000 -0600 @@ -93,7 +93,7 @@ virtual double thermo(int) {return 0.0;} virtual int dof(int) {return 0;} - virtual void dilate(int, double, double, double, double) {} + virtual void deform(int) {} virtual int modify_param(int, char **) {return 0;} diff -Naur lammps-20Jul07/src/fix_nph.cpp lammps-3Aug07/src/fix_nph.cpp --- lammps-20Jul07/src/fix_nph.cpp 2007-07-03 13:52:58.000000000 -0600 +++ lammps-3Aug07/src/fix_nph.cpp 2007-08-03 13:40:30.000000000 -0600 @@ -22,7 +22,7 @@ #include "atom.h" #include "force.h" #include "comm.h" -#include "output.h" +#include "domain.h" #include "modify.h" #include "compute.h" #include "kspace.h" @@ -110,7 +110,7 @@ // process extra keywords drag = 0.0; - dilate_partial = 0; + allremap = 1; int iarg; if (press_couple == 0) iarg = 7; @@ -123,8 +123,8 @@ iarg += 2; } else if (strcmp(arg[iarg],"dilate") == 0) { if (iarg+2 > narg) error->all("Illegal fix nph command"); - if (strcmp(arg[iarg+1],"all") == 0) dilate_partial = 0; - else if (strcmp(arg[iarg+1],"partial") == 0) dilate_partial = 1; + if (strcmp(arg[iarg+1],"all") == 0) allremap = 1; + else if (strcmp(arg[iarg+1],"partial") == 0) allremap = 0; else error->all("Illegal fix nph command"); iarg += 2; } else error->all("Illegal fix nph command"); @@ -283,7 +283,7 @@ step_respa = ((Respa *) update->integrate)->step; } - // detect if any rigid fixes exist so rigid bodies move when box is dilated + // detect if any rigid fixes exist so rigid bodies move when box is remapped // rfix[] = indices to each fix rigid delete [] rfix; @@ -306,7 +306,7 @@ void FixNPH::setup() { - p_target[0] = p_start[0]; // used by thermo() + p_target[0] = p_start[0]; // needed by thermo() method p_target[1] = p_start[1]; p_target[2] = p_start[2]; @@ -371,9 +371,9 @@ } } - // rescale simulation box and all owned atoms by 1/2 step + // remap simulation box and all owned atoms by 1/2 step - box_dilate(0); + remap(0); // x update by full step only for atoms in NPH group @@ -385,10 +385,10 @@ } } - // rescale simulation box and all owned atoms by 1/2 step + // remap simulation box and all owned atoms by 1/2 step // redo KSpace coeffs since volume has changed - box_dilate(0); + remap(0); if (kspace_flag) force->kspace->setup(); } @@ -452,11 +452,11 @@ void FixNPH::initial_integrate_respa(int ilevel, int flag) { // if flag = 1, then is 2nd call at outermost level from rRESPA - // perform 2nd half of box dilate on own + ghost atoms and return + // perform 2nd half of box remap on own + ghost atoms and return // redo KSpace coeffs since volume has changed if (flag == 1) { - box_dilate(1); + remap(1); if (kspace_flag) force->kspace->setup(); return; } @@ -478,7 +478,7 @@ int *mask = atom->mask; int nlocal = atom->nlocal; - // outermost level - update omega_dot, apply to v, dilate box + // outermost level - update omega_dot, apply to v, remap box // all other levels - NVE update of v // x,v updates only performed for atoms in NPH group @@ -516,9 +516,9 @@ } } - // rescale simulation box and all owned atoms by 1/2 step + // remap simulation box and all owned atoms by 1/2 step - box_dilate(0); + remap(0); } else { @@ -613,89 +613,62 @@ } /* ---------------------------------------------------------------------- - dilate the box around center of box + change box size + remap owned or owned+ghost atoms depending on flag + remap all atoms or fix group atoms depending on allremap flag + if rigid bodies exist, scale rigid body centers-of-mass ------------------------------------------------------------------------- */ -void FixNPH::box_dilate(int flag) +void FixNPH::remap(int flag) { int i,n; - - // ctr = geometric center of box in a dimension - // scale owned or owned+ghost atoms depending on flag - // re-define simulation box via xprd/yprd/zprd - // scale atom coords for all atoms or only for fix group atoms - // if fix rigid exists, scale rigid body centers-of-mass - // don't do anything if non-periodic or press style is constant volume + double oldlo,oldhi,ctr; double **x = atom->x; int *mask = atom->mask; if (flag) n = atom->nlocal + atom->nghost; else n = atom->nlocal; - double oldlo,oldhi,ctr; + // convert pertinent atoms and rigid bodies to lamda coords - if (domain->xperiodic && p_flag[0]) { - oldlo = domain->boxlo[0]; - oldhi = domain->boxhi[0]; - ctr = 0.5 * (oldlo + oldhi); - domain->boxlo[0] = (oldlo-ctr)*dilation[0] + ctr; - domain->boxhi[0] = (oldhi-ctr)*dilation[0] + ctr; - domain->prd[0] = domain->xprd = domain->boxhi[0] - domain->boxlo[0]; - if (dilate_partial) { - for (i = 0; i < n; i++) - if (mask[i] & groupbit) - x[i][0] = ctr + (x[i][0]-ctr)*dilation[0]; - } else { - for (i = 0; i < n; i++) - x[i][0] = ctr + (x[i][0]-ctr)*dilation[0]; - } - if (nrigid) - for (i = 0; i < nrigid; i++) - modify->fix[rfix[i]]-> - dilate(0,oldlo,oldhi,domain->boxlo[0],domain->boxhi[0]); - } - - if (domain->yperiodic && p_flag[1]) { - oldlo = domain->boxlo[1]; - oldhi = domain->boxhi[1]; - ctr = 0.5 * (oldlo + oldhi); - domain->boxlo[1] = (oldlo-ctr)*dilation[1] + ctr; - domain->boxhi[1] = (oldhi-ctr)*dilation[1] + ctr; - domain->prd[1] = domain->yprd = domain->boxhi[1] - domain->boxlo[1]; - if (dilate_partial) { - for (i = 0; i < n; i++) - if (mask[i] & groupbit) - x[i][1] = ctr + (x[i][1]-ctr)*dilation[1]; - } else { - for (i = 0; i < n; i++) - x[i][1] = ctr + (x[i][1]-ctr)*dilation[1]; - } - if (nrigid) - for (i = 0; i < nrigid; i++) - modify->fix[rfix[i]]-> - dilate(1,oldlo,oldhi,domain->boxlo[1],domain->boxhi[1]); - } - - if (domain->zperiodic && p_flag[2]) { - oldlo = domain->boxlo[2]; - oldhi = domain->boxhi[2]; - ctr = 0.5 * (oldlo + oldhi); - domain->boxlo[2] = (oldlo-ctr)*dilation[2] + ctr; - domain->boxhi[2] = (oldhi-ctr)*dilation[2] + ctr; - domain->prd[2] = domain->zprd = domain->boxhi[2] - domain->boxlo[2]; - if (dilate_partial) { - for (i = 0; i < n; i++) - if (mask[i] & groupbit) - x[i][2] = ctr + (x[i][2]-ctr)*dilation[2]; - } else { - for (i = 0; i < n; i++) - x[i][2] = ctr + (x[i][2]-ctr)*dilation[2]; + if (allremap) domain->x2lamda(n); + else { + for (i = 0; i < n; i++) + if (mask[i] & groupbit) + domain->x2lamda(x[i],x[i]); + } + + if (nrigid) + for (i = 0; i < nrigid; i++) + modify->fix[rfix[i]]->deform(0); + + // reset global and local box to new size/shape + + for (i = 0; i < 3; i++) { + if (p_flag[i]) { + oldlo = domain->boxlo[i]; + oldhi = domain->boxhi[i]; + ctr = 0.5 * (oldlo + oldhi); + domain->boxlo[i] = (oldlo-ctr)*dilation[i] + ctr; + domain->boxhi[i] = (oldhi-ctr)*dilation[i] + ctr; } - if (nrigid) - for (i = 0; i < nrigid; i++) - modify->fix[rfix[i]]-> - dilate(2,oldlo,oldhi,domain->boxlo[2],domain->boxhi[2]); } + + domain->set_global_box(); + domain->set_local_box(); + + // convert pertinent atoms and rigid bodies back to box coords + + if (allremap) domain->lamda2x(n); + else { + for (i = 0; i < n; i++) + if (mask[i] & groupbit) + domain->lamda2x(x[i],x[i]); + } + + if (nrigid) + for (i = 0; i < nrigid; i++) + modify->fix[rfix[i]]->deform(1); } /* ---------------------------------------------------------------------- diff -Naur lammps-20Jul07/src/fix_nph.h lammps-3Aug07/src/fix_nph.h --- lammps-20Jul07/src/fix_nph.h 2007-07-03 13:52:58.000000000 -0600 +++ lammps-3Aug07/src/fix_nph.h 2007-08-03 13:36:35.000000000 -0600 @@ -40,7 +40,7 @@ double boltz,nktv2p; double vol0,nkt; - int press_couple,dilate_partial; + int press_couple,allremap; int p_flag[3]; // 1 if control P on this dim, 0 if not double p_start[3],p_stop[3]; double p_freq[3],p_target[3]; @@ -60,7 +60,7 @@ int tflag,pflag; void couple(); - void box_dilate(int); + void remap(int); }; } diff -Naur lammps-20Jul07/src/fix_npt.cpp lammps-3Aug07/src/fix_npt.cpp --- lammps-20Jul07/src/fix_npt.cpp 2007-07-03 13:52:58.000000000 -0600 +++ lammps-3Aug07/src/fix_npt.cpp 2007-08-03 13:36:35.000000000 -0600 @@ -22,7 +22,6 @@ #include "atom.h" #include "force.h" #include "comm.h" -#include "output.h" #include "modify.h" #include "compute.h" #include "kspace.h" @@ -117,7 +116,7 @@ // process extra keywords drag = 0.0; - dilate_partial = 0; + allremap = 1; int iarg; if (press_couple == 0) iarg = 10; @@ -130,8 +129,8 @@ iarg += 2; } else if (strcmp(arg[iarg],"dilate") == 0) { if (iarg+2 > narg) error->all("Illegal fix npt command"); - if (strcmp(arg[iarg+1],"all") == 0) dilate_partial = 0; - else if (strcmp(arg[iarg+1],"partial") == 0) dilate_partial = 1; + if (strcmp(arg[iarg+1],"all") == 0) allremap = 1; + else if (strcmp(arg[iarg+1],"partial") == 0) allremap = 0; else error->all("Illegal fix npt command"); iarg += 2; } else error->all("Illegal fix npt command"); @@ -284,7 +283,7 @@ step_respa = ((Respa *) update->integrate)->step; } - // detect if any rigid fixes exist so rigid bodies move when box is dilated + // detect if any rigid fixes exist so rigid bodies move when box is remapped // rfix[] = indices to each fix rigid delete [] rfix; @@ -381,9 +380,9 @@ } } - // rescale simulation box and all owned atoms by 1/2 step + // remap simulation box and all owned atoms by 1/2 step - box_dilate(0); + remap(0); // x update by full step only for atoms in NPT group @@ -395,10 +394,10 @@ } } - // rescale simulation box and all owned atoms by 1/2 step + // remap simulation box and all owned atoms by 1/2 step // redo KSpace coeffs since volume has changed - box_dilate(0); + remap(0); if (kspace_flag) force->kspace->setup(); } @@ -468,11 +467,11 @@ void FixNPT::initial_integrate_respa(int ilevel, int flag) { // if flag = 1, then is 2nd call at outermost level from rRESPA - // perform 2nd half of box dilate on own + ghost atoms and return + // perform 2nd half of box remap on own + ghost atoms and return // redo KSpace coeffs since volume has changed if (flag == 1) { - box_dilate(1); + remap(1); if (kspace_flag) force->kspace->setup(); return; } @@ -494,7 +493,7 @@ int *mask = atom->mask; int nlocal = atom->nlocal; - // outermost level - update eta_dot and omega_dot, apply to v, dilate box + // outermost level - update eta_dot and omega_dot, apply to v, remap box // all other levels - NVE update of v // x,v updates only performed for atoms in NPT group @@ -540,9 +539,9 @@ } } - // rescale simulation box and all owned atoms by 1/2 step + // remap simulation box and all owned atoms by 1/2 step - box_dilate(0); + remap(0); } else { @@ -638,89 +637,62 @@ } /* ---------------------------------------------------------------------- - dilate the box around center of box + change box size + remap owned or owned+ghost atoms depending on flag + remap all atoms or fix group atoms depending on allremap flag + if rigid bodies exist, scale rigid body centers-of-mass ------------------------------------------------------------------------- */ -void FixNPT::box_dilate(int flag) +void FixNPT::remap(int flag) { int i,n; - - // ctr = geometric center of box in a dimension - // scale owned or owned+ghost atoms depending on flag - // re-define simulation box via xprd/yprd/zprd - // scale atom coords for all atoms or only for fix group atoms - // if fix rigid exists, scale rigid body centers-of-mass - // don't do anything if non-periodic or press style is constant volume + double oldlo,oldhi,ctr; double **x = atom->x; int *mask = atom->mask; if (flag) n = atom->nlocal + atom->nghost; else n = atom->nlocal; - double oldlo,oldhi,ctr; + // convert pertinent atoms and rigid bodies to lamda coords - if (domain->xperiodic && p_flag[0]) { - oldlo = domain->boxlo[0]; - oldhi = domain->boxhi[0]; - ctr = 0.5 * (oldlo + oldhi); - domain->boxlo[0] = (oldlo-ctr)*dilation[0] + ctr; - domain->boxhi[0] = (oldhi-ctr)*dilation[0] + ctr; - domain->prd[0] = domain->xprd = domain->boxhi[0] - domain->boxlo[0]; - if (dilate_partial) { - for (i = 0; i < n; i++) - if (mask[i] & groupbit) - x[i][0] = ctr + (x[i][0]-ctr)*dilation[0]; - } else { - for (i = 0; i < n; i++) - x[i][0] = ctr + (x[i][0]-ctr)*dilation[0]; - } - if (nrigid) - for (i = 0; i < nrigid; i++) - modify->fix[rfix[i]]-> - dilate(0,oldlo,oldhi,domain->boxlo[0],domain->boxhi[0]); - } - - if (domain->yperiodic && p_flag[1]) { - oldlo = domain->boxlo[1]; - oldhi = domain->boxhi[1]; - ctr = 0.5 * (oldlo + oldhi); - domain->boxlo[1] = (oldlo-ctr)*dilation[1] + ctr; - domain->boxhi[1] = (oldhi-ctr)*dilation[1] + ctr; - domain->prd[1] = domain->yprd = domain->boxhi[1] - domain->boxlo[1]; - if (dilate_partial) { - for (i = 0; i < n; i++) - if (mask[i] & groupbit) - x[i][1] = ctr + (x[i][1]-ctr)*dilation[1]; - } else { - for (i = 0; i < n; i++) - x[i][1] = ctr + (x[i][1]-ctr)*dilation[1]; - } - if (nrigid) - for (i = 0; i < nrigid; i++) - modify->fix[rfix[i]]-> - dilate(1,oldlo,oldhi,domain->boxlo[1],domain->boxhi[1]); - } - - if (domain->zperiodic && p_flag[2]) { - oldlo = domain->boxlo[2]; - oldhi = domain->boxhi[2]; - ctr = 0.5 * (oldlo + oldhi); - domain->boxlo[2] = (oldlo-ctr)*dilation[2] + ctr; - domain->boxhi[2] = (oldhi-ctr)*dilation[2] + ctr; - domain->prd[2] = domain->zprd = domain->boxhi[2] - domain->boxlo[2]; - if (dilate_partial) { - for (i = 0; i < n; i++) - if (mask[i] & groupbit) - x[i][2] = ctr + (x[i][2]-ctr)*dilation[2]; - } else { - for (i = 0; i < n; i++) - x[i][2] = ctr + (x[i][2]-ctr)*dilation[2]; + if (allremap) domain->x2lamda(n); + else { + for (i = 0; i < n; i++) + if (mask[i] & groupbit) + domain->x2lamda(x[i],x[i]); + } + + if (nrigid) + for (i = 0; i < nrigid; i++) + modify->fix[rfix[i]]->deform(0); + + // reset global and local box to new size/shape + + for (i = 0; i < 3; i++) { + if (p_flag[i]) { + oldlo = domain->boxlo[i]; + oldhi = domain->boxhi[i]; + ctr = 0.5 * (oldlo + oldhi); + domain->boxlo[i] = (oldlo-ctr)*dilation[i] + ctr; + domain->boxhi[i] = (oldhi-ctr)*dilation[i] + ctr; } - if (nrigid) - for (i = 0; i < nrigid; i++) - modify->fix[rfix[i]]-> - dilate(2,oldlo,oldhi,domain->boxlo[2],domain->boxhi[2]); } + + domain->set_global_box(); + domain->set_local_box(); + + // convert pertinent atoms and rigid bodies back to box coords + + if (allremap) domain->lamda2x(n); + else { + for (i = 0; i < n; i++) + if (mask[i] & groupbit) + domain->lamda2x(x[i],x[i]); + } + + if (nrigid) + for (i = 0; i < nrigid; i++) + modify->fix[rfix[i]]->deform(1); } /* ---------------------------------------------------------------------- diff -Naur lammps-20Jul07/src/fix_npt.h lammps-3Aug07/src/fix_npt.h --- lammps-20Jul07/src/fix_npt.h 2007-07-03 13:52:58.000000000 -0600 +++ lammps-3Aug07/src/fix_npt.h 2007-08-03 13:36:35.000000000 -0600 @@ -45,7 +45,7 @@ double t_freq; double f_eta,eta_dot,eta; - int press_couple,dilate_partial; + int press_couple,allremap; int p_flag[3]; // 1 if control P on this dim, 0 if not double p_start[3],p_stop[3]; double p_freq[3],p_target[3]; @@ -65,7 +65,7 @@ int tflag,pflag; void couple(); - void box_dilate(int); + void remap(int); }; } diff -Naur lammps-20Jul07/src/fix_nvt.cpp lammps-3Aug07/src/fix_nvt.cpp --- lammps-20Jul07/src/fix_nvt.cpp 2007-06-22 10:59:17.000000000 -0600 +++ lammps-3Aug07/src/fix_nvt.cpp 2007-08-03 13:43:35.000000000 -0600 @@ -129,7 +129,7 @@ void FixNVT::setup() { - t_target = t_start; // used by thermo() + t_target = t_start; // needed by thermo() method t_current = temperature->compute_scalar(); } diff -Naur lammps-20Jul07/src/fix_rigid.cpp lammps-3Aug07/src/fix_rigid.cpp --- lammps-20Jul07/src/fix_rigid.cpp 2007-06-27 10:45:46.000000000 -0600 +++ lammps-3Aug07/src/fix_rigid.cpp 2007-08-03 13:36:35.000000000 -0600 @@ -847,18 +847,19 @@ } /* ---------------------------------------------------------------------- - adjust xcm of each rigid body due to box dilation in idim - called by various fixes that change box + adjust xcm of each rigid body due to box deformation + called by various fixes that change box size/shape + flag = 0/1 means map from box to lamda coords or vice versa ------------------------------------------------------------------------- */ -void FixRigid::dilate(int idim, - double oldlo, double oldhi, double newlo, double newhi) -{ - double ratio; - for (int ibody = 0; ibody < nbody; ibody++) { - ratio = (xcm[ibody][idim] - oldlo) / (oldhi - oldlo); - xcm[ibody][idim] = newlo + ratio*(newhi - newlo); - } +void FixRigid::deform(int flag) +{ + if (flag == 0) + for (int ibody = 0; ibody < nbody; ibody++) + domain->x2lamda(xcm[ibody],xcm[ibody]); + else + for (int ibody = 0; ibody < nbody; ibody++) + domain->lamda2x(xcm[ibody],xcm[ibody]); } /* ---------------------------------------------------------------------- diff -Naur lammps-20Jul07/src/fix_rigid.h lammps-3Aug07/src/fix_rigid.h --- lammps-20Jul07/src/fix_rigid.h 2007-01-29 17:22:05.000000000 -0700 +++ lammps-3Aug07/src/fix_rigid.h 2007-08-03 13:40:30.000000000 -0600 @@ -37,7 +37,7 @@ int unpack_exchange(int, double *); int dof(int); - void dilate(int, double, double, double, double); + void deform(int); private: double dtv,dtf,dtq; diff -Naur lammps-20Jul07/src/fix_set_force.cpp lammps-3Aug07/src/fix_set_force.cpp --- lammps-20Jul07/src/fix_set_force.cpp 2007-03-22 10:35:57.000000000 -0600 +++ lammps-3Aug07/src/fix_set_force.cpp 2007-08-03 13:38:10.000000000 -0600 @@ -131,7 +131,10 @@ post_force(vflag); } -/* ---------------------------------------------------------------------- */ +/* ---------------------------------------------------------------------- + allow setforce stats to be printed with thermo output + n = 1,2,3 = components of total force on fix group before reset +------------------------------------------------------------------------- */ double FixSetForce::thermo(int n) { diff -Naur lammps-20Jul07/src/input.cpp lammps-3Aug07/src/input.cpp --- lammps-20Jul07/src/input.cpp 2007-06-22 10:59:17.000000000 -0600 +++ lammps-3Aug07/src/input.cpp 2007-08-03 13:37:42.000000000 -0600 @@ -452,8 +452,12 @@ else flag = 0; + // return if command was listed above + if (flag) return 0; + // check if command is added via style.h + if (0) return 0; // dummy line to enable else-if macro expansion #define CommandClass @@ -466,6 +470,8 @@ #include "style.h" #undef CommandClass + // unrecognized command + return -1; } diff -Naur lammps-20Jul07/src/modify.cpp lammps-3Aug07/src/modify.cpp --- lammps-20Jul07/src/modify.cpp 2007-02-20 17:15:29.000000000 -0700 +++ lammps-3Aug07/src/modify.cpp 2007-08-03 13:37:42.000000000 -0600 @@ -234,7 +234,7 @@ /* ---------------------------------------------------------------------- thermo energy call only for relevant fixes - called by Thermo clas + called by Thermo class arg to Fix thermo() is 0, so fix will return its energy contribution ------------------------------------------------------------------------- */ diff -Naur lammps-20Jul07/src/neigh_stencil.cpp lammps-3Aug07/src/neigh_stencil.cpp --- lammps-20Jul07/src/neigh_stencil.cpp 2007-03-26 15:30:26.000000000 -0600 +++ lammps-3Aug07/src/neigh_stencil.cpp 2007-08-03 13:43:14.000000000 -0600 @@ -49,6 +49,8 @@ { int i; + // for 2d, sz = 0 + int nmax = (2*sz+1) * (2*sy+1) * (2*sx+1); if (half) { diff -Naur lammps-20Jul07/src/read_restart.cpp lammps-3Aug07/src/read_restart.cpp --- lammps-20Jul07/src/read_restart.cpp 2007-07-12 07:07:06.000000000 -0600 +++ lammps-3Aug07/src/read_restart.cpp 2007-08-03 13:40:47.000000000 -0600 @@ -269,7 +269,7 @@ } /* ---------------------------------------------------------------------- - search for all files in dir matching infile which contains "*" + search for all files in dir matching infile which contains a "*" replace "*" with latest timestep value to create outfile name if infile also contains "%", need to use "base" when search directory ------------------------------------------------------------------------- */ diff -Naur lammps-20Jul07/src/respa.cpp lammps-3Aug07/src/respa.cpp --- lammps-20Jul07/src/respa.cpp 2007-06-20 07:17:59.000000000 -0600 +++ lammps-3Aug07/src/respa.cpp 2007-08-03 13:42:43.000000000 -0600 @@ -454,8 +454,8 @@ if (ilevel) recurse(ilevel-1); // at outermost level: - // call initial_integrate w/ flag = 1 so fix NPT,NPH perform 2nd dilate - // at correct symmetric position in rRESPA hierarchy + // call initial_integrate w/ flag = 1 so fix NPT,NPH performs + // 2nd box deform at correct symmetric position in rRESPA hierarchy // check on rebuilding neighbor list // at innermost level, communicate // at middle levels, do nothing