diff --git a/src/MC/fix_atom_swap.cpp b/src/MC/fix_atom_swap.cpp
index 610f1f099e..b3420e423e 100644
--- a/src/MC/fix_atom_swap.cpp
+++ b/src/MC/fix_atom_swap.cpp
@@ -144,7 +144,7 @@ void FixAtomSwap::options(int narg, char **arg)
   if (narg < 0) error->all(FLERR,"Illegal fix atom/swap command");
 
   regionflag = 0;
-  conserve_ke_flag = 1;
+  ke_flag = 1;
   semi_grand_flag = 0;
   nswaptypes = 0;
   nmutypes = 0;
@@ -162,7 +162,7 @@ void FixAtomSwap::options(int narg, char **arg)
       iarg += 2;
     } else if (strcmp(arg[iarg],"ke") == 0) {
       if (iarg+2 > narg) error->all(FLERR,"Illegal fix atom/swap command");
-      conserve_ke_flag = utils::logical(FLERR,arg[iarg+1],false,lmp);
+      ke_flag = utils::logical(FLERR,arg[iarg+1],false,lmp);
       iarg += 2;
     } else if (strcmp(arg[iarg],"semi-grand") == 0) {
       if (iarg+2 > narg) error->all(FLERR,"Illegal fix atom/swap command");
@@ -303,6 +303,8 @@ void FixAtomSwap::pre_exchange()
 
   if (next_reneighbor != update->ntimestep) return;
 
+  // insure current system is ready to compute energy
+
   if (domain->triclinic) domain->x2lamda(atom->nlocal);
   domain->pbc();
   comm->exchange();
@@ -311,8 +313,13 @@ void FixAtomSwap::pre_exchange()
   if (modify->n_pre_neighbor) modify->pre_neighbor();
   neighbor->build(1);
 
+  // energy_stored = energy of current state
+  // will be updated after accepted swaps
+
   energy_stored = energy_full();
 
+  // attempt Ncycle atom swaps
+
   int nsuccess = 0;
   if (semi_grand_flag) {
     update_semi_grand_atoms_list();
@@ -322,23 +329,30 @@ void FixAtomSwap::pre_exchange()
     for (int i = 0; i < ncycles; i++) nsuccess += attempt_swap();
   }
 
+  // udpate MC stats
+
   nswap_attempts += ncycles;
   nswap_successes += nsuccess;
 
-  energy_full();
   next_reneighbor = update->ntimestep + nevery;
 }
 
 /* ----------------------------------------------------------------------
-Note: atom charges are assumed equal and so are not updated
+   attempt a semd-grand swap of a single atom
+   compare before/after energy and accept/reject the swap
+   NOTE: atom charges are assumed equal and so are not updated
 ------------------------------------------------------------------------- */
 
 int FixAtomSwap::attempt_semi_grand()
 {
   if (nswap == 0) return 0;
 
+  // pre-swap energy
+
   double energy_before = energy_stored;
 
+  // pick a random atom and perform swap
+
   int itype,jtype,jswaptype;
   int i = pick_semi_grand_atom();
   if (i >= 0) {
@@ -352,9 +366,12 @@ int FixAtomSwap::attempt_semi_grand()
     atom->type[i] = jtype;
   }
 
+  // if unequal_cutoffs, call comm->borders() and rebuild neighbor list
+  // else communicate ghost atoms
+  // call to comm->exchange() is a no-op but clears ghost atoms
+
   if (unequal_cutoffs) {
     if (domain->triclinic) domain->x2lamda(atom->nlocal);
-    domain->pbc();
     comm->exchange();
     comm->borders();
     if (domain->triclinic) domain->lamda2x(atom->nlocal+atom->nghost);
@@ -364,6 +381,8 @@ int FixAtomSwap::attempt_semi_grand()
     comm->forward_comm_fix(this);
   }
 
+  // post-swap energy
+
   if (force->kspace) force->kspace->qsum_qsq();
   double energy_after = energy_full();
 
@@ -376,10 +395,12 @@ int FixAtomSwap::attempt_semi_grand()
   int success_all = 0;
   MPI_Allreduce(&success,&success_all,1,MPI_INT,MPI_MAX,world);
 
+  // swap accepted, return 1
+
   if (success_all) {
     update_semi_grand_atoms_list();
     energy_stored = energy_after;
-    if (conserve_ke_flag) {
+    if (ke_flag) {
       if (i >= 0) {
         atom->v[i][0] *= sqrt_mass_ratio[itype][jtype];
         atom->v[i][1] *= sqrt_mass_ratio[itype][jtype];
@@ -387,38 +408,35 @@ int FixAtomSwap::attempt_semi_grand()
       }
     }
     return 1;
-  } else {
-    if (i >= 0) {
-      atom->type[i] = itype;
-    }
-    if (force->kspace) force->kspace->qsum_qsq();
-    energy_stored = energy_before;
-
-    if (unequal_cutoffs) {
-      if (domain->triclinic) domain->x2lamda(atom->nlocal);
-      domain->pbc();
-      comm->exchange();
-      comm->borders();
-      if (domain->triclinic) domain->lamda2x(atom->nlocal+atom->nghost);
-      if (modify->n_pre_neighbor) modify->pre_neighbor();
-      neighbor->build(1);
-    } else {
-      comm->forward_comm_fix(this);
-    }
   }
+
+  // swap not accepted, return 0
+  // restore the swapped atom
+  // do not need to re-call comm->borders() and rebuild neighbor list
+  //   since will be done on next cycle or in Verlet when this fix finishes
+
+  if (i >= 0) atom->type[i] = itype;
+  if (force->kspace) force->kspace->qsum_qsq();
+
   return 0;
 }
 
-
 /* ----------------------------------------------------------------------
+   attempt a swap of a pair of atoms
+   compare before/after energy and accept/reject the swap
 ------------------------------------------------------------------------- */
 
 int FixAtomSwap::attempt_swap()
 {
   if ((niswap == 0) || (njswap == 0)) return 0;
 
+  // pre-swap energy
+
   double energy_before = energy_stored;
 
+  // pick a random pair of atoms
+  // swap their properties
+
   int i = pick_i_swap_atom();
   int j = pick_j_swap_atom();
   int itype = type_list[0];
@@ -433,6 +451,10 @@ int FixAtomSwap::attempt_swap()
     if (atom->q_flag) atom->q[j] = qtype[0];
   }
 
+  // if unequal_cutoffs, call comm->borders() and rebuild neighbor list
+  // else communicate ghost atoms
+  // call to comm->exchange() is a no-op but clears ghost atoms
+
   if (unequal_cutoffs) {
     if (domain->triclinic) domain->x2lamda(atom->nlocal);
     domain->pbc();
@@ -445,13 +467,17 @@ int FixAtomSwap::attempt_swap()
     comm->forward_comm_fix(this);
   }
 
+  // post-swap energy
+
   double energy_after = energy_full();
 
+  // swap accepted, return 1
+  // if ke_flag, rescale atom velocities
+
   if (random_equal->uniform() <
       exp(beta*(energy_before - energy_after))) {
     update_swap_atoms_list();
-    energy_stored = energy_after;
-    if (conserve_ke_flag) {
+    if (ke_flag) {
       if (i >= 0) {
         atom->v[i][0] *= sqrt_mass_ratio[itype][jtype];
         atom->v[i][1] *= sqrt_mass_ratio[itype][jtype];
@@ -463,30 +489,24 @@ int FixAtomSwap::attempt_swap()
         atom->v[j][2] *= sqrt_mass_ratio[jtype][itype];
       }
     }
+    energy_stored = energy_after;
     return 1;
-  } else {
-    if (i >= 0) {
-      atom->type[i] =  type_list[0];
-      if (atom->q_flag) atom->q[i] = qtype[0];
-    }
-    if (j >= 0) {
-      atom->type[j] =  type_list[1];
-      if (atom->q_flag) atom->q[j] = qtype[1];
-    }
-    energy_stored = energy_before;
-
-    if (unequal_cutoffs) {
-      if (domain->triclinic) domain->x2lamda(atom->nlocal);
-      domain->pbc();
-      comm->exchange();
-      comm->borders();
-      if (domain->triclinic) domain->lamda2x(atom->nlocal+atom->nghost);
-      if (modify->n_pre_neighbor) modify->pre_neighbor();
-      neighbor->build(1);
-    } else {
-      comm->forward_comm_fix(this);
-    }
   }
+
+  // swap not accepted, return 0
+  // restore the swapped itype & jtype atoms
+  // do not need to re-call comm->borders() and rebuild neighbor list
+  //   since will be done on next cycle or in Verlet when this fix finishes
+
+  if (i >= 0) {
+    atom->type[i] = type_list[0];
+    if (atom->q_flag) atom->q[i] = qtype[0];
+  }
+  if (j >= 0) {
+    atom->type[j] = type_list[1];
+    if (atom->q_flag) atom->q[j] = qtype[1];
+  }
+
   return 0;
 }
 
@@ -499,7 +519,6 @@ double FixAtomSwap::energy_full()
   int eflag = 1;
   int vflag = 0;
 
-  if (modify->n_pre_neighbor) modify->pre_neighbor();
   if (modify->n_pre_force) modify->pre_force(vflag);
 
   if (force->pair) force->pair->compute(eflag,vflag);
diff --git a/src/MC/fix_atom_swap.h b/src/MC/fix_atom_swap.h
index 65d811e2b9..2511fe6544 100644
--- a/src/MC/fix_atom_swap.h
+++ b/src/MC/fix_atom_swap.h
@@ -40,7 +40,7 @@ class FixAtomSwap : public Fix {
 
  private:
   int nevery, seed;
-  int conserve_ke_flag;    // yes = conserve ke, no = do not conserve ke
+  int ke_flag;              // yes = conserve ke, no = do not conserve ke
   int semi_grand_flag;     // yes = semi-grand canonical, no = constant composition
   int ncycles;
   int niswap, njswap;                  // # of i,j swap atoms on all procs
diff --git a/src/MC/fix_mol_swap.cpp b/src/MC/fix_mol_swap.cpp
index 842a76593a..cd4d9d3cb8 100644
--- a/src/MC/fix_mol_swap.cpp
+++ b/src/MC/fix_mol_swap.cpp
@@ -35,12 +35,14 @@
 using namespace LAMMPS_NS;
 using namespace FixConst;
 
+#define BIG 1.0e20
+
 /* ---------------------------------------------------------------------- */
 
 FixMolSwap::FixMolSwap(LAMMPS *lmp, int narg, char **arg) :
   Fix(lmp, narg, arg)
 {
-  if (narg != 9) error->all(FLERR,"Illegal fix mol/swap command");
+  if (narg < 9) error->all(FLERR,"Illegal fix mol/swap command");
 
   vector_flag = 1;
   size_vector = 2;
@@ -58,13 +60,32 @@ FixMolSwap::FixMolSwap(LAMMPS *lmp, int narg, char **arg) :
   seed = utils::inumeric(FLERR,arg[7],false,lmp);
   double temperature = utils::numeric(FLERR,arg[8],false,lmp);
 
+  // optional args
+
+  ke_flag = 1;
+
+  int iarg = 9;
+  while (iarg < narg) {
+    if (strcmp(arg[iarg],"ke") == 0) {
+      if (iarg+2 > narg) error->all(FLERR,"Illegal fix mol/swap command");
+      ke_flag = utils::logical(FLERR,arg[iarg+1],false,lmp);
+      iarg += 2;
+    } else error->all(FLERR,"Illegal fix mol/swap command");
+  }
+
+  // error check
+
   if (nevery <= 0) error->all(FLERR,"Illegal fix mol/swap command");
   if (ncycles < 0) error->all(FLERR,"Illegal fix mol/swap command");
+  if (itype == jtype) error->all(FLERR,"Illegal fix mol/swap command");
   if (itype <= 0 || itype > atom->ntypes ||
       jtype <= 0 || jtype > atom->ntypes)
     error->all(FLERR,"Fix mol/swap atom types are invalid");
   if (seed <= 0) error->all(FLERR,"Illegal fix mol/swap command");
   if (temperature <= 0.0) error->all(FLERR,"Illegal fix mol/swap command");
+  if (ke_flag && atom->rmass) 
+    error->all(FLERR,"Cannot conserve kinetic energy with fix mol/swap "
+               "unless per-type masses");
 
   beta = 1.0/(force->boltz*temperature);
 
@@ -79,12 +100,17 @@ FixMolSwap::FixMolSwap(LAMMPS *lmp, int narg, char **arg) :
 
   // zero out counters
 
-  nswap_attempts = 0.0;
-  nswap_successes = 0.0;
+  nswap_attempt = 0.0;
+  nswap_accept = 0.0;
+
+  // qflag = 1 if charges are defined
+
+  if (atom->q_flag) qflag = 1;
+  else qflag = 0;
 
   // set comm size needed by this Fix
 
-  if (atom->q_flag) comm_forward = 2;
+  if (qflag) comm_forward = 2;
   else comm_forward = 1;
 }
 
@@ -141,6 +167,54 @@ void FixMolSwap::init()
   MPI_Allreduce(&minmol_me,&minmol,1,MPI_LMP_TAGINT,MPI_MIN,world);
   MPI_Allreduce(&maxmol_me,&maxmol,1,MPI_LMP_TAGINT,MPI_MAX,world);
 
+  // if ke_flag, check if itype/jtype masses are different
+  // if yes, pre-calcuate velocity scale factors that keep KE constant
+  // if no, turn ke_flag off
+
+  if (ke_flag) {
+    double *mass = atom->mass;
+    if (mass[itype] != mass[jtype]) {
+      i2j_vscale = sqrt(mass[itype]/mass[jtype]);
+      j2i_vscale = sqrt(mass[jtype]/mass[itype]);
+    } else ke_flag = 0;
+  }
+
+  // if qflag, check if all charges for itype are identical, ditto for jtype
+  // if not, cannot swap charges, issue warning
+
+  if (qflag) {
+    double *q = atom->q;
+
+    double iqone,jqone;
+    iqone = jqone = -BIG;
+
+    for (int i = 0; i < nlocal; i++) {
+      if (molecule[i] == 0) continue;
+      if (!(mask[i] & groupbit)) continue;
+      if (type[i] == itype) iqone = q[i];
+      if (type[i] == jtype) jqone = q[i];
+    }
+
+    MPI_Allreduce(&iqone,&iq,1,MPI_DOUBLE,MPI_MAX,world);
+    MPI_Allreduce(&jqone,&jq,1,MPI_DOUBLE,MPI_MAX,world);
+
+    int flag = 0;
+
+    for (int i = 0; i < nlocal; i++) {
+      if (molecule[i] == 0) continue;
+      if (!(mask[i] & groupbit)) continue;
+      if (type[i] == itype && q[i] != iq) flag = 1;
+      if (type[i] == jtype && q[i] != jq) flag = 1;
+    }
+
+    int flagall;
+    MPI_Allreduce(&flag,&flagall,1,MPI_INT,MPI_MAX,world);
+    if (flagall) qflag = 0;
+
+    if (!qflag && comm->me == 0) 
+      error->warning(FLERR,"Cannot swap charges in fix mol/swap");
+  }
+
   // check if all cutoffs involving itype and jtype are the same
   // if not, reneighboring will be needed after swaps
 
@@ -152,7 +226,7 @@ void FixMolSwap::init()
 }
 
 /* ----------------------------------------------------------------------
-   perform Ncycle Monte Carlo swaps
+   perform ncycle Monte Carlo swaps
 ------------------------------------------------------------------------- */
 
 void FixMolSwap::pre_exchange()
@@ -178,11 +252,13 @@ void FixMolSwap::pre_exchange()
 
   // attempt Ncycle molecule swaps
 
-  int nsuccess = 0;
-  for (int m = 0; m < ncycles; m++) nsuccess += attempt_swap();
+  int naccept = 0;
+  for (int m = 0; m < ncycles; m++) naccept += attempt_swap();
 
-  nswap_attempts += ncycles;
-  nswap_successes += nsuccess;
+  // udpate MC stats
+
+  nswap_attempt += ncycles;
+  nswap_accept += naccept;
 
   next_reneighbor = update->ntimestep + nevery;
 }
