diff --git a/src/region.h b/src/region.h
index 8f7959025..0a1d888a9 100644
--- a/src/region.h
+++ b/src/region.h
@@ -1,116 +1,115 @@
 /* ----------------------------------------------------------------------
    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.
 ------------------------------------------------------------------------- */
 
 #ifndef LMP_REGION_H
 #define LMP_REGION_H
 
 #include "pointers.h"
 
 namespace LAMMPS_NS {
 
 class Region : protected Pointers {
  public:
   char *id,*style;
   int interior;                     // 1 for interior, 0 for exterior
   int scaleflag;                    // 1 for lattice, 0 for box
   double xscale,yscale,zscale;      // scale factors for box/lattice units
   double extent_xlo,extent_xhi;     // bounding box on region
   double extent_ylo,extent_yhi;
   double extent_zlo,extent_zhi;
   int bboxflag;                     // 1 if bounding box is computable
+  int varshape;                     // 1 if region shape changes over time
 
   // contact = particle near region surface
 
   struct Contact {
     double r;                 // distance between particle & surf, r > 0.0
     double delx,dely,delz;    // vector from surface pt to particle
   };
   Contact *contact;           // list of contacts
   int cmax;                   // max # of contacts possible with region
 
   Region(class LAMMPS *, int, char **);
   virtual ~Region();
   virtual void init();
   virtual int dynamic_check();
 
   // called by other classes to check point versus region
 
   int match(double, double, double);
   int surface(double, double, double, double);
 
   // implemented by each region, not called by other classes
 
   virtual int inside(double, double, double) = 0;
   virtual int surface_interior(double *, double) = 0;
   virtual int surface_exterior(double *, double) = 0;
   virtual void shape_update() {}
 
  protected:
-  int varshape;       // 1 if region shape changes over time
-
   void add_contact(int, double *, double, double, double);
   void options(int, char **);
 
  private:
   int dynamic;        // 1 if region position/orientation changes over time
   int moveflag,rotateflag;   // 1 if position/orientation changes
 
   double point[3],axis[3],runit[3];
   char *xstr,*ystr,*zstr,*tstr;
   int xvar,yvar,zvar,tvar;
   double dx,dy,dz,theta;
   bigint lastshape,lastdynamic;
 
   void forward_transform(double &, double &, double &);
   void inverse_transform(double &, double &, double &);
   void rotate(double &, double &, double &, double);
 };
 
 }
 
 #endif
 
 /* ERROR/WARNING messages:
 
 E: Variable name for region does not exist
 
 Self-explanatory.
 
 E: Variable for region is invalid style
 
 Only equal-style variables can be used.
 
 E: Variable for region is not equal style
 
 Self-explanatory.
 
 E: Illegal ... command
 
 Self-explanatory.  Check the input script syntax and compare to the
 documentation for the command.  You can use -echo screen as a
 command-line option when running LAMMPS to see the offending line.
 
 E: Region union or intersect cannot be dynamic
 
 The sub-regions can be dynamic, but not the combined region.
 
 E: Use of region with undefined lattice
 
 If units = lattice (the default) for the region command, then a
 lattice must first be defined via the lattice command.
 
 E: Region cannot have 0 length rotation vector
 
 Self-explanatory.
 
 */
