diff --git a/src/USER-MISC/pair_lj_sf.cpp b/src/USER-MISC/pair_lj_sf.cpp
index d378cfab4..2e1c92c93 100644
--- a/src/USER-MISC/pair_lj_sf.cpp
+++ b/src/USER-MISC/pair_lj_sf.cpp
@@ -1,351 +1,356 @@
 /* ----------------------------------------------------------------------
    LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
    http://lammps.sandia.gov, Sandia National Laboratories
    Steve Plimpton, sjplimp@sandia.gov
 
    Copyright (2003) Sandia Corporation.  Under the terms of Contract
    DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
    certain rights in this software.  This software is distributed under
    the GNU General Public License.
 
    See the README file in the top-level LAMMPS directory.
 ------------------------------------------------------------------------- */
 
 /* ----------------------------------------------------------------------
    Contributing author: Laurent Joly (U Lyon), ljoly.ulyon@gmail.com
 ------------------------------------------------------------------------- */
 
 #include "math.h"
 #include "stdio.h"
 #include "stdlib.h"
 #include "pair_lj_sf.h"
 #include "atom.h"
 #include "comm.h"
 #include "force.h"
 #include "neigh_list.h"
 #include "memory.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
 
 /* ---------------------------------------------------------------------- */
 
-PairLJShiftedForce::PairLJShiftedForce(LAMMPS *lmp) : Pair(lmp) {}
+PairLJShiftedForce::PairLJShiftedForce(LAMMPS *lmp) : Pair(lmp)
+{
+  respa_enable = 0;
+}
 
 /* ---------------------------------------------------------------------- */
 
 PairLJShiftedForce::~PairLJShiftedForce()
 {
   if (allocated) {
     memory->destroy(setflag);
     memory->destroy(cutsq);
 
     memory->destroy(cut);
     memory->destroy(epsilon);
     memory->destroy(sigma);
     memory->destroy(lj1);
     memory->destroy(lj2);
     memory->destroy(lj3);
     memory->destroy(lj4);
     memory->destroy(foffset);
     memory->destroy(offset);
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void PairLJShiftedForce::compute(int eflag, int vflag)
 {
   int i,j,ii,jj,inum,jnum,itype,jtype;
   double xtmp,ytmp,ztmp,delx,dely,delz,evdwl,fpair;
   double rsq,r2inv,r6inv,forcelj,factor_lj;
   double r,t;
   int *ilist,*jlist,*numneigh,**firstneigh;
 
   evdwl = 0.0;
   if (eflag || vflag) ev_setup(eflag,vflag);
   else evflag = vflag_fdotr = 0;
 
   double **x = atom->x;
   double **f = atom->f;
   int *type = atom->type;
   int nlocal = atom->nlocal;
   int nall = nlocal + atom->nghost;
   double *special_lj = force->special_lj;
   int newton_pair = force->newton_pair;
 
   inum = list->inum;
   ilist = list->ilist;
   numneigh = list->numneigh;
   firstneigh = list->firstneigh;
 
   // loop over neighbors of my atoms
 
   for (ii = 0; ii < inum; ii++) {
     i = ilist[ii];
     xtmp = x[i][0];
     ytmp = x[i][1];
     ztmp = x[i][2];
     itype = type[i];
     jlist = firstneigh[i];
     jnum = numneigh[i];
 
     for (jj = 0; jj < jnum; jj++) {
       j = jlist[jj];
       factor_lj = special_lj[sbmask(j)];
       j &= NEIGHMASK;
 
       delx = xtmp - x[j][0];
       dely = ytmp - x[j][1];
       delz = ztmp - x[j][2];
       rsq = delx*delx + dely*dely + delz*delz;
       jtype = type[j];
 
       if (rsq < cutsq[itype][jtype]) {
         r2inv = 1.0/rsq;
         r6inv = r2inv*r2inv*r2inv;
         r = sqrt(rsq);
         t = r/cut[itype][jtype];
         forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]) -
           t*foffset[itype][jtype];
         fpair = factor_lj*forcelj*r2inv;
 
         f[i][0] += delx*fpair;
         f[i][1] += dely*fpair;
         f[i][2] += delz*fpair;
         if (newton_pair || j < nlocal) {
           f[j][0] -= delx*fpair;
           f[j][1] -= dely*fpair;
           f[j][2] -= delz*fpair;
         }
 
         if (eflag) {
           evdwl = r6inv*(lj3[itype][jtype]*r6inv-lj4[itype][jtype]) +
             (t-1.0)*foffset[itype][jtype] - offset[itype][jtype];
           evdwl *= factor_lj;
         }
 
         if (evflag) ev_tally(i,j,nlocal,newton_pair,
                              evdwl,0.0,fpair,delx,dely,delz);
       }
     }
   }
 
   if (vflag_fdotr) virial_fdotr_compute();
 }
 
 /* ----------------------------------------------------------------------
    allocate all arrays
 ------------------------------------------------------------------------- */
 
 void PairLJShiftedForce::allocate()
 {
   allocated = 1;
   int n = atom->ntypes;
 
   memory->create(setflag,n+1,n+1,"pair:setflag");
   for (int i = 1; i <= n; i++)
     for (int j = i; j <= n; j++)
       setflag[i][j] = 0;
 
   memory->create(cutsq,n+1,n+1,"pair:cutsq");
+
   memory->create(cut,n+1,n+1,"pair:cut");
   memory->create(epsilon,n+1,n+1,"pair:epsilon");
   memory->create(sigma,n+1,n+1,"pair:sigma");
   memory->create(lj1,n+1,n+1,"pair:lj1");
-  memory->create(lj3,n+1,n+1,"pair:lj2");
+  memory->create(lj2,n+1,n+1,"pair:lj2");
   memory->create(lj3,n+1,n+1,"pair:lj3");
   memory->create(lj4,n+1,n+1,"pair:lj4");
   memory->create(foffset,n+1,n+1,"pair:foffset");
   memory->create(offset,n+1,n+1,"pair:offset");
 }
 
 /* ----------------------------------------------------------------------
    global settings
 ------------------------------------------------------------------------- */
 
 void PairLJShiftedForce::settings(int narg, char **arg)
 {
   if (narg != 1) error->all(FLERR,"Illegal pair_style command");
 
   cut_global = force->numeric(arg[0]);
 
   if (cut_global <= 0.0)
     error->all(FLERR,"Illegal pair_style command");
 
   // reset cutoffs that have been explicitly set
 
   if (allocated) {
     int i,j;
     for (i = 1; i <= atom->ntypes; i++)
       for (j = i+1; j <= atom->ntypes; j++)
         if (setflag[i][j]) cut[i][j] = cut_global;
   }
 }
 
 /* ----------------------------------------------------------------------
    set coeffs for one or more type pairs
 ------------------------------------------------------------------------- */
 
 void PairLJShiftedForce::coeff(int narg, char **arg)
 {
-  if (narg < 4 || narg > 5) error->all(FLERR,"Incorrect args for pair coefficients");
+  if (narg < 4 || narg > 5)
+    error->all(FLERR,"Incorrect args for pair coefficients");
   if (!allocated) allocate();
 
   int ilo,ihi,jlo,jhi;
   force->bounds(arg[0],atom->ntypes,ilo,ihi);
   force->bounds(arg[1],atom->ntypes,jlo,jhi);
 
   double epsilon_one = force->numeric(arg[2]);
   double sigma_one = force->numeric(arg[3]);
 
   double cut_one = cut_global;
   if (narg == 5) cut_one = force->numeric(arg[4]);
 
   if (cut_one <= 0.0)
     error->all(FLERR,"Incorrect args for pair coefficients");
 
   int count = 0;
   for (int i = ilo; i <= ihi; i++) {
     for (int j = MAX(jlo,i); j <= jhi; j++) {
       epsilon[i][j] = epsilon_one;
       sigma[i][j] = sigma_one;
       cut[i][j] = cut_one;
       setflag[i][j] = 1;
       count++;
     }
   }
 
   if (count == 0) error->all(FLERR,"Incorrect args for pair coefficients");
 }
 
 /* ----------------------------------------------------------------------
    init for one type pair i,j and corresponding j,i
 ------------------------------------------------------------------------- */
 
 double PairLJShiftedForce::init_one(int i, int j)
 {
   if (setflag[i][j] == 0) {
     epsilon[i][j] = mix_energy(epsilon[i][i],epsilon[j][j],
                                sigma[i][i],sigma[j][j]);
     sigma[i][j] = mix_distance(sigma[i][i],sigma[j][j]);
     cut[i][j] = mix_distance(cut[i][i],cut[j][j]);
   }
 
   lj1[i][j] = 48.0 * epsilon[i][j] * pow(sigma[i][j],12.0);
   lj2[i][j] = 24.0 * epsilon[i][j] * pow(sigma[i][j],6.0);
   lj3[i][j] = 4.0 * epsilon[i][j] * pow(sigma[i][j],12.0);
   lj4[i][j] = 4.0 * epsilon[i][j] * pow(sigma[i][j],6.0);
 
   double ratio = sigma[i][j] / cut[i][j];
   foffset[i][j] = 4.0 * epsilon[i][j] * (12.0 * pow(ratio,12.0) -
                                          6.0 * pow(ratio,6.0));
   offset[i][j] = 4.0 * epsilon[i][j] * (pow(ratio,12.0) - pow(ratio,6.0));
 
   cut[j][i] = cut[i][j];
   lj1[j][i] = lj1[i][j];
   lj2[j][i] = lj2[i][j];
   lj3[j][i] = lj3[i][j];
   lj4[j][i] = lj4[i][j];
   foffset[j][i] = foffset[i][j];
   offset[j][i] = offset[i][j];
 
   return cut[i][j];
 }
 
 /* ----------------------------------------------------------------------
    proc 0 writes to restart file
 ------------------------------------------------------------------------- */
 
 void PairLJShiftedForce::write_restart(FILE *fp)
 {
   write_restart_settings(fp);
 
   int i,j;
   for (i = 1; i <= atom->ntypes; i++)
     for (j = i; j <= atom->ntypes; j++) {
       fwrite(&setflag[i][j],sizeof(int),1,fp);
       if (setflag[i][j]) {
         fwrite(&epsilon[i][j],sizeof(double),1,fp);
         fwrite(&sigma[i][j],sizeof(double),1,fp);
         fwrite(&cut[i][j],sizeof(double),1,fp);
       }
     }
 }
 
 /* ----------------------------------------------------------------------
    proc 0 reads from restart file, bcasts
 ------------------------------------------------------------------------- */
 
 void PairLJShiftedForce::read_restart(FILE *fp)
 {
   read_restart_settings(fp);
   allocate();
 
   int i,j;
   int me = comm->me;
   for (i = 1; i <= atom->ntypes; i++)
     for (j = i; j <= atom->ntypes; j++) {
       if (me == 0) fread(&setflag[i][j],sizeof(int),1,fp);
       MPI_Bcast(&setflag[i][j],1,MPI_INT,0,world);
       if (setflag[i][j]) {
         if (me == 0) {
           fread(&epsilon[i][j],sizeof(double),1,fp);
           fread(&sigma[i][j],sizeof(double),1,fp);
           fread(&cut[i][j],sizeof(double),1,fp);
         }
         MPI_Bcast(&epsilon[i][j],1,MPI_DOUBLE,0,world);
         MPI_Bcast(&sigma[i][j],1,MPI_DOUBLE,0,world);
         MPI_Bcast(&cut[i][j],1,MPI_DOUBLE,0,world);
       }
     }
 }
 
 /* ----------------------------------------------------------------------
    proc 0 writes to restart file
 ------------------------------------------------------------------------- */
 
 void PairLJShiftedForce::write_restart_settings(FILE *fp)
 {
   fwrite(&cut_global,sizeof(double),1,fp);
   fwrite(&offset_flag,sizeof(int),1,fp);
   fwrite(&mix_flag,sizeof(int),1,fp);
 }
 
 /* ----------------------------------------------------------------------
    proc 0 reads from restart file, bcasts
 ------------------------------------------------------------------------- */
 
 void PairLJShiftedForce::read_restart_settings(FILE *fp)
 {
   int me = comm->me;
   if (me == 0) {
     fread(&cut_global,sizeof(double),1,fp);
     fread(&offset_flag,sizeof(int),1,fp);
     fread(&mix_flag,sizeof(int),1,fp);
   }
   MPI_Bcast(&cut_global,1,MPI_DOUBLE,0,world);
   MPI_Bcast(&offset_flag,1,MPI_INT,0,world);
   MPI_Bcast(&mix_flag,1,MPI_INT,0,world);
 }
 
 /* ---------------------------------------------------------------------- */
 
 double PairLJShiftedForce::single(int i, int j, int itype, int jtype, double rsq,
                          double factor_coul, double factor_lj,
                          double &fforce)
 {
   double r2inv,r6inv,forcelj,philj,r,t;
 
   r2inv = 1.0/rsq;
   r6inv = r2inv*r2inv*r2inv;
   r = sqrt(rsq);
   t = r/cut[itype][jtype];
   forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]) -
     t*foffset[itype][jtype];
   fforce = factor_lj*forcelj*r2inv;
 
   philj = r6inv*(lj3[itype][jtype]*r6inv-lj4[itype][jtype]) +
     (t-1.0)*foffset[itype][jtype] - offset[itype][jtype];
   return factor_lj*philj;
 }
