diff --git a/src/COLLOID/fix_wall_colloid.cpp b/src/COLLOID/fix_wall_colloid.cpp
index 0c03ad57a..7627fca70 100644
--- a/src/COLLOID/fix_wall_colloid.cpp
+++ b/src/COLLOID/fix_wall_colloid.cpp
@@ -1,160 +1,158 @@
 /* ----------------------------------------------------------------------
    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 authors: Jeremy Lechman (SNL)
 ------------------------------------------------------------------------- */
 
 #include "math.h"
 #include "string.h"
 #include "fix_wall_colloid.h"
 #include "atom.h"
 #include "atom_vec.h"
 #include "update.h"
 #include "respa.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
 
 /* ---------------------------------------------------------------------- */
 
 FixWallColloid::FixWallColloid(LAMMPS *lmp, int narg, char **arg) : 
   FixWall(lmp, narg, arg) {}
 
 /* ---------------------------------------------------------------------- */
 
 void FixWallColloid::init()
 {
   if (!atom->sphere_flag) 
     error->all("Fix wall/colloid requires atom style sphere");
 
   // insure all particles in group are extended particles
 
   double *radius = atom->radius;
   int *mask = atom->mask;
   int nlocal = atom->nlocal;
 
   int flag = 0;
   for (int i = 0; i < nlocal; i++)
     if (mask[i] & groupbit)
       if (radius[i] == 0.0) flag = 1;
 
   int flagall;
   MPI_Allreduce(&flag,&flagall,1,MPI_INT,MPI_SUM,world);
   if (flagall) error->all("Fix wall/colloid requires extended particles");
 
   FixWall::init();
 }
 
 /* ---------------------------------------------------------------------- */
 
 void FixWallColloid::precompute(int m)
 {
   coeff1[m] = 4.0/315.0 * epsilon[m] * pow(sigma[m],6.0);
   coeff2[m] = 2.0/3.0 * epsilon[m];
   coeff3[m] = epsilon[m] * pow(sigma[m],6.0)/7560.0;
   coeff4[m] = epsilon[m]/6.0;
 }
 
 /* ----------------------------------------------------------------------
    interaction of all particles in group with a wall
    m = index of wall coeffs
    which = xlo,xhi,ylo,yhi,zlo,zhi
    error if any finite-size particle is touching or penetrating wall
 ------------------------------------------------------------------------- */
 
 void FixWallColloid::wall_particle(int m, int which, double coord)
 {
   double delta,delta2,rinv,r2inv,r4inv,r8inv,fwall;
-  double r2,rinv2,r2inv2,r4inv2,r6inv2;
+  double r2,rinv2,r2inv2,r4inv2;
   double r3,rinv3,r2inv3,r4inv3,r6inv3;
   double rad,rad2,rad4,rad8,diam,new_coeff2;
   double eoffset;
 
   double **x = atom->x;
   double **f = atom->f;
   double *radius = atom->radius;
   int *mask = atom->mask;
   int nlocal = atom->nlocal;
 
   int dim = which  / 2;
   int side = which % 2;
   if (side == 0) side = -1;
 
   int onflag = 0;
 
   for (int i = 0; i < nlocal; i++)
     if (mask[i] & groupbit) {
       if (side < 0) delta = x[i][dim] - coord;
       else delta = coord - x[i][dim];
       if (delta >= cutoff[m]) continue;
       rad = radius[i];
       if (rad >= delta) {
 	onflag = 1;
 	continue;
       }
 
       new_coeff2 = coeff2[m]*rad*rad*rad;
       diam = 2.0*rad;
       rad2 = rad*rad;
       rad4 = rad2*rad2;
       rad8 = rad4*rad4;
       delta2 = rad2 - delta*delta;
       rinv = 1.0/delta2;
       r2inv = rinv*rinv;
       r4inv = r2inv*r2inv;
       r8inv = r4inv*r4inv;
       fwall = side * (coeff1[m]*(rad8*rad + 27.0*rad4*rad2*rad*pow(delta,2.0)
 				 + 63.0*rad4*rad*pow(delta,4.0)
 				 + 21.0*rad2*rad*pow(delta,6.0))*r8inv - 
 		      new_coeff2*r2inv);
       f[i][dim] -= fwall;
 
       r2 = rad - delta;
       rinv2 = 1.0/r2;
       r2inv2 = rinv2*rinv2;
       r4inv2 = r2inv2*r2inv2;
-      r6inv2 = r4inv2*r2inv2;
       r3 = delta + rad;
       rinv3 = 1.0/r3;
       r2inv3 = rinv3*rinv3;
       r4inv3 = r2inv3*r2inv3;
       r6inv3 = r4inv3*r2inv3;
       ewall[0] += coeff3[m]*((-3.5*diam+delta)*r4inv2*r2inv2*rinv2
 			     + (3.5*diam+delta)*r4inv3*r2inv3*rinv3) -
 	coeff4[m]*((diam*delta-r2*r3*(log(-r2)-log(r3)))*
 		   (-rinv2)*rinv3);
 
       // offset depends on particle size
 
       r2 = rad - cutoff[m];
       rinv2 = 1.0/r2;
       r2inv2 = rinv2*rinv2;
       r4inv2 = r2inv2*r2inv2;
-      r6inv2 = r4inv2*r2inv2;
       r3 = cutoff[m] + rad;
       rinv3 = 1.0/r3;
       r2inv3 = rinv3*rinv3;
       r4inv3 = r2inv3*r2inv3;
       r6inv3 = r4inv3*r2inv3;
       eoffset = coeff3[m]*((-3.5*diam+cutoff[m])*r4inv2*r2inv2*rinv2
 			   + (3.5*diam+cutoff[m])*r4inv3*r2inv3*rinv3) -
 	coeff4[m]*((diam*cutoff[m]-r2*r3*(log(-r2)-log(r3)))*
 		   (-rinv2)*rinv3);
       ewall[0] -= eoffset;
 
       ewall[m+1] += fwall;
     }
 
   if (onflag) error->one("Particle on or inside fix wall surface");
 }