@@ -199,13 +275,15 @@ int FixMolSwap::attempt_swap()
   double energy_before = energy_stored;
 
   // pick a random molecule
-  // swap all of its eligible itype & jtype atoms
+  // swap properties of all its eligible itype & jtype atoms
 
   tagint molID = 
     minmol + static_cast<tagint> (random->uniform() * (maxmol-minmol+1));
   if (molID > maxmol) molID = maxmol;
 
   int *mask = atom->mask;
+  double **v = atom->v;
+  double *q = atom->q;
   int *type = atom->type;
   tagint *molecule = atom->molecule;
   int nlocal = atom->nlocal;
@@ -213,8 +291,23 @@ int FixMolSwap::attempt_swap()
   for (int i = 0; i < nlocal; i++) {
     if (molecule[i] != molID) continue;
     if (!(mask[i] & groupbit)) continue;
-    if (type[i] == itype) type[i] = jtype;
-    else if (type[i] == jtype) type[i] = itype;
+    if (type[i] == itype) {
+      type[i] = jtype;
+      if (qflag) q[i] = jq;
+      if (ke_flag) {
+        v[i][0] *= i2j_vscale;
+        v[i][1] *= i2j_vscale;
+        v[i][2] *= i2j_vscale;
+      }
+    } else if (type[i] == jtype) {
+      type[i] = itype;
+      if (qflag) q[i] = iq;
+      if (ke_flag) {
+        v[i][0] *= j2i_vscale;
+        v[i][1] *= j2i_vscale;
+        v[i][2] *= j2i_vscale;
+      }
+    }
   }
 
   // if unequal_cutoffs, call comm->borders() and rebuild neighbor list
