diff --git a/src/fix_nh.cpp b/src/fix_nh.cpp
index fcba0fceaf..35803f4fbe 100644
--- a/src/fix_nh.cpp
+++ b/src/fix_nh.cpp
@@ -69,6 +69,10 @@ FixNH::FixNH(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
   deviatoric_flag = 0;
   nreset_h0 = 0;
 
+  // Used by FixNVTSllod to preserve non-default value  
+
+  mtchain_default_flag = 1;
+
   tstat_flag = 0;
   double t_period = 0.0;
 
@@ -220,6 +224,8 @@ FixNH::FixNH(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
     } else if (strcmp(arg[iarg],"tchain") == 0) {
       if (iarg+2 > narg) error->all("Illegal fix nvt/npt/nph command");
       mtchain = atoi(arg[iarg+1]);
+      // Used by FixNVTSllod to preserve non-default value  
+      mtchain_default_flag = 0;
       if (mtchain < 1) error->all("Illegal fix nvt/npt/nph command");
       iarg += 2;
     } else if (strcmp(arg[iarg],"pchain") == 0) {
diff --git a/src/fix_nh.h b/src/fix_nh.h
index 3edcf33302..9ce1ed395f 100644
--- a/src/fix_nh.h
+++ b/src/fix_nh.h
@@ -73,6 +73,7 @@ class FixNH : public Fix {
   double *eta_dotdot;
   double *eta_mass;
   int mtchain;                     // length of chain
+  int mtchain_default_flag;        // 1 = mtchain is default   
                                    
   double *etap;                    // chain thermostat for barostat
   double *etap_dot;
diff --git a/src/fix_nvt_sllod.cpp b/src/fix_nvt_sllod.cpp
index 5740243305..2a4bc99aff 100644
--- a/src/fix_nvt_sllod.cpp
+++ b/src/fix_nvt_sllod.cpp
@@ -41,6 +41,11 @@ FixNVTSllod::FixNVTSllod(LAMMPS *lmp, int narg, char **arg) :
   if (pstat_flag)
     error->all("Pressure control can not be used with fix nvt/sllod");
 
+  // default values
+
+  if (mtchain_default_flag) mtchain = 1;
+
+
   // create a new compute temp style
   // id = fix-ID + temp
 
@@ -110,10 +115,10 @@ void FixNVTSllod::nh_v_temp()
 
   for (int i = 0; i < nlocal; i++) {
     if (mask[i] & groupbit) {
-      temperature->remove_bias(i,v[i]);
       vdelu[0] = h_two[0]*v[i][0] + h_two[5]*v[i][1] + h_two[4]*v[i][2];
       vdelu[1] = h_two[1]*v[i][1] + h_two[3]*v[i][2];
       vdelu[2] = h_two[2]*v[i][2];
+      temperature->remove_bias(i,v[i]);
       v[i][0] = v[i][0]*factor_eta - dthalf*vdelu[0];
       v[i][1] = v[i][1]*factor_eta - dthalf*vdelu[1];
       v[i][2] = v[i][2]*factor_eta - dthalf*vdelu[2];