diff --git a/src/compute_property_local.cpp b/src/compute_property_local.cpp
index 5aa765246..a0e615bc4 100644
--- a/src/compute_property_local.cpp
+++ b/src/compute_property_local.cpp
@@ -1,830 +1,827 @@
 /* ----------------------------------------------------------------------
    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 "string.h"
 #include "compute_property_local.h"
 #include "atom.h"
 #include "atom_vec.h"
 #include "update.h"
 #include "force.h"
 #include "pair.h"
 #include "neighbor.h"
 #include "neigh_request.h"
 #include "neigh_list.h"
 #include "memory.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
 
 enum{NONE,NEIGH,PAIR,BOND,ANGLE,DIHEDRAL,IMPROPER};
 
 #define DELTA 10000
 
 /* ---------------------------------------------------------------------- */
 
 ComputePropertyLocal::ComputePropertyLocal(LAMMPS *lmp, int narg, char **arg) :
   Compute(lmp, narg, arg)
 {
   if (narg < 4) error->all("Illegal compute property/local command");
 
   local_flag = 1;
   nvalues = narg - 3;
   if (nvalues == 1) size_local_cols = 0;
   else size_local_cols = nvalues;
 
   pack_choice = new FnPtrPack[nvalues];
 
   kindflag = NONE;
 
   int i;
   for (int iarg = 3; iarg < narg; iarg++) {
     i = iarg-3;
 
     if (strcmp(arg[iarg],"natom1") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_patom1;
       if (kindflag != NONE && kindflag != NEIGH) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = NEIGH;
     } else if (strcmp(arg[iarg],"natom2") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_patom2;
       if (kindflag != NONE && kindflag != NEIGH) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = NEIGH;
 
     } else if (strcmp(arg[iarg],"patom1") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_patom1;
       if (kindflag != NONE && kindflag != PAIR) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = PAIR;
     } else if (strcmp(arg[iarg],"patom2") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_patom2;
       if (kindflag != NONE && kindflag != PAIR) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = PAIR;
 
     } else if (strcmp(arg[iarg],"batom1") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_batom1;
       if (kindflag != NONE && kindflag != BOND) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = BOND;
     } else if (strcmp(arg[iarg],"batom2") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_batom2;
       if (kindflag != NONE && kindflag != BOND) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = BOND;
     } else if (strcmp(arg[iarg],"btype") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_btype;
       if (kindflag != NONE && kindflag != BOND) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = BOND;
 
     } else if (strcmp(arg[iarg],"aatom1") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_aatom1;
       if (kindflag != NONE && kindflag != ANGLE) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = ANGLE;
     } else if (strcmp(arg[iarg],"aatom2") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_aatom2;
       if (kindflag != NONE && kindflag != ANGLE) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = ANGLE;
     } else if (strcmp(arg[iarg],"aatom3") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_aatom3;
       if (kindflag != NONE && kindflag != ANGLE) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = ANGLE;
     } else if (strcmp(arg[iarg],"atype") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_atype;
       if (kindflag != NONE && kindflag != ANGLE) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = ANGLE;
 
     } else if (strcmp(arg[iarg],"datom1") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_datom1;
       if (kindflag != NONE && kindflag != DIHEDRAL) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = DIHEDRAL;
     } else if (strcmp(arg[iarg],"datom2") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_datom2;
       if (kindflag != NONE && kindflag != DIHEDRAL) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = DIHEDRAL;
     } else if (strcmp(arg[iarg],"datom3") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_datom3;
       if (kindflag != NONE && kindflag != DIHEDRAL) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = DIHEDRAL;
     } else if (strcmp(arg[iarg],"datom4") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_datom4;
       if (kindflag != NONE && kindflag != DIHEDRAL) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = DIHEDRAL;
     } else if (strcmp(arg[iarg],"dtype") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_dtype;
       if (kindflag != NONE && kindflag != DIHEDRAL) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = DIHEDRAL;
 
     } else if (strcmp(arg[iarg],"iatom1") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_iatom1;
       if (kindflag != NONE && kindflag != IMPROPER) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = IMPROPER;
     } else if (strcmp(arg[iarg],"iatom2") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_iatom2;
       if (kindflag != NONE && kindflag != IMPROPER) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = IMPROPER;
     } else if (strcmp(arg[iarg],"iatom3") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_iatom3;
       if (kindflag != NONE && kindflag != IMPROPER) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = IMPROPER;
     } else if (strcmp(arg[iarg],"iatom4") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_iatom4;
       if (kindflag != NONE && kindflag != IMPROPER) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = IMPROPER;
     } else if (strcmp(arg[iarg],"itype") == 0) {
       pack_choice[i] = &ComputePropertyLocal::pack_itype;
       if (kindflag != NONE && kindflag != IMPROPER) 
 	error->all("Compute property/local cannot use these inputs together");
       kindflag = IMPROPER;
 
     } else error->all("Invalid keyword in compute property/local command");
   }
 
   // error check
 
   if (kindflag == BOND && atom->avec->bonds_allow == 0)
     error->all("Compute property/local for property that isn't allocated");
   if (kindflag == ANGLE && atom->avec->angles_allow == 0)
     error->all("Compute property/local for property that isn't allocated");
   if (kindflag == DIHEDRAL && atom->avec->dihedrals_allow == 0)
     error->all("Compute property/local for property that isn't allocated");
   if (kindflag == IMPROPER && atom->avec->impropers_allow == 0)
     error->all("Compute property/local for property that isn't allocated");
 
   nmax = 0;
   vector = NULL;
   array = NULL;
   indices = NULL;
 }
 
 /* ---------------------------------------------------------------------- */
 
 ComputePropertyLocal::~ComputePropertyLocal()
 {
   delete [] pack_choice;
   memory->destroy(vector);
   memory->destroy(array);
   memory->destroy(indices);
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::init()
 {
   if (kindflag == NEIGH || kindflag == PAIR) {
     if (force->pair == NULL) 
       error->all("No pair style is defined for compute property/local");
     if (force->pair->single_enable == 0)
       error->all("Pair style does not support compute property/local");
   }
 
   // for NEIGH/PAIR need an occasional half neighbor list
 
   if (kindflag == NEIGH || kindflag == PAIR) {
     int irequest = neighbor->request((void *) this);
     neighbor->requests[irequest]->pair = 0;
     neighbor->requests[irequest]->compute = 1;
     neighbor->requests[irequest]->occasional = 1;
   }
 
   // do initial memory allocation so that memory_usage() is correct
   // cannot be done yet for NEIGH/PAIR, since neigh list does not exist
 
   if (kindflag == NEIGH) ncount = 0;
   else if (kindflag == PAIR) ncount = 0;
   else if (kindflag == BOND) ncount = count_bonds(0);
   else if (kindflag == ANGLE) ncount = count_angles(0);
   else if (kindflag == DIHEDRAL) ncount = count_dihedrals(0);
   else if (kindflag == IMPROPER) ncount = count_impropers(0);
 
   if (ncount > nmax) reallocate(ncount);
   size_local_rows = ncount;
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::init_list(int id, NeighList *ptr)
 {
   list = ptr;
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::compute_local()
 {
   invoked_local = update->ntimestep;
 
   // count local entries and generate list of indices
 
   if (kindflag == NEIGH) ncount = count_pairs(0,0);
   else if (kindflag == PAIR) ncount = count_pairs(0,1);
   else if (kindflag == BOND) ncount = count_bonds(0);
   else if (kindflag == ANGLE) ncount = count_angles(0);
   else if (kindflag == DIHEDRAL) ncount = count_dihedrals(0);
   else if (kindflag == IMPROPER) ncount = count_impropers(0);
 
   if (ncount > nmax) reallocate(ncount);
   size_local_rows = ncount;
 
   if (kindflag == NEIGH) ncount = count_pairs(1,0);
   else if (kindflag == PAIR) ncount = count_pairs(1,1);
   else if (kindflag == BOND) ncount = count_bonds(1);
   else if (kindflag == ANGLE) ncount = count_angles(1);
   else if (kindflag == DIHEDRAL) ncount = count_dihedrals(1);
   else if (kindflag == IMPROPER) ncount = count_impropers(1);
 
   // fill vector or array with local values
 
   if (nvalues == 1) {
     buf = vector;
     (this->*pack_choice[0])(0);
   } else {
     if (array) buf = &array[0][0];
     for (int n = 0; n < nvalues; n++)
       (this->*pack_choice[n])(n);
   }
 }
 
 /* ----------------------------------------------------------------------
    count pairs and compute pair info on this proc
    only count pair once if newton_pair is off
    both atom I,J must be in group
    if allflag is set, compute requested info about pair
    if forceflag = 1, pair must be within force cutoff, else neighbor cutoff
 ------------------------------------------------------------------------- */
 
 int ComputePropertyLocal::count_pairs(int allflag, int forceflag)
 {
   int i,j,m,n,ii,jj,inum,jnum,itype,jtype;
-  double xtmp,ytmp,ztmp,delx,dely,delz;
-  double rsq,factor_coul,factor_lj;
+  double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
   int *ilist,*jlist,*numneigh,**firstneigh;
 
   double **x = atom->x;
   int *tag = atom->tag;
   int *type = atom->type;
   int *mask = atom->mask;
   int nlocal = atom->nlocal;
   double *special_coul = force->special_coul;
   double *special_lj = force->special_lj;
   int newton_pair = force->newton_pair;
 
   // invoke half neighbor list (will copy or build if necessary)
 
   if (allflag == 0) neighbor->build_one(list->index);
 
   inum = list->inum;
   ilist = list->ilist;
   numneigh = list->numneigh;
   firstneigh = list->firstneigh;
 
   // loop over neighbors of my atoms
   // skip if I or J are not in group
 
   double **cutsq = force->pair->cutsq;
 
   m = n = 0;
   for (ii = 0; ii < inum; ii++) {
     i = ilist[ii];
     if (!(mask[i] & groupbit)) continue;
 
     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)];
-      factor_coul = special_coul[sbmask(j)];
       j &= NEIGHMASK;
 
       if (!(mask[j] & groupbit)) continue;
       if (newton_pair == 0 && j >= nlocal) continue;
 
       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 (forceflag && rsq >= cutsq[itype][jtype]) continue;
 
       if (allflag) {
 	indices[m][0] = tag[i];
 	indices[m][1] = tag[j];
       }
       m++;
     }
   }
 
   return m;
 }
 
 /* ----------------------------------------------------------------------
    count bonds on this proc
    only count bond once if newton_bond is off
    all atoms in interaction must be in group
    all atoms in interaction must be known to proc
    if bond is deleted (type = 0), do not count
    if bond is turned off (type < 0), still count
 ------------------------------------------------------------------------- */
 
 int ComputePropertyLocal::count_bonds(int flag)
 {
   int i,atom1,atom2;
 
   int *num_bond = atom->num_bond;
   int **bond_atom = atom->bond_atom;
   int **bond_type = atom->bond_type;
   int *tag = atom->tag;
   int *mask = atom->mask;
   int nlocal = atom->nlocal;
   int newton_bond = force->newton_bond;
 
   int m = 0;
   for (atom1 = 0; atom1 < nlocal; atom1++) {
     if (!(mask[atom1] & groupbit)) continue;
     for (i = 0; i < num_bond[atom1]; i++) {
       atom2 = atom->map(bond_atom[atom1][i]);
       if (atom2 < 0 || !(mask[atom2] & groupbit)) continue;
       if (newton_bond == 0 && tag[atom1] > tag[atom2]) continue;
       if (bond_type[atom1][i] == 0) continue;
 
       if (flag) {
 	indices[m][0] = atom1;
 	indices[m][1] = i;
       }
       m++;
     }
   }
 
   return m;
 }
 
 /* ----------------------------------------------------------------------
    count angles on this proc
    only count if 2nd atom is the one storing the angle
    all atoms in interaction must be in group
    all atoms in interaction must be known to proc
    if angle is deleted (type = 0), do not count
    if angle is turned off (type < 0), still count
 ------------------------------------------------------------------------- */
 
 int ComputePropertyLocal::count_angles(int flag)
 {
   int i,atom1,atom2,atom3;
 
   int *num_angle = atom->num_angle;
   int **angle_atom1 = atom->angle_atom1;
   int **angle_atom2 = atom->angle_atom2;
   int **angle_atom3 = atom->angle_atom3;
   int **angle_type = atom->angle_type;
   int *tag = atom->tag;
   int *mask = atom->mask;
   int nlocal = atom->nlocal;
 
   int m = 0;
   for (atom2 = 0; atom2 < nlocal; atom2++) {
     if (!(mask[atom2] & groupbit)) continue;
     for (i = 0; i < num_angle[atom2]; i++) {
       if (tag[atom2] != angle_atom2[atom2][i]) continue;
       atom1 = atom->map(angle_atom1[atom2][i]);
       if (atom1 < 0 || !(mask[atom1] & groupbit)) continue;
       atom3 = atom->map(angle_atom3[atom2][i]);
       if (atom3 < 0 || !(mask[atom3] & groupbit)) continue;
       if (angle_type[atom2][i] == 0) continue;
 
       if (flag) {
 	indices[m][0] = atom2;
 	indices[m][1] = i;
       }
       m++;
     }
   }
 
   return m;
 }
 
 /* ----------------------------------------------------------------------
    count dihedrals on this proc
    only count if 2nd atom is the one storing the dihedral
    all atoms in interaction must be in group
    all atoms in interaction must be known to proc
 ------------------------------------------------------------------------- */
 
 int ComputePropertyLocal::count_dihedrals(int flag)
 {
   int i,atom1,atom2,atom3,atom4;
 
   int *num_dihedral = atom->num_dihedral;
   int **dihedral_atom1 = atom->dihedral_atom1;
   int **dihedral_atom2 = atom->dihedral_atom2;
   int **dihedral_atom3 = atom->dihedral_atom3;
   int **dihedral_atom4 = atom->dihedral_atom4;
   int *tag = atom->tag;
   int *mask = atom->mask;
   int nlocal = atom->nlocal;
 
   int m = 0;
   for (atom2 = 0; atom2 < nlocal; atom2++) {
     if (!(mask[atom2] & groupbit)) continue;
     for (i = 0; i < num_dihedral[atom2]; i++) {
       if (tag[atom2] != dihedral_atom2[atom2][i]) continue;
       atom1 = atom->map(dihedral_atom1[atom2][i]);
       if (atom1 < 0 || !(mask[atom1] & groupbit)) continue;
       atom3 = atom->map(dihedral_atom3[atom2][i]);
       if (atom3 < 0 || !(mask[atom3] & groupbit)) continue;
       atom4 = atom->map(dihedral_atom4[atom2][i]);
       if (atom4 < 0 || !(mask[atom4] & groupbit)) continue;
 
       if (flag) {
 	indices[m][0] = atom2;
 	indices[m][1] = i;
       }
       m++;
     }
   }
 
   return m;
 }
 
 /* ----------------------------------------------------------------------
    count impropers on this proc
    only count if 2nd atom is the one storing the improper
    all atoms in interaction must be in group
    all atoms in interaction must be known to proc
 ------------------------------------------------------------------------- */
 
 int ComputePropertyLocal::count_impropers(int flag)
 {
   int i,atom1,atom2,atom3,atom4;
 
   int *num_improper = atom->num_improper;
   int **improper_atom1 = atom->improper_atom1;
   int **improper_atom2 = atom->improper_atom2;
   int **improper_atom3 = atom->improper_atom3;
   int **improper_atom4 = atom->improper_atom4;
   int *tag = atom->tag;
   int *mask = atom->mask;
   int nlocal = atom->nlocal;
 
   int m = 0;
   for (atom2 = 0; atom2 < nlocal; atom2++) {
     if (!(mask[atom2] & groupbit)) continue;
     for (i = 0; i < num_improper[atom2]; i++) {
       if (tag[atom2] != improper_atom2[atom2][i]) continue;
       atom1 = atom->map(improper_atom1[atom2][i]);
       if (atom1 < 0 || !(mask[atom1] & groupbit)) continue;
       atom3 = atom->map(improper_atom3[atom2][i]);
       if (atom3 < 0 || !(mask[atom3] & groupbit)) continue;
       atom4 = atom->map(improper_atom4[atom2][i]);
       if (atom4 < 0 || !(mask[atom4] & groupbit)) continue;
 
       if (flag) {
 	indices[m][0] = atom2;
 	indices[m][1] = i;
       }
       m++;
     }
   }
 
   return m;
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::reallocate(int n)
 {
   // grow vector or array and indices array
 
   while (nmax < n) nmax += DELTA;
   if (nvalues == 1) {
     memory->destroy(vector);
     memory->create(vector,nmax,"property/local:vector");
     vector_local = vector;
   } else {
     memory->destroy(array);
     memory->create(array,nmax,nvalues,"property/local:array");
     array_local = array;
   }
 
   memory->destroy(indices);
   memory->create(indices,nmax,2,"property/local:indices");
 }
 
 /* ----------------------------------------------------------------------
    memory usage of local data
 ------------------------------------------------------------------------- */
 
 double ComputePropertyLocal::memory_usage()
 {
   double bytes = nmax*nvalues * sizeof(double);
   bytes += nmax*2 * sizeof(int);
   return bytes;
 }
 
 /* ----------------------------------------------------------------------
    one method for every keyword compute property/local can output
    the atom property is packed into buf starting at n with stride nvalues
    customize a new keyword by adding a method
 ------------------------------------------------------------------------- */
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_patom1(int n)
 {
   for (int m = 0; m < ncount; m++) {
     buf[n] = indices[m][0];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_patom2(int n)
 {
   for (int m = 0; m < ncount; m++) {
     buf[n] = indices[m][1];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_batom1(int n)
 {
   int i;
   int *tag = atom->tag;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     buf[n] = tag[i];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_batom2(int n)
 {
   int i,j;
   int **bond_atom = atom->bond_atom;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = bond_atom[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_btype(int n)
 {
   int i,j;
   int **bond_type = atom->bond_type;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = bond_type[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_aatom1(int n)
 {
   int i,j;
   int **angle_atom1 = atom->angle_atom1;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = angle_atom1[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_aatom2(int n)
 {
   int i,j;
   int **angle_atom2 = atom->angle_atom2;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = angle_atom2[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_aatom3(int n)
 {
   int i,j;
   int **angle_atom3 = atom->angle_atom3;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = angle_atom3[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_atype(int n)
 {
   int i,j;
   int **angle_type = atom->angle_type;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = angle_type[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_datom1(int n)
 {
   int i,j;
   int **dihedral_atom1 = atom->dihedral_atom1;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = dihedral_atom1[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_datom2(int n)
 {
   int i,j;
   int **dihedral_atom2 = atom->dihedral_atom2;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = dihedral_atom2[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_datom3(int n)
 {
   int i,j;
   int **dihedral_atom3 = atom->dihedral_atom3;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = dihedral_atom3[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_datom4(int n)
 {
   int i,j;
   int **dihedral_atom4 = atom->dihedral_atom4;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = dihedral_atom4[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_dtype(int n)
 {
   int i,j;
   int **dihedral_type = atom->dihedral_type;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = dihedral_type[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_iatom1(int n)
 {
   int i,j;
   int **improper_atom1 = atom->improper_atom1;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = improper_atom1[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_iatom2(int n)
 {
   int i,j;
   int **improper_atom2 = atom->improper_atom2;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = improper_atom2[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_iatom3(int n)
 {
   int i,j;
   int **improper_atom3 = atom->improper_atom3;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = improper_atom3[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_iatom4(int n)
 {
   int i,j;
   int **improper_atom4 = atom->improper_atom4;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = improper_atom4[i][j];
     n += nvalues;
   }
 }
 
 /* ---------------------------------------------------------------------- */
 
 void ComputePropertyLocal::pack_itype(int n)
 {
   int i,j;
   int **improper_type = atom->improper_type;
 
   for (int m = 0; m < ncount; m++) {
     i = indices[m][0];
     j = indices[m][1];
     buf[n] = improper_type[i][j];
     n += nvalues;
   }
 }
diff --git a/src/memory.cpp b/src/memory.cpp
index 541a7101d..6af16e68e 100644
--- a/src/memory.cpp
+++ b/src/memory.cpp
@@ -1,86 +1,85 @@
 /* ----------------------------------------------------------------------
    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 "lmptype.h"
-#include "mpi.h"
 #include "stdio.h"
 #include "stdlib.h"
 #include "string.h"
 #include "memory.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
 
 /* ---------------------------------------------------------------------- */
 
 Memory::Memory(LAMMPS *lmp) : Pointers(lmp) {}
 
 /* ----------------------------------------------------------------------
    safe malloc 
 ------------------------------------------------------------------------- */
 
 void *Memory::smalloc(bigint nbytes, const char *name)
 {
   if (nbytes == 0) return NULL;
 
   void *ptr = malloc(nbytes);
   if (ptr == NULL) {
     char str[128];
     sprintf(str,"Failed to allocate " BIGINT_FORMAT " bytes for array %s",
 	    nbytes,name);
     error->one(str);
   }
   return ptr;
 }
 
 /* ----------------------------------------------------------------------
    safe realloc 
 ------------------------------------------------------------------------- */
 
 void *Memory::srealloc(void *ptr, bigint nbytes, const char *name)
 {
   if (nbytes == 0) {
     destroy(ptr);
     return NULL;
   }
 
   ptr = realloc(ptr,nbytes);
   if (ptr == NULL) {
     char str[128];
     sprintf(str,"Failed to reallocate " BIGINT_FORMAT " bytes for array %s",
 	    nbytes,name);
     error->one(str);
   }
   return ptr;
 }
 
 /* ----------------------------------------------------------------------
    safe free 
 ------------------------------------------------------------------------- */
 
 void Memory::sfree(void *ptr)
 {
   if (ptr == NULL) return;
   free(ptr);
 }
 
 /* ----------------------------------------------------------------------
    erroneous usage of templated create/grow functions
 ------------------------------------------------------------------------- */
 
 void Memory::fail(const char *name)
 {
   char str[128];
   sprintf(str,"Cannot create/grow a vector/array of pointers for %s",name);
   error->one(str);
 }