diff --git a/src/region_intersect.cpp b/src/region_intersect.cpp
index 736ec3c6e..bae73952d 100644
--- a/src/region_intersect.cpp
+++ b/src/region_intersect.cpp
@@ -1,208 +1,233 @@
 /* ----------------------------------------------------------------------
    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 "stdlib.h"
 #include "string.h"
 #include "region_intersect.h"
 #include "domain.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
 
 /* ---------------------------------------------------------------------- */
 
 RegIntersect::RegIntersect(LAMMPS *lmp, int narg, char **arg) :
   Region(lmp, narg, arg)
 {
   if (narg < 5) error->all(FLERR,"Illegal region command");
   int n = atoi(arg[2]);
   if (n < 2) error->all(FLERR,"Illegal region command");
   options(narg-(n+3),&arg[n+3]);
 
   // build list of regions to intersect
 
   list = new int[n];
   nregion = 0;
 
   int iregion;
   for (int iarg = 0; iarg < n; iarg++) {
     iregion = domain->find_region(arg[iarg+3]);
     if (iregion == -1) 
       error->all(FLERR,"Region intersect region ID does not exist");
     list[nregion++] = iregion;
   }
 
-  // extent of intersection of regions
-  // has bounding box if interior and any sub-region has bounding box
+  // this region is variable shape if any of sub-regions are
 
   Region **regions = domain->regions;
+  for (int ilist = 0; ilist < nregion; ilist++)
+    if (regions[list[ilist]]->varshape) varshape = 1;
+
+  // extent of intersection of regions
+  // has bounding box if interior and any sub-region has bounding box
 
   bboxflag = 0;
   for (int ilist = 0; ilist < nregion; ilist++)
     if (regions[list[ilist]]->bboxflag == 1) bboxflag = 1;
   if (!interior) bboxflag = 0;
 
   if (bboxflag) {
     int first = 1;
     for (int ilist = 0; ilist < nregion; ilist++) {
       if (regions[list[ilist]]->bboxflag == 0) continue;
       if (first) {
         extent_xlo = regions[list[ilist]]->extent_xlo;
         extent_ylo = regions[list[ilist]]->extent_ylo;
         extent_zlo = regions[list[ilist]]->extent_zlo;
         extent_xhi = regions[list[ilist]]->extent_xhi;
         extent_yhi = regions[list[ilist]]->extent_yhi;
         extent_zhi = regions[list[ilist]]->extent_zhi;
         first = 0;
       }
 
       extent_xlo = MAX(extent_xlo,regions[list[ilist]]->extent_xlo);
       extent_ylo = MAX(extent_ylo,regions[list[ilist]]->extent_ylo);
       extent_zlo = MAX(extent_zlo,regions[list[ilist]]->extent_zlo);
       extent_xhi = MIN(extent_xhi,regions[list[ilist]]->extent_xhi);
       extent_yhi = MIN(extent_yhi,regions[list[ilist]]->extent_yhi);
       extent_zhi = MIN(extent_zhi,regions[list[ilist]]->extent_zhi);
     }
   }
 
   // possible contacts = sum of possible contacts in all sub-regions
 
   cmax = 0;
   for (int ilist = 0; ilist < nregion; ilist++)
     cmax += regions[list[ilist]]->cmax;
   contact = new Contact[cmax];
 }
 
 /* ---------------------------------------------------------------------- */
 
 RegIntersect::~RegIntersect()
 {
   delete [] list;
   delete [] contact;
 }
 
