Merge pull request #285 from akohlmey/fix-ipi-update

update for fix ipi from michele ceriotti

doc/src/fix_ipi.txt

@@ -10,18 +10,19 @@ fix ipi command :h3
 
 [Syntax:]
 
-fix ID group-ID ipi address port \[unix\] :pre
+fix ID group-ID ipi address port \[unix\] \[reset\] :pre
 
 ID, group-ID are documented in "fix"_fix.html command
 ipi = style name of this fix command
 address = internet address (FQDN or IP), or UNIX socket name
 port = port number (ignored for UNIX sockets)
-optional keyword = {unix}, if present uses a unix socket :ul
+optional keyword = {unix}, if present uses a unix socket
+optional keyword = {reset}, if present reset electrostatics at each call :ul
 
 [Examples:]
 
 fix 1 all ipi my.server.com 12345
-fix 1 all ipi mysocket 666 unix
+fix 1 all ipi mysocket 666 unix reset
 
 [Description:]
 
@@ -57,6 +58,15 @@ input are listed in the same order as in the data file of LAMMPS. The
 initial configuration is ignored, as it will be substituted with the
 coordinates received from i-PI before forces are ever evaluated.
 
+A note of caution when using potentials that contain long-range 
+electrostatics, or that contain parameters that depend on box size:
+all of these options will be initialized based on the cell size in the
+LAMMPS-side initial configuration and kept constant during the run. 
+This is required to e.g. obtain reproducible and conserved forces. 
+If the cell varies too wildly, it may be advisable to reinitialize 
+these interactions at each call. This behavior can be requested by 
+setting the {reset} switch. 
+
 [Restart, fix_modify, output, run start/stop, minimize info:]
 
 There is no restart information associated with this fix, since all

src/USER-MISC/fix_ipi.cpp

@@ -30,6 +30,7 @@
 #include "compute.h"
 #include "comm.h"
 #include "neighbor.h"
+#include "irregular.h"
 #include "domain.h"
 #include "compute_pressure.h"
 #include <errno.h>
@@ -171,7 +172,7 @@ static void readbuffer(int sockfd, char *data, int len, Error* error)
 /* ---------------------------------------------------------------------- */
 
 FixIPI::FixIPI(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg), irregular(NULL)
 {
   /* format for fix:
    *  fix  num  group_id ipi host port [unix]
@@ -191,8 +192,11 @@ FixIPI::FixIPI(LAMMPS *lmp, int narg, char **arg) :
   host = strdup(arg[3]);
   port = force->inumeric(FLERR,arg[4]);
 
-  inet   = ((narg > 5) && (strcmp(arg[5],"unix") ==0) ) ? 0 : 1;
+  inet   = ((narg > 5) && (strcmp(arg[5],"unix") == 0) ) ? 0 : 1;
   master = (comm->me==0) ? 1 : 0;
+  // check if forces should be reinitialized and set flag
+  reset_flag = ((narg > 6 && (strcmp(arg[5],"reset") == 0 )) || ((narg > 5) && (strcmp(arg[5],"reset") == 0)) ) ? 1 : 0;
+
   hasdata = bsize = 0;
 
   // creates a temperature compute for all atoms
@@ -212,6 +216,9 @@ FixIPI::FixIPI(LAMMPS *lmp, int narg, char **arg) :
   newarg[4] = (char *) "virial";
   modify->add_compute(5,newarg);
   delete [] newarg;
+
+  // create instance of Irregular class
+  irregular = new Irregular(lmp);
 }
 
 /* ---------------------------------------------------------------------- */
@@ -222,6 +229,7 @@ FixIPI::~FixIPI()
   free(host);
   modify->delete_compute("IPI_TEMP");
   modify->delete_compute("IPI_PRESS");
+  delete irregular;
 }
 
 
@@ -331,10 +339,12 @@ void FixIPI::initial_integrate(int vflag)
   domain->xz = cellh[2]*posconv;
   domain->yz = cellh[5]*posconv;
 
-  // signal that the box has (or may have) changed
+  // do error checks on simulation box and set small for triclinic boxes
+  domain->set_initial_box();
+  // reset global and local box using the new box dimensions
   domain->reset_box();
+  // signal that the box has (or may have) changed
   domain->box_change = 1;
-  if (kspace_flag) force->kspace->setup();
 
   // picks local atoms from the buffer
   double **x = atom->x;
@@ -349,6 +359,36 @@ void FixIPI::initial_integrate(int vflag)
     }
   }
 
+  // insure atoms are in current box & update box via shrink-wrap
+  // has to be be done before invoking Irregular::migrate_atoms()
+  //   since it requires atoms be inside simulation box
+
+  if (domain->triclinic) domain->x2lamda(atom->nlocal);
+  domain->pbc();
+  domain->reset_box();
+  if (domain->triclinic) domain->lamda2x(atom->nlocal);
+
+  // move atoms to new processors via irregular()
+  // only needed if migrate_check() says an atom moves to far
+  if (domain->triclinic) domain->x2lamda(atom->nlocal);
+  if (irregular->migrate_check()) irregular->migrate_atoms();
+  if (domain->triclinic) domain->lamda2x(atom->nlocal);
+
+  // check if kspace solver is used
+  if (reset_flag && kspace_flag) {
+    // reset kspace, pair, angles, ... b/c simulation box might have changed.
+    //   kspace->setup() is in some cases not enough since, e.g., g_ewald needs
+    //   to be reestimated due to changes in box dimensions.
+    force->init();
+    // setup_grid() is necessary for pppm since init() is not calling
+    //   setup() nor setup_grid() upon calling init().
+    if (force->kspace->pppmflag) force->kspace->setup_grid();
+    // other kspace styles might need too another setup()?
+  } else if (!reset_flag && kspace_flag) {
+    // original version
+    force->kspace->setup();
+  }
+
   // compute PE. makes sure that it will be evaluated at next step
   modify->compute[modify->find_compute("thermo_pe")]->invoked_scalar = -1;
   modify->addstep_compute_all(update->ntimestep+1);

src/USER-MISC/fix_ipi.h

@@ -37,6 +37,10 @@ class FixIPI : public Fix {
   char *host; int port; int inet, master, hasdata;
   int ipisock, me; double *buffer; long bsize;
   int kspace_flag;
+  int reset_flag;
+
+ private:
+  class Irregular *irregular;
 };
 
 }