@@ -236,23 +329,37 @@ int FixMolSwap::attempt_swap()
 
   double energy_after = energy_full();
 
-  // swap accepted
+  // swap accepted, return 1
 
   if (random->uniform() < exp(beta*(energy_before - energy_after))) {
     energy_stored = energy_after;
     return 1;
+  }
 
-  // swap not accepted
+  // swap not accepted, return 0
   // restore the swapped itype & jtype atoms
   // do not need to re-call comm->borders() and rebuild neighbor list
-  //   since will be done on next cycle or in Verlet when this fix is done
+  //   since will be done on next cycle or in Verlet when this fix finishes
 
-  } else {
-    for (int i = 0; i < nlocal; i++) {
-      if (molecule[i] != molID) continue;
-      if (!(mask[i] & groupbit)) continue;
-      if (type[i] == itype) type[i] = jtype;
-      else if (type[i] == jtype) type[i] = itype;
+  for (int i = 0; i < nlocal; i++) {
+    if (molecule[i] != molID) continue;
+    if (!(mask[i] & groupbit)) continue;
+    if (type[i] == itype) {
+      type[i] = jtype;
+      if (qflag) q[i] = jq;
+      if (ke_flag) {
+        v[i][0] *= i2j_vscale;
+        v[i][1] *= i2j_vscale;
+        v[i][2] *= i2j_vscale;
+      }
+    } else if (type[i] == jtype) {
+      type[i] = itype;
+      if (qflag) q[i] = iq;
+      if (ke_flag) {
+        v[i][0] *= j2i_vscale;
+        v[i][1] *= j2i_vscale;
+        v[i][2] *= j2i_vscale;
+      }
     }
   }
 
@@ -300,7 +407,7 @@ int FixMolSwap::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*
 
   m = 0;
 
-  if (atom->q_flag) {
+  if (qflag) {
     for (i = 0; i < n; i++) {
       j = list[i];
       buf[m++] = type[j];
@@ -328,7 +435,7 @@ void FixMolSwap::unpack_forward_comm(int n, int first, double *buf)
   m = 0;
   last = first + n;
 
-  if (atom->q_flag) {
+  if (qflag) {
     for (i = first; i < last; i++) {
       type[i] = static_cast<int> (buf[m++]);
       q[i] = buf[m++];
@@ -345,8 +452,8 @@ void FixMolSwap::unpack_forward_comm(int n, int first, double *buf)
 
 double FixMolSwap::compute_vector(int n)
 {
-  if (n == 0) return nswap_attempts;
-  if (n == 1) return nswap_successes;
+  if (n == 0) return nswap_attempt;
+  if (n == 1) return nswap_accept;
   return 0.0;
 }
 
@@ -360,8 +467,8 @@ void FixMolSwap::write_restart(FILE *fp)
   double list[6];
   list[n++] = random->state();
   list[n++] = ubuf(next_reneighbor).d;
-  list[n++] = nswap_attempts;
-  list[n++] = nswap_successes;
+  list[n++] = nswap_attempt;
+  list[n++] = nswap_accept;
   list[n++] = ubuf(update->ntimestep).d;
 
   if (comm->me == 0) {
@@ -385,8 +492,8 @@ void FixMolSwap::restart(char *buf)
 
   next_reneighbor = (bigint) ubuf(list[n++]).i;
 
-  nswap_attempts = static_cast<int>(list[n++]);
-  nswap_successes = static_cast<int>(list[n++]);
+  nswap_attempt = static_cast<int>(list[n++]);
+  nswap_accept = static_cast<int>(list[n++]);
 
   bigint ntimestep_restart = (bigint) ubuf(list[n++]).i;
   if (ntimestep_restart != update->ntimestep)
diff --git a/src/MC/fix_mol_swap.h b/src/MC/fix_mol_swap.h
index c84ed815c4..6c91d27b91 100644
--- a/src/MC/fix_mol_swap.h
+++ b/src/MC/fix_mol_swap.h
@@ -42,12 +42,21 @@ class FixMolSwap : public Fix {
   int itype, jtype;
   double temperature;
 
-  bool unequal_cutoffs;
-  tagint minmol,maxmol;
-  double nswap_attempts;
-  double nswap_successes;
-  double beta;
-  double energy_stored;
+  int ke_flag;            // 1 if kinetic energy is also swapped
+  double i2j_vscale;      // scale factors for velocity to keep KE constant
+  double j2i_vscale; 
+
+  int qflag;              // 1 if charge is also swapped
+  double iq,jq;           // charge values for all itype,jtype atoms
+
+  bool unequal_cutoffs;   // 1 if itype and jtype have any different cutoffs
+  tagint minmol,maxmol;   // range of mol IDs selected for swaps
+
+  double nswap_attempt;   // cummulative stats on MC attempts and accepts
+  double nswap_accept;
+
+  double beta;            // 1/kT
+  double energy_stored;   // energy of current state as swaps are accepted
 
   class RanPark *random;
   class Compute *c_pe;