+/* ---------------------------------------------------------------------- */
+
+void RegIntersect::init()
+{
+  Region::init();
+  Region **regions = domain->regions;
+  for (int ilist = 0; ilist < nregion; ilist++)
+    regions[list[ilist]]->init();
+}
+
 /* ----------------------------------------------------------------------
    return 1 if region is dynamic, 0 if static
    dynamic if any sub-region is dynamic, else static
 ------------------------------------------------------------------------- */
 
 int RegIntersect::dynamic_check()
 {
   Region **regions = domain->regions;
   for (int ilist = 0; ilist < nregion; ilist++)
     if (regions[list[ilist]]->dynamic_check()) return 1;
   return 0;
 }
 
 /* ----------------------------------------------------------------------
    inside = 1 if x,y,z is match() with all sub-regions
    else inside = 0
 ------------------------------------------------------------------------- */
 
 int RegIntersect::inside(double x, double y, double z)
 {
   int ilist;
   Region **regions = domain->regions;
   for (ilist = 0; ilist < nregion; ilist++)
     if (!regions[list[ilist]]->match(x,y,z)) break;
 
   if (ilist == nregion) return 1;
   return 0;
 }
 
 /* ----------------------------------------------------------------------
    compute contacts with interior of intersection of sub-regions
    (1) compute contacts in each sub-region
    (2) only keep a contact if surface point is match() to all other regions
 ------------------------------------------------------------------------- */
 
 int RegIntersect::surface_interior(double *x, double cutoff)
 {
   int m,ilist,jlist,iregion,jregion,ncontacts;
   double xs,ys,zs;
 
   Region **regions = domain->regions;
   int n = 0;
 
   for (ilist = 0; ilist < nregion; ilist++) {
     iregion = list[ilist];
     ncontacts = regions[iregion]->surface(x[0],x[1],x[2],cutoff);
     for (m = 0; m < ncontacts; m++) {
       xs = x[0] - regions[iregion]->contact[m].delx;
       ys = x[1] - regions[iregion]->contact[m].dely;
       zs = x[2] - regions[iregion]->contact[m].delz;
       for (jlist = 0; jlist < nregion; jlist++) {
         if (jlist == ilist) continue;
         jregion = list[jlist];
         if (!regions[jregion]->match(xs,ys,zs)) break;
       }
       if (jlist == nregion) {
         contact[n].r = regions[iregion]->contact[m].r;
         contact[n].delx = regions[iregion]->contact[m].delx;
         contact[n].dely = regions[iregion]->contact[m].dely;
         contact[n].delz = regions[iregion]->contact[m].delz;
         n++;
       }
     }
   }
 
   return n;
 }
 
 /* ----------------------------------------------------------------------
    compute contacts with interior of intersection of sub-regions
    (1) flip interior/exterior flag of each sub-region
    (2) compute contacts in each sub-region
    (3) only keep a contact if surface point is not match() to all other regions
    (4) flip interior/exterior flags back to original settings
    this is effectively same algorithm as surface_interior() for RegUnion
 ------------------------------------------------------------------------- */
 
 int RegIntersect::surface_exterior(double *x, double cutoff)
 {
   int m,ilist,jlist,iregion,jregion,ncontacts;
   double xs,ys,zs;
 
   Region **regions = domain->regions;
   int n = 0;
 
   for (ilist = 0; ilist < nregion; ilist++)
     regions[list[ilist]]->interior ^= 1;
 
   for (ilist = 0; ilist < nregion; ilist++) {
     iregion = list[ilist];
     ncontacts = regions[iregion]->surface(x[0],x[1],x[2],cutoff);
     for (m = 0; m < ncontacts; m++) {
       xs = x[0] - regions[iregion]->contact[m].delx;
       ys = x[1] - regions[iregion]->contact[m].dely;
       zs = x[2] - regions[iregion]->contact[m].delz;
       for (jlist = 0; jlist < nregion; jlist++) {
         if (jlist == ilist) continue;
         jregion = list[jlist];
         if (regions[jregion]->match(xs,ys,zs)) break;
       }
       if (jlist == nregion) {
         contact[n].r = regions[iregion]->contact[m].r;
         contact[n].delx = regions[iregion]->contact[m].delx;
         contact[n].dely = regions[iregion]->contact[m].dely;
         contact[n].delz = regions[iregion]->contact[m].delz;
         n++;
       }
     }
   }
 
   for (ilist = 0; ilist < nregion; ilist++)
     regions[list[ilist]]->interior ^= 1;
 
   return n;
 }
+
+/* ----------------------------------------------------------------------
+   change region shape of all sub-regions
+------------------------------------------------------------------------- */
+
+void RegIntersect::shape_update()
+{
+  Region **regions = domain->regions;
+  for (int ilist = 0; ilist < nregion; ilist++)
+    regions[list[ilist]]->shape_update();
+}
diff --git a/src/region_intersect.h b/src/region_intersect.h
index a8bd460f3..62946a65c 100644
--- a/src/region_intersect.h
+++ b/src/region_intersect.h
@@ -1,58 +1,60 @@
 /* ----------------------------------------------------------------------
    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.
 ------------------------------------------------------------------------- */
 
 #ifdef REGION_CLASS
 
 RegionStyle(intersect,RegIntersect)
 
 #else
 
 #ifndef LMP_REGION_INTERSECT_H
 #define LMP_REGION_INTERSECT_H
 
 #include "region.h"
 
 namespace LAMMPS_NS {
 
 class RegIntersect : public Region {
  public:
   RegIntersect(class LAMMPS *, int, char **);
   ~RegIntersect();
+  void init();
   int dynamic_check();
   int inside(double, double, double);
   int surface_interior(double *, double);
   int surface_exterior(double *, double);
+  void shape_update();
 
  private:
   int nregion;
   int *list;
 };
 
 }
 
 #endif
 #endif
 
 /* ERROR/WARNING messages:
 
 E: Illegal ... command
 
 Self-explanatory.  Check the input script syntax and compare to the
 documentation for the command.  You can use -echo screen as a
 command-line option when running LAMMPS to see the offending line.
 
 E: Region intersect region ID does not exist
 
 Self-explanatory.
 
 */
