Added flag to make only fix_nphug use constant eta_mass

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

src/SHOCK/fix_nphug.cpp

@@ -36,6 +36,10 @@ FixNPHug::FixNPHug(LAMMPS *lmp, int narg, char **arg) :
   FixNH(lmp, narg, arg)
 {
 
+  // Prevent eta_mass from being updated every timestep
+
+  eta_mass_flag = 0;
+
   // extend vector of base-class computes
 
   size_vector += 3;

src/fix_nh.cpp

@@ -72,6 +72,7 @@ FixNH::FixNH(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
   mtk_flag = 1;
   deviatoric_flag = 0;
   nreset_h0 = 0;
+  eta_mass_flag = 1;
 
   // turn on tilt factor scaling, whenever applicable
 
@@ -731,6 +732,11 @@ void FixNH::initial_integrate(int vflag)
 
   if (tstat_flag) {
     compute_temp_target();
+    if (eta_mass_flag) { 
+      eta_mass[0] = tdof * boltz * t_target / (t_freq*t_freq);
+      for (int ich = 1; ich < mtchain; ich++)
+	eta_mass[ich] = boltz * t_target / (t_freq*t_freq);
+    }
     nhc_temp_integrate();
   }
 
@@ -827,6 +833,11 @@ void FixNH::initial_integrate_respa(int vflag, int ilevel, int iloop)
 
     if (tstat_flag) {
       compute_temp_target();
+      if (eta_mass_flag) { 
+	eta_mass[0] = tdof * boltz * t_target / (t_freq*t_freq);
+	for (int ich = 1; ich < mtchain; ich++)
+	  eta_mass[ich] = boltz * t_target / (t_freq*t_freq);
+      }
       nhc_temp_integrate();
     }
 
@@ -1570,48 +1581,11 @@ void FixNH::reset_dt()
   
   if (strstr(update->integrate_style,"respa"))
     dto = 0.5*step_respa[0];
-  
-  p_freq_max = 0.0;
-  if (pstat_flag) {
-    p_freq_max = MAX(p_freq[0],p_freq[1]);
-    p_freq_max = MAX(p_freq_max,p_freq[2]);
-    if (pstyle == TRICLINIC) {
-      p_freq_max = MAX(p_freq_max,p_freq[3]);
-      p_freq_max = MAX(p_freq_max,p_freq[4]);
-      p_freq_max = MAX(p_freq_max,p_freq[5]);
-    }
-
-    double kt = boltz * t_target;
-    double nkt = atom->natoms * kt;
-
-    for (int i = 0; i < 3; i++)
-      if (p_flag[i])
-	omega_mass[i] = nkt/(p_freq[i]*p_freq[i]);
-
-    if (pstyle == TRICLINIC) {
-      for (int i = 3; i < 6; i++)
-	if (p_flag[i]) omega_mass[i] = nkt/(p_freq[i]*p_freq[i]);
-    }
-
-  // masses and initial forces on barostat thermostat variables
-
-    if (mpchain) {
-      etap_mass[0] = boltz * t_target / (p_freq_max*p_freq_max);
-      for (int ich = 1; ich < mpchain; ich++)
-	etap_mass[ich] = boltz * t_target / (p_freq_max*p_freq_max);
-      for (int ich = 1; ich < mpchain; ich++)
-	etap_dotdot[ich] = 
-	  (etap_mass[ich-1]*etap_dot[ich-1]*etap_dot[ich-1] -
-	   boltz * t_target) / etap_mass[ich];
-    }
 
+  if (pstat_flag)
     pdrag_factor = 1.0 - (update->dt * p_freq_max * drag / nc_pchain);
-  }
-
+  
   if (tstat_flag)
-    eta_mass[0] = tdof * boltz * t_target / (t_freq*t_freq);
-    for (int ich = 1; ich < mtchain; ich++)
-      eta_mass[ich] = boltz * t_target / (t_freq*t_freq);
     tdrag_factor = 1.0 - (update->dt * t_freq * drag / nc_tchain);
 }
 

src/fix_nh.h

@@ -104,6 +104,8 @@ class FixNH : public Fix {
 
   double mtk_term1,mtk_term2;      // Martyna-Tobias-Klein corrections
 
+  int eta_mass_flag;               // 1 if eta_mass updated, 0 if not.
+
   int scaleyz;                     // 1 if yz scaled with lz 
   int scalexz;                     // 1 if xz scaled with lz 
   int scalexy;                     // 1 if xy scaled with ly