diff --git a/src/fix_adapt.cpp b/src/fix_adapt.cpp
index 2f0c735ef..49684020c 100644
--- a/src/fix_adapt.cpp
+++ b/src/fix_adapt.cpp
@@ -1,393 +1,389 @@
 /* ----------------------------------------------------------------------
    LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
    http://lammps.sandia.gov, Sandia National Laboratories
    Steve Plimpton, sjplimp@sandia.gov
 
    Copyright (2003) Sandia Corporation.  Under the terms of Contract
    DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
    certain rights in this software.  This software is distributed under
    the GNU General Public License.
 
    See the README file in the top-level LAMMPS directory.
 ------------------------------------------------------------------------- */
 
 #include "math.h"
 #include "string.h"
 #include "stdlib.h"
 #include "fix_adapt.h"
 #include "atom.h"
 #include "update.h"
 #include "modify.h"
 #include "force.h"
 #include "pair.h"
 #include "pair_hybrid.h"
 #include "kspace.h"
 #include "input.h"
 #include "variable.h"
 #include "math_const.h"
 #include "memory.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
 using namespace FixConst;
 using namespace MathConst;
 
 enum{PAIR,KSPACE,ATOM};
 enum{DIAMETER,CHARGE};
 
 /* ---------------------------------------------------------------------- */
 
 FixAdapt::FixAdapt(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
 {
   if (narg < 5) error->all(FLERR,"Illegal fix adapt command");
   nevery = atoi(arg[3]);
   if (nevery < 0) error->all(FLERR,"Illegal fix adapt command");
 
   // count # of adaptations
 
   nadapt = 0;
 
   int iarg = 4;
   while (iarg < narg) {
     if (strcmp(arg[iarg],"pair") == 0) {
       if (iarg+6 > narg) error->all(FLERR,"Illegal fix adapt command");
       nadapt++;
       iarg += 6;
     } else if (strcmp(arg[iarg],"kspace") == 0) {
       if (iarg+2 > narg) error->all(FLERR,"Illegal fix adapt command");
       nadapt++;
       iarg += 2;
     } else if (strcmp(arg[iarg],"atom") == 0) {
       if (iarg+3 > narg) error->all(FLERR,"Illegal fix adapt command");
       nadapt++;
       iarg += 3;
     } else break;
   }
 
   if (nadapt == 0) error->all(FLERR,"Illegal fix adapt command");
   adapt = new Adapt[nadapt];
 
   // parse keywords
 
   nadapt = 0;
   diamflag = 0;
   chgflag = 0;
 
   iarg = 4;
   while (iarg < narg) {
     if (strcmp(arg[iarg],"pair") == 0) {
       if (iarg+6 > narg) error->all(FLERR,"Illegal fix adapt command");
       adapt[nadapt].which = PAIR;
       int n = strlen(arg[iarg+1]) + 1;
       adapt[nadapt].pstyle = new char[n];
       strcpy(adapt[nadapt].pstyle,arg[iarg+1]);
       n = strlen(arg[iarg+2]) + 1;
       adapt[nadapt].pparam = new char[n];
       strcpy(adapt[nadapt].pparam,arg[iarg+2]);
       force->bounds(arg[iarg+3],atom->ntypes,
                     adapt[nadapt].ilo,adapt[nadapt].ihi);
       force->bounds(arg[iarg+4],atom->ntypes,
                     adapt[nadapt].jlo,adapt[nadapt].jhi);
       if (strstr(arg[iarg+5],"v_") == arg[iarg+5]) {
         n = strlen(&arg[iarg+5][2]) + 1;
         adapt[nadapt].var = new char[n];
         strcpy(adapt[nadapt].var,&arg[iarg+5][2]);
       } else error->all(FLERR,"Illegal fix adapt command");
       nadapt++;
       iarg += 6;
     } else if (strcmp(arg[iarg],"kspace") == 0) {
       if (iarg+2 > narg) error->all(FLERR,"Illegal fix adapt command");
       adapt[nadapt].which = KSPACE;
       if (strstr(arg[iarg+1],"v_") == arg[iarg+1]) {
         int n = strlen(&arg[iarg+1][2]) + 1;
         adapt[nadapt].var = new char[n];
         strcpy(adapt[nadapt].var,&arg[iarg+1][2]);
       } else error->all(FLERR,"Illegal fix adapt command");
       nadapt++;
       iarg += 2;
     } else if (strcmp(arg[iarg],"atom") == 0) {
       if (iarg+3 > narg) error->all(FLERR,"Illegal fix adapt command");
       adapt[nadapt].which = ATOM;
       if (strcmp(arg[iarg+1],"diameter") == 0) {
         adapt[nadapt].aparam = DIAMETER;
         diamflag = 1;
       } else if (strcmp(arg[iarg+1],"charge") == 0) {
         adapt[nadapt].aparam = CHARGE; 
         chgflag = 1; 
       } else error->all(FLERR,"Illegal fix adapt command");
       if (strstr(arg[iarg+2],"v_") == arg[iarg+2]) {
         int n = strlen(&arg[iarg+2][2]) + 1;
         adapt[nadapt].var = new char[n];
         strcpy(adapt[nadapt].var,&arg[iarg+2][2]);
       } else error->all(FLERR,"Illegal fix adapt command");
       nadapt++;
       iarg += 3;
     } else break;
   }
 
   // optional keywords
 
   resetflag = 0;
   scaleflag = 0;
 
   while (iarg < narg) {
     if (strcmp(arg[iarg],"reset") == 0) {
       if (iarg+2 > narg) error->all(FLERR,"Illegal fix adapt command");
       if (strcmp(arg[iarg+1],"no") == 0) resetflag = 0;
       else if (strcmp(arg[iarg+1],"yes") == 0) resetflag = 1;
       else error->all(FLERR,"Illegal fix adapt command");
       iarg += 2;
     } else if (strcmp(arg[iarg],"scale") == 0) {
       if (iarg+2 > narg) error->all(FLERR,"Illegal fix adapt command");
       if (strcmp(arg[iarg+1],"no") == 0) scaleflag = 0;
       else if (strcmp(arg[iarg+1],"yes") == 0) scaleflag = 1;
       else error->all(FLERR,"Illegal fix adapt command");
       iarg += 2;
     } else error->all(FLERR,"Illegal fix adapt command");
   }
 
   // allocate pair style arrays
 
   int n = atom->ntypes;
   for (int m = 0; m < nadapt; m++)
     if (adapt[m].which == PAIR)
       memory->create(adapt[m].array_orig,n+1,n+1,"adapt:array_orig");
 }
 
 /* ---------------------------------------------------------------------- */
 
 FixAdapt::~FixAdapt()
 {
   for (int m = 0; m < nadapt; m++) {
     delete [] adapt[m].var;
     if (adapt[m].which == PAIR) {
       delete [] adapt[m].pstyle;
       delete [] adapt[m].pparam;
       memory->destroy(adapt[m].array_orig);
     }
   }
   delete [] adapt;
 }
 
 /* ---------------------------------------------------------------------- */
 
 int FixAdapt::setmask()
 {
   int mask = 0;
   mask |= PRE_FORCE;
   mask |= POST_RUN;
   return mask;
 }
 
 /* ---------------------------------------------------------------------- */
 
 void FixAdapt::init()
 {
   int i,j;
 
   // setup and error checks
 
   anypair = 0;
 
   for (int m = 0; m < nadapt; m++) {
     Adapt *ad = &adapt[m];
 
     ad->ivar = input->variable->find(ad->var);
     if (ad->ivar < 0)
       error->all(FLERR,"Variable name for fix adapt does not exist");
     if (!input->variable->equalstyle(ad->ivar))
       error->all(FLERR,"Variable for fix adapt is invalid style");
 
     if (ad->which == PAIR) {
       anypair = 1;
 
       Pair *pair = force->pair_match(ad->pstyle,1);
       if (pair == NULL) error->all(FLERR,"Fix adapt pair style does not exist");
       void *ptr = pair->extract(ad->pparam,ad->pdim);
       if (ptr == NULL) 
         error->all(FLERR,"Fix adapt pair style param not supported");
 
       ad->pdim = 2;
       if (ad->pdim == 0) ad->scalar = (double *) ptr;
       if (ad->pdim == 2) ad->array = (double **) ptr;
 
       // if pair hybrid, test that ilo,ihi,jlo,jhi are valid for sub-style
 
       if (ad->pdim == 2 && (strcmp(force->pair_style,"hybrid") == 0 ||
                             strcmp(force->pair_style,"hybrid/overlay") == 0)) {
         PairHybrid *pair = (PairHybrid *) force->pair;
         for (i = ad->ilo; i <= ad->ihi; i++)
           for (j = MAX(ad->jlo,i); j <= ad->jhi; j++)
             if (!pair->check_ijtype(i,j,ad->pstyle))
               error->all(FLERR,"Fix adapt type pair range is not valid for "
                          "pair hybrid sub-style");
       }
 
     } else if (ad->which == KSPACE) {
       if (force->kspace == NULL)
         error->all(FLERR,"Fix adapt kspace style does not exist");
       kspace_scale = (double *) force->kspace->extract("scale");
 
     } else if (ad->which == ATOM) {
       if (ad->aparam == DIAMETER) {
         if (!atom->radius_flag)
           error->all(FLERR,"Fix adapt requires atom attribute diameter");
       }
       if (ad->aparam == CHARGE) {
 	if (!atom->q_flag)
 	  error->all(FLERR,"Fix adapt requires atom attribute charge");
       }
     }
   }
 
   // make copy of original pair array values
 
   for (int m = 0; m < nadapt; m++) {
     Adapt *ad = &adapt[m];
     if (ad->which == PAIR && ad->pdim == 2) {
       for (i = ad->ilo; i <= ad->ihi; i++)
         for (j = MAX(ad->jlo,i); j <= ad->jhi; j++)
           ad->array_orig[i][j] = ad->array[i][j];
     }
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void FixAdapt::setup_pre_force(int vflag)
 {
   change_settings();
 }
 
 /* ---------------------------------------------------------------------- */
 
 void FixAdapt::pre_force(int vflag)
 {
   if (nevery == 0) return;
   if (update->ntimestep % nevery) return;
   change_settings();
 }
 
 /* ---------------------------------------------------------------------- */
 
 void FixAdapt::post_run()
 {
   if (resetflag) restore_settings();
 }
 
 /* ----------------------------------------------------------------------
    change pair,kspace,atom parameters based on variable evaluation
 ------------------------------------------------------------------------- */
 
 void FixAdapt::change_settings()
 {
   int i,j;
 
   // variable evaluation may invoke computes so wrap with clear/add
 
   modify->clearstep_compute();
 
   for (int m = 0; m < nadapt; m++) {
     Adapt *ad = &adapt[m];
     double value = input->variable->compute_equal(ad->ivar);
 
     // set global scalar or type pair array values
 
     if (ad->which == PAIR) {
       if (ad->pdim == 0) {
         if (scaleflag) *ad->scalar = value * ad->scalar_orig;
         else *ad->scalar = value;
       } else if (ad->pdim == 2) {
         if (scaleflag)
           for (i = ad->ilo; i <= ad->ihi; i++)
             for (j = MAX(ad->jlo,i); j <= ad->jhi; j++)
               ad->array[i][j] = value*ad->array_orig[i][j];
         else
           for (i = ad->ilo; i <= ad->ihi; i++)
             for (j = MAX(ad->jlo,i); j <= ad->jhi; j++)
               ad->array[i][j] = value;
       }
 
     // set kspace scale factor
 
     } else if (ad->which == KSPACE) {
       *kspace_scale = value;
 
     } else if (ad->which == ATOM) {
 
       // set radius from diameter
       // also scale rmass to new value
 
       if (ad->aparam == DIAMETER) {
         int mflag = 0;
         if (atom->rmass_flag) mflag = 1;
         double density;
 
         double *radius = atom->radius;
         double *rmass = atom->rmass;
         int *mask = atom->mask;
         int nlocal = atom->nlocal;
 
         if (mflag == 0) {
           for (i = 0; i < nlocal; i++)
             if (mask[i] & groupbit)
               radius[i] = 0.5*value;
         } else {
           for (i = 0; i < nlocal; i++)
             if (mask[i] & groupbit) {
               density = rmass[i] / (4.0*MY_PI/3.0 *
                                     radius[i]*radius[i]*radius[i]);
               radius[i] = 0.5*value;
               rmass[i] = 4.0*MY_PI/3.0 *
                 radius[i]*radius[i]*radius[i] * density;
             }
         }
       } else if (ad->aparam == CHARGE) {
         double *q = atom->q; 
 	int *mask = atom->mask;
 	int nlocal = atom->nlocal;
 	for (i = 0; i < nlocal; i++)
-	  if (mask[i] & groupbit) {
-//              cout<<"Fix adapt called now!\n";
-//            cout<<" Charge of atom: "<<i<<" is: "<<value; 
-            q[i] = value; 
-	  }
+	  if (mask[i] & groupbit) q[i] = value; 
       }
     }
   }
 
   modify->addstep_compute(update->ntimestep + nevery);
 
   // re-initialize pair styles if any PAIR settings were changed
   // this resets other coeffs that may depend on changed values,
   // and also offset and tail corrections
 
   if (anypair) force->pair->reinit();
 }
 
 /* ----------------------------------------------------------------------
    restore pair,kspace.atom parameters to original values
 ------------------------------------------------------------------------- */
 
 void FixAdapt::restore_settings()
 {
   for (int m = 0; m < nadapt; m++) {
     Adapt *ad = &adapt[m];
     if (ad->which == PAIR) {
       if (ad->pdim == 0) *ad->scalar = ad->scalar_orig;
       else if (ad->pdim == 2) {
         for (int i = ad->ilo; i <= ad->ihi; i++)
           for (int j = MAX(ad->jlo,i); j <= ad->jhi; j++)
             ad->array[i][j] = ad->array_orig[i][j];
       }
 
     } else if (ad->which == KSPACE) {
       *kspace_scale = 1.0;
 
     } else if (ad->which == ATOM) {
 
     }
   }
 
   if (anypair) force->pair->reinit();
 }