diff --git a/src/region_union.cpp b/src/region_union.cpp
index ca1e693f3..3741e8fdd 100644
--- a/src/region_union.cpp
+++ b/src/region_union.cpp
@@ -1,199 +1,225 @@
 /* ----------------------------------------------------------------------
    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 "stdlib.h"
 #include "string.h"
 #include "region_union.h"
 #include "domain.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
 
 #define BIG 1.0e20
 
 /* ---------------------------------------------------------------------- */
 
 RegUnion::RegUnion(LAMMPS *lmp, int narg, char **arg) : Region(lmp, narg, arg)
 {
   if (narg < 5) error->all(FLERR,"Illegal region command");
   int n = atoi(arg[2]);
   if (n < 2) error->all(FLERR,"Illegal region command");
   options(narg-(n+3),&arg[n+3]);
 
   // build list of regions to union
 
   list = new int[n];
   nregion = 0;
 
   int iregion;
   for (int iarg = 0; iarg < n; iarg++) {
     iregion = domain->find_region(arg[iarg+3]);
-    if (iregion == -1) error->all(FLERR,"Region union region ID does not exist");
+    if (iregion == -1) 
+      error->all(FLERR,"Region union region ID does not exist");
     list[nregion++] = iregion;
   }
 
-  // extent of union of regions
-  // has bounding box if interior and all sub-regions have bounding box
+  // this region is variable shape if any of sub-regions are
 
   Region **regions = domain->regions;
+  for (int ilist = 0; ilist < nregion; ilist++)
+    if (regions[list[ilist]]->varshape) varshape = 1;
+
+  // extent of union of regions
+  // has bounding box if interior and all sub-regions have bounding box
 
   bboxflag = 1;
   for (int ilist = 0; ilist < nregion; ilist++)
     if (regions[list[ilist]]->bboxflag == 0) bboxflag = 0;
   if (!interior) bboxflag = 0;
 
   if (bboxflag) {
     extent_xlo = extent_ylo = extent_zlo = BIG;
     extent_xhi = extent_yhi = extent_zhi = -BIG;
 
     for (int ilist = 0; ilist < nregion; ilist++) {
       extent_xlo = MIN(extent_xlo,regions[list[ilist]]->extent_xlo);
       extent_ylo = MIN(extent_ylo,regions[list[ilist]]->extent_ylo);
       extent_zlo = MIN(extent_zlo,regions[list[ilist]]->extent_zlo);
       extent_xhi = MAX(extent_xhi,regions[list[ilist]]->extent_xhi);
       extent_yhi = MAX(extent_yhi,regions[list[ilist]]->extent_yhi);
       extent_zhi = MAX(extent_zhi,regions[list[ilist]]->extent_zhi);
     }
   }
 
   // possible contacts = sum of possible contacts in all sub-regions
 
   cmax = 0;
   for (int ilist = 0; ilist < nregion; ilist++)
     cmax += regions[list[ilist]]->cmax;
   contact = new Contact[cmax];
 }
 
 /* ---------------------------------------------------------------------- */
 
 RegUnion::~RegUnion()
 {
   delete [] list;
   delete [] contact;
 }
 
