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

src/Make.py

@@ -590,7 +590,7 @@ Syntax: Make.py switch args ...
   switch for actions:
     -a lib-all, lib-dir, clean, file, exe or machine
     list one or more actions, in any order
-    machine is a Makefile.machine suffix, must be last if used
+    machine is a Makefile.machine suffix
   one-letter switches:
     -d (dir), -j (jmake), -m (makefile), -o (output),
     -p (packages), -r (redo), -s (settings), -v (verbose)

src/RIGID/fix_shake.cpp

@@ -202,6 +202,10 @@ FixShake::FixShake(LAMMPS *lmp, int narg, char **arg) :
     a_min_all = new double[na];
   }
 
+  // SHAKE vs RATTLE
+
+  rattle = 0;
+
   // identify all SHAKE clusters
 
   find_clusters();
@@ -428,12 +432,13 @@ void FixShake::setup(int vflag)
   // half timestep constraint on pre-step, full timestep thereafter
 
   if (strstr(update->integrate_style,"verlet")) {
-    dtv = update->dt;
-    dtfsq = 0.5 * update->dt * update->dt * force->ftm2v;
-    post_force(vflag);
-    dtfsq = update->dt * update->dt * force->ftm2v;
+    dtv     = update->dt;
+    dtfsq   = 0.5 * update->dt * update->dt * force->ftm2v;
+    FixShake::post_force(vflag);
+    if (!rattle) dtfsq = update->dt * update->dt * force->ftm2v;
+
   } else {
-    dtv = step_respa[0];
+    dtv = step_respa[0]; 
     dtf_innerhalf = 0.5 * step_respa[0] * force->ftm2v;
     dtf_inner = dtf_innerhalf;
 
@@ -441,11 +446,10 @@ void FixShake::setup(int vflag)
 
     for (int ilevel = 0; ilevel < nlevels_respa; ilevel++) {
       ((Respa *) update->integrate)->copy_flevel_f(ilevel);
-      post_force_respa(vflag,ilevel,loop_respa[ilevel]-1);
+      FixShake::post_force_respa(vflag,ilevel,loop_respa[ilevel]-1);
       ((Respa *) update->integrate)->copy_f_flevel(ilevel);
     }
-
-    dtf_inner = step_respa[0] * force->ftm2v;
+    if (!rattle) dtf_inner = step_respa[0] * force->ftm2v;
   }
 }
 
@@ -652,8 +656,11 @@ void FixShake::find_clusters()
   tagint tagprev;
   double massone;
   tagint *buf;
-
-  if (me == 0 && screen) fprintf(screen,"Finding SHAKE clusters ...\n");
+  
+  if (me == 0 && screen) {
+    if (!rattle) fprintf(screen,"Finding SHAKE clusters ...\n");
+    else fprintf(screen,"Finding RATTLE clusters ...\n");
+  }
 
   atommols = atom->avec->onemols;
 
@@ -1625,7 +1632,7 @@ void FixShake::shake3(int m)
 
 void FixShake::shake4(int m)
 {
-  int nlist,list[4];
+ int nlist,list[4];
   double v[6];
   double invmass0,invmass1,invmass2,invmass3;
 
@@ -2636,11 +2643,13 @@ void FixShake::reset_dt()
 {
   if (strstr(update->integrate_style,"verlet")) {
     dtv = update->dt;
-    dtfsq = update->dt * update->dt * force->ftm2v;
+    if (rattle) dtfsq   = 0.5 * update->dt * update->dt * force->ftm2v;
+    else dtfsq = update->dt * update->dt * force->ftm2v;
   } else {
     dtv = step_respa[0];
     dtf_innerhalf = 0.5 * step_respa[0] * force->ftm2v;
-    dtf_inner = step_respa[0] * force->ftm2v;
+    if (rattle) dtf_inner = dtf_innerhalf;
+    else dtf_inner = step_respa[0] * force->ftm2v;
   }
 }
 

src/RIGID/fix_shake.h

@@ -27,32 +27,33 @@ namespace LAMMPS_NS {
 class FixShake : public Fix {
  public:
   FixShake(class LAMMPS *, int, char **);
-  ~FixShake();
-  int setmask();
-  void init();
+  virtual ~FixShake();
+  virtual int setmask();
+  virtual void init();
   void setup(int);
   void pre_neighbor();
-  void post_force(int);
-  void post_force_respa(int, int, int);
+  virtual void post_force(int);
+  virtual void post_force_respa(int, int, int);
 
-  double memory_usage();
-  void grow_arrays(int);
-  void copy_arrays(int, int, int);
+  virtual double memory_usage();
+  virtual void grow_arrays(int);
+  virtual void copy_arrays(int, int, int);
   void set_arrays(int);
-  void update_arrays(int, int);
+  virtual void update_arrays(int, int);
   void set_molecule(int, tagint, int, double *, double *, double *);
 
-  int pack_exchange(int, double *);
-  int unpack_exchange(int, double *);
-  int pack_forward_comm(int, int *, double *, int, int *);
-  void unpack_forward_comm(int, int, double *);
+  virtual int pack_exchange(int, double *);
+  virtual int unpack_exchange(int, double *);
+  virtual int pack_forward_comm(int, int *, double *, int, int *);
+  virtual void unpack_forward_comm(int, int, double *);
 
   int dof(int);
-  void reset_dt();
+  virtual void reset_dt();
   void *extract(const char *, int &);
 
- private:
+ protected:
   int me,nprocs;
+  int rattle;                            // 0 = SHAKE, 1 = RATTLE
   double tolerance;                      // SHAKE tolerance
   int max_iter;                          // max # of SHAKE iterations
   int output_every;                      // SHAKE stat output every so often

src/fix_respa.cpp

@@ -14,9 +14,8 @@
 #include "stdlib.h"
 #include "fix_respa.h"
 #include "atom.h"
-#include "memory.h"
-#include "error.h"
 #include "force.h"
+#include "memory.h"
 
 using namespace LAMMPS_NS;
 using namespace FixConst;

src/fix_respa.h

@@ -27,7 +27,7 @@ namespace LAMMPS_NS {
 class FixRespa : public Fix {
   friend class Respa;
   friend class FixShake;
-  friend class FixShake2;
+  friend class FixRattle;
 
  public:
   FixRespa(class LAMMPS *, int, char **);
@@ -50,6 +50,7 @@ class FixRespa : public Fix {
 
 #endif
 #endif
+
 /* ERROR/WARNING messages:
 
 */