+/* ---------------------------------------------------------------------- */
+
+void RegUnion::init()
+{
+  Region::init();
+  Region **regions = domain->regions;
+  for (int ilist = 0; ilist < nregion; ilist++)
+    regions[list[ilist]]->init();
+}
+
 /* ----------------------------------------------------------------------
    return 1 if region is dynamic, 0 if static
    dynamic if any sub-region is dynamic, else static
 ------------------------------------------------------------------------- */
 
 int RegUnion::dynamic_check()
 {
   Region **regions = domain->regions;
   for (int ilist = 0; ilist < nregion; ilist++)
     if (regions[list[ilist]]->dynamic_check()) return 1;
   return 0;
 }
 
 /* ----------------------------------------------------------------------
    inside = 1 if x,y,z is match() with any sub-region
    else inside = 0
 ------------------------------------------------------------------------- */
 
 int RegUnion::inside(double x, double y, double z)
 {
   int ilist;
   Region **regions = domain->regions;
   for (ilist = 0; ilist < nregion; ilist++)
     if (regions[list[ilist]]->match(x,y,z)) break;
 
   if (ilist == nregion) return 0;
   return 1;
 }
 
 /* ----------------------------------------------------------------------
    compute contacts with interior of union of sub-regions
    (1) compute contacts in each sub-region
    (2) only keep a contact if surface point is not match() to all other regions
 ------------------------------------------------------------------------- */
 
 int RegUnion::surface_interior(double *x, double cutoff)
 {
   int m,ilist,jlist,iregion,jregion,ncontacts;
   double xs,ys,zs;
 
   Region **regions = domain->regions;
   int n = 0;
 
   for (ilist = 0; ilist < nregion; ilist++) {
     iregion = list[ilist];
     ncontacts = regions[iregion]->surface(x[0],x[1],x[2],cutoff);
     for (m = 0; m < ncontacts; m++) {
       xs = x[0] - regions[iregion]->contact[m].delx;
       ys = x[1] - regions[iregion]->contact[m].dely;
       zs = x[2] - regions[iregion]->contact[m].delz;
       for (jlist = 0; jlist < nregion; jlist++) {
         if (jlist == ilist) continue;
         jregion = list[jlist];
         if (regions[jregion]->match(xs,ys,zs)) break;
       }
       if (jlist == nregion) {
         contact[n].r = regions[iregion]->contact[m].r;
         contact[n].delx = regions[iregion]->contact[m].delx;
         contact[n].dely = regions[iregion]->contact[m].dely;
         contact[n].delz = regions[iregion]->contact[m].delz;
         n++;
       }
     }
   }
 
   return n;
 }
 
 /* ----------------------------------------------------------------------
    compute contacts with exterior of union of sub-regions
    (1) flip interior/exterior flag of each sub-region
    (2) compute contacts in each sub-region
    (3) only keep a contact if surface point is match() to all other regions
    (4) flip interior/exterior flags back to original settings
    this is effectively same algorithm as surface_interior() for RegIntersect
 ------------------------------------------------------------------------- */
 
 int RegUnion::surface_exterior(double *x, double cutoff)
 {
   int m,ilist,jlist,iregion,jregion,ncontacts;
   double xs,ys,zs;
 
   Region **regions = domain->regions;
   int n = 0;
 
   for (ilist = 0; ilist < nregion; ilist++)
     regions[list[ilist]]->interior ^= 1;
 
   for (ilist = 0; ilist < nregion; ilist++) {
     iregion = list[ilist];
     ncontacts = regions[iregion]->surface(x[0],x[1],x[2],cutoff);
     for (m = 0; m < ncontacts; m++) {
       xs = x[0] - regions[iregion]->contact[m].delx;
       ys = x[1] - regions[iregion]->contact[m].dely;
       zs = x[2] - regions[iregion]->contact[m].delz;
       for (jlist = 0; jlist < nregion; jlist++) {
         if (jlist == ilist) continue;
         jregion = list[jlist];
         if (!regions[jregion]->match(xs,ys,zs)) break;
       }
       if (jlist == nregion) {
         contact[n].r = regions[iregion]->contact[m].r;
         contact[n].delx = regions[iregion]->contact[m].delx;
         contact[n].dely = regions[iregion]->contact[m].dely;
         contact[n].delz = regions[iregion]->contact[m].delz;
         n++;
       }
     }
   }
 
   for (ilist = 0; ilist < nregion; ilist++)
     regions[list[ilist]]->interior ^= 1;
 
   return n;
 }
+
+/* ----------------------------------------------------------------------
+   change region shape of all sub-regions
+------------------------------------------------------------------------- */
+
+void RegUnion::shape_update()
+{
+  Region **regions = domain->regions;
+  for (int ilist = 0; ilist < nregion; ilist++)
+    regions[list[ilist]]->shape_update();
+}
diff --git a/src/region_union.h b/src/region_union.h
index 44d99375c..6d0ec063c 100644
--- a/src/region_union.h
+++ b/src/region_union.h
@@ -1,59 +1,61 @@
 /* ----------------------------------------------------------------------
    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.
 ------------------------------------------------------------------------- */
 
 #ifdef REGION_CLASS
 
 RegionStyle(union,RegUnion)
 
 #else
 
 #ifndef LMP_REGION_UNION_H
 #define LMP_REGION_UNION_H
 
 #include "region.h"
 
 namespace LAMMPS_NS {
 
 class RegUnion : public Region {
  public:
   RegUnion(class LAMMPS *, int, char **);
   ~RegUnion();
+  void init();
   int dynamic_check();
   int inside(double, double, double);
   int surface_interior(double *, double);
   int surface_exterior(double *, double);
+  void shape_update();
 
  private:
   int nregion;
   int *list;
 };
 
 }
 
 #endif
 #endif
 
 /* ERROR/WARNING messages:
 
 E: Illegal ... command
 
 Self-explanatory.  Check the input script syntax and compare to the
 documentation for the command.  You can use -echo screen as a
 command-line option when running LAMMPS to see the offending line.
 
 E: Region union region ID does not exist
 
 One or more of the region IDs specified by the region union command
 does not exist.
 
 */