diff --git a/src/MANYBODY/pair_airebo.cpp b/src/MANYBODY/pair_airebo.cpp
index fa01ac15ac95c2340503f452a324c8de3b4a5c5e..f6555063dc42ca3b261cd257e111842214a92c97 100644
--- a/src/MANYBODY/pair_airebo.cpp
+++ b/src/MANYBODY/pair_airebo.cpp
@@ -31,6 +31,7 @@
 #include "neighbor.h"
 #include "neigh_list.h"
 #include "neigh_request.h"
+#include "my_page.h"
 #include "math_const.h"
 #include "math_special.h"
 #include "memory.h"
@@ -55,8 +56,9 @@ PairAIREBO::PairAIREBO(LAMMPS *lmp) : Pair(lmp)
   maxlocal = 0;
   REBO_numneigh = NULL;
   REBO_firstneigh = NULL;
-  maxpage = 0;
-  pages = NULL;
+  ipage = NULL;
+  pgsize = oneatom = 0;
+
   nC = nH = NULL;
   manybody_flag = 1;
 }
@@ -69,8 +71,7 @@ PairAIREBO::~PairAIREBO()
 {
   memory->destroy(REBO_numneigh);
   memory->sfree(REBO_firstneigh);
-  for (int i = 0; i < maxpage; i++) memory->destroy(pages[i]);
-  memory->sfree(pages);
+  delete [] ipage;
   memory->destroy(nC);
   memory->destroy(nH);
 
@@ -221,11 +222,24 @@ void PairAIREBO::init_style()
   neighbor->requests[irequest]->full = 1;
   neighbor->requests[irequest]->ghost = 1;
 
-  // local REBO neighbor list memory
+  // local REBO neighbor list
+  // create pages if first time or if neighbor pgsize/oneatom has changed
+
+  int create = 0;
+  if (ipage == NULL) create = 1;
+  if (pgsize != neighbor->pgsize) create = 1;
+  if (oneatom != neighbor->oneatom) create = 1;
 
-  pgsize = neighbor->pgsize;
-  oneatom = neighbor->oneatom;
-  if (maxpage == 0) add_pages();
+  if (create) {
+    delete [] ipage;
+    pgsize = neighbor->pgsize;
+    oneatom = neighbor->oneatom;
+
+    int nmypage= comm->nthreads;
+    ipage = new MyPage<int>[nmypage];
+    for (int i = 0; i < nmypage; i++)
+      ipage[i].init(oneatom,pgsize,PGDELTA);
+  }
 }
 
 /* ----------------------------------------------------------------------
@@ -329,19 +343,13 @@ void PairAIREBO::REBO_neigh()
   // store all REBO neighs of owned and ghost atoms
   // scan full neighbor list of I
 
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (ii = 0; ii < allnum; ii++) {
     i = ilist[ii];
 
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == maxpage) add_pages();
-    }
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     xtmp = x[i][0];
     ytmp = x[i][1];
@@ -371,10 +379,9 @@ void PairAIREBO::REBO_neigh()
 
     REBO_firstneigh[i] = neighptr;
     REBO_numneigh[i] = n;
-    npnt += n;
-    if (npnt >= pgsize)
-      error->one(FLERR,
-                 "Neighbor list overflow, boost neigh_modify one or page");
+    ipage->vgot(n);
+    if (ipage->status())
+      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 }
 
@@ -3285,21 +3292,6 @@ double PairAIREBO::TijSpline(double Nij, double Nji,
   return Tijf;
 }
 
-/* ----------------------------------------------------------------------
-   add pages to REBO neighbor list
-------------------------------------------------------------------------- */
-
-void PairAIREBO::add_pages(int howmany)
-{
-  int toppage = maxpage;
-  maxpage += howmany*PGDELTA;
-
-  pages = (int **)
-    memory->srealloc(pages,maxpage*sizeof(int *),"AIREBO:pages");
-  for (int i = toppage; i < maxpage; i++)
-    memory->create(pages[i],pgsize,"AIREBO:pages[i]");
-}
-
 /* ----------------------------------------------------------------------
    read AIREBO potential file
 ------------------------------------------------------------------------- */
@@ -4199,7 +4191,10 @@ double PairAIREBO::memory_usage()
   double bytes = 0.0;
   bytes += maxlocal * sizeof(int);
   bytes += maxlocal * sizeof(int *);
-  bytes += maxpage * neighbor->pgsize * sizeof(int);
-  bytes += 3 * maxlocal * sizeof(double);
+
+  for (int i = 0; i < comm->nthreads; i++)
+    bytes += ipage[i].size();
+  
+  bytes += 2*maxlocal * sizeof(double);
   return bytes;
 }
diff --git a/src/MANYBODY/pair_airebo.h b/src/MANYBODY/pair_airebo.h
index fc8378b98be312fc2b9279714ec0b0cb897b747a..8536d6ab012e03772e569f5010bc47b140c712bf 100644
--- a/src/MANYBODY/pair_airebo.h
+++ b/src/MANYBODY/pair_airebo.h
@@ -21,6 +21,7 @@ PairStyle(airebo,PairAIREBO)
 #define LMP_PAIR_AIREBO_H
 
 #include "pair.h"
+#include "my_page.h"
 #include "math.h"
 #include "math_const.h"
 
@@ -38,15 +39,10 @@ class PairAIREBO : public Pair {
   double memory_usage();
 
  protected:
-  int **pages;                     // neighbor list pages
   int *map;                        // 0 (C), 1 (H), or -1 (NULL) for each type
 
   int me;
   int ljflag,torflag;              // 0/1 if LJ,torsion terms included
-  int maxlocal;                    // size of numneigh, firstneigh arrays
-  int maxpage;                     // # of pages currently allocated
-  int pgsize;                      // size of neighbor page
-  int oneatom;                     // max # of neighbors for one atom
 
   double cutlj;                    // user-specified LJ cutoff
   double cutljrebosq;              // cut for when to compute
@@ -56,8 +52,13 @@ class PairAIREBO : public Pair {
   double **lj1,**lj2,**lj3,**lj4;  // pre-computed LJ coeffs for C,H types
   double cut3rebo;                 // maximum distance for 3rd REBO neigh
 
+  int maxlocal;                    // size of numneigh, firstneigh arrays
+  int pgsize;                      // size of neighbor page
+  int oneatom;                     // max # of neighbors for one atom
+  MyPage<int> *ipage;              // neighbor list pages
   int *REBO_numneigh;              // # of pair neighbors for each atom
   int **REBO_firstneigh;           // ptr to 1st neighbor of each atom
+
   double *closestdistsq;           // closest owned atom dist to each ghost
   double *nC,*nH;                  // sum of weighting fns with REBO neighs
 
@@ -101,7 +102,6 @@ class PairAIREBO : public Pair {
   double piRCSpline(double, double, double, int, int, double *);
   double TijSpline(double, double, double, double *);
 
-  void add_pages(int howmany = 1);
   void read_file(char *);
 
   double Sp5th(double, double *, double *);
diff --git a/src/MANYBODY/pair_comb.cpp b/src/MANYBODY/pair_comb.cpp
index d5982caf0eb068dd5ce8497fececf77802ba58f7..2da2fc93d96f3144e52e5ccacd71e71db70736c0 100644
--- a/src/MANYBODY/pair_comb.cpp
+++ b/src/MANYBODY/pair_comb.cpp
@@ -31,6 +31,7 @@
 #include "neigh_request.h"
 #include "group.h"
 #include "update.h"
+#include "my_page.h"
 #include "math_const.h"
 #include "memory.h"
 #include "error.h"
@@ -73,8 +74,9 @@ PairComb::PairComb(LAMMPS *lmp) : Pair(lmp)
 
   sht_num = NULL;
   sht_first = NULL;
-  maxpage = 0;
-  pages = NULL;
+
+  ipage = NULL;
+  pgsize = oneatom = 0;
 
   // set comm size needed by this Pair
 
@@ -108,6 +110,8 @@ PairComb::~PairComb()
   memory->destroy(sht_num);
   memory->destroy(sht_first);
 
+  delete [] ipage;
+
   if (allocated) {
     memory->destroy(setflag);
     memory->destroy(cutsq);
@@ -541,9 +545,24 @@ void PairComb::init_style()
   neighbor->requests[irequest]->half = 0;
   neighbor->requests[irequest]->full = 1;
 
-  pgsize = neighbor->pgsize;
-  oneatom = neighbor->oneatom;
-  if (maxpage == 0) add_pages();
+  // local Comb neighbor list
+  // create pages if first time or if neighbor pgsize/oneatom has changed
+
+  int create = 0;
+  if (ipage == NULL) create = 1;
+  if (pgsize != neighbor->pgsize) create = 1;
+  if (oneatom != neighbor->oneatom) create = 1;
+
+  if (create) {
+    delete [] ipage;
+    pgsize = neighbor->pgsize;
+    oneatom = neighbor->oneatom;
+
+    int nmypage = comm->nthreads;
+    ipage = new MyPage<int>[nmypage];
+    for (int i = 0; i < nmypage; i++)
+      ipage[i].init(oneatom,pgsize);
+  }
 }
 
 /* ----------------------------------------------------------------------
@@ -1364,8 +1383,7 @@ void PairComb::sm_table()
   memory->create(NCo,nmax,"pair:NCo");
   memory->create(bbij,nmax,MAXNEIGH,"pair:bbij");
   memory->create(sht_num,nmax,"pair:sht_num");
-  sht_first = (int **) memory->smalloc(nmax*sizeof(int *),
-            "pair:sht_first");
+  sht_first = (int **) memory->smalloc(nmax*sizeof(int *),"pair:sht_first");
 
   // set interaction number: 0-0=0, 1-1=1, 0-1=1-0=2
 
@@ -2023,18 +2041,6 @@ void PairComb::unpack_reverse_comm(int n, int *list, double *buf)
   }
 }
 
-/* ----------------------------------------------------------------------
-   memory usage of local atom-based arrays
-------------------------------------------------------------------------- */
-
-double PairComb::memory_usage()
-{
-  double bytes = maxeatom * sizeof(double);
-  bytes += maxvatom*6 * sizeof(double);
-  bytes += nmax * sizeof(int);
-  bytes += MAXNEIGH * nmax * sizeof(int);
-  return bytes;
-}
 /* ---------------------------------------------------------------------- */
 
 void PairComb::Short_neigh()
@@ -2064,21 +2070,17 @@ void PairComb::Short_neigh()
   ilist = list->ilist;
   numneigh = list->numneigh;
   firstneigh = list->firstneigh;
-  int npntj = 0;
-  int npage = 0;
+
+  // create Comb neighbor list
+
+  ipage->reset();
 
   for (ii = 0; ii < inum; ii++) {
     i = ilist[ii];
     itype = type[i];
 
-    if (pgsize - npntj < oneatom) {
-      npntj = 0;
-      npage++;
-      if (npage == maxpage) add_pages();
-    }
-
     nj = 0;
-    neighptrj = &pages[npage][npntj];
+    neighptrj = ipage->vget();
 
     xtmp = x[i][0];
     ytmp = x[i][1];
@@ -2103,19 +2105,27 @@ void PairComb::Short_neigh()
 
     sht_first[i] = neighptrj;
     sht_num[i] = nj;
-    npntj += nj;
+    ipage->vgot(nj);
+    if (ipage->status())
+      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 }
 
-/* ---------------------------------------------------------------------- */
+/* ----------------------------------------------------------------------
+   memory usage of local atom-based arrays
+------------------------------------------------------------------------- */
 
-void PairComb::add_pages(int howmany)
+double PairComb::memory_usage()
 {
-  int toppage = maxpage;
-  maxpage += howmany*PGDELTA;
+  double bytes = maxeatom * sizeof(double);
+  bytes += maxvatom*6 * sizeof(double);
+  bytes += nmax * sizeof(int);
+  bytes += nmax * sizeof(int *);
+
+  for (int i = 0; i < comm->nthreads; i++)
+    bytes += ipage[i].size();
 
-  pages = (int **)
-    memory->srealloc(pages,maxpage*sizeof(int *),"pair:pages");
-  for (int i = toppage; i < maxpage; i++)
-    memory->create(pages[i],pgsize,"pair:pages[i]");
+  bytes += nmax * sizeof(int);
+  bytes += MAXNEIGH*nmax * sizeof(double);
+  return bytes;
 }
diff --git a/src/MANYBODY/pair_comb.h b/src/MANYBODY/pair_comb.h
index 3131c95b1b82dd9081e443a2c194efc6ae25a9f6..cd3947a3d7e1d07944af75e0182570b933ea5d89 100644
--- a/src/MANYBODY/pair_comb.h
+++ b/src/MANYBODY/pair_comb.h
@@ -21,6 +21,7 @@ PairStyle(comb,PairComb)
 #define LMP_PAIR_COMB_H
 
 #include "pair.h"
+#include "my_page.h"
 
 namespace LAMMPS_NS {
 
@@ -77,6 +78,12 @@ class PairComb : public Pair {
   int *NCo, cor_flag, cuo_flag, cuo_flag1, cuo_flag2;
   double **bbij;
 
+  int pgsize;                      // size of neighbor page
+  int oneatom;                     // max # of neighbors for one atom
+  int *sht_num,**sht_first;        // short-range neighbor list
+  MyPage<int> *ipage;              // neighbor list pages
+  double cutmin;
+
   void allocate();
   virtual void read_file(char *);
   void setup();
@@ -145,13 +152,7 @@ class PairComb : public Pair {
   int pack_comm(int , int *, double *, int, int *);
   void unpack_comm(int , int , double *);
 
-  // short range neighbor list
-
-  void add_pages(int howmany = 1);
   void Short_neigh();
-  int maxpage, pgsize, oneatom, **pages;
-  int *sht_num, **sht_first;        // short-range neighbor list
-  double cutmin;
 
   // vector functions, inline for efficiency
 
diff --git a/src/MANYBODY/pair_lcbop.cpp b/src/MANYBODY/pair_lcbop.cpp
index 5e7c6974c20cf66e06c60064107d8e2e235671cf..cc511d59a40064ee8fe1a164534b2a94717c4816 100644
--- a/src/MANYBODY/pair_lcbop.cpp
+++ b/src/MANYBODY/pair_lcbop.cpp
@@ -30,6 +30,7 @@
 #include "neighbor.h"
 #include "neigh_list.h"
 #include "neigh_request.h"
+#include "my_page.h"
 #include "math_const.h"
 #include "memory.h"
 #include "error.h"
@@ -52,8 +53,9 @@ PairLCBOP::PairLCBOP(LAMMPS *lmp) : Pair(lmp) {
   maxlocal = 0;
   SR_numneigh = NULL;
   SR_firstneigh = NULL;
-  maxpage = 0;
-  pages = NULL;
+  ipage = NULL;
+  pgsize = oneatom = 0;
+
   N = NULL;
   M = NULL;
 }
@@ -62,11 +64,11 @@ PairLCBOP::PairLCBOP(LAMMPS *lmp) : Pair(lmp) {
    check if allocated, since class can be destructed when incomplete
 ------------------------------------------------------------------------- */
 
-PairLCBOP::~PairLCBOP() {
+PairLCBOP::~PairLCBOP()
+{
   memory->destroy(SR_numneigh);
   memory->sfree(SR_firstneigh);
-  for (int i = 0; i < maxpage; i++) memory->destroy(pages[i]);
-  memory->sfree(pages);
+  delete [] ipage;
   memory->destroy(N);
   memory->destroy(M);
 
@@ -191,11 +193,24 @@ void PairLCBOP::init_style()
   neighbor->requests[irequest]->full = 1;
   neighbor->requests[irequest]->ghost = 1;
 
-  // local SR neighbor list memory
+  // local SR neighbor list
+  // create pages if first time or if neighbor pgsize/oneatom has changed
+
+  int create = 0;
+  if (ipage == NULL) create = 1;
+  if (pgsize != neighbor->pgsize) create = 1;
+  if (oneatom != neighbor->oneatom) create = 1;
 
-  pgsize = neighbor->pgsize;
-  oneatom = neighbor->oneatom;
-  if (maxpage == 0) add_pages();
+  if (create) {
+    delete [] ipage;
+    pgsize = neighbor->pgsize;
+    oneatom = neighbor->oneatom;
+
+    int nmypage = comm->nthreads;
+    ipage = new MyPage<int>[nmypage];
+    for (int i = 0; i < nmypage; i++)
+      ipage[i].init(oneatom,pgsize,PGDELTA);
+  }
 }
 
 /* ----------------------------------------------------------------------
@@ -238,9 +253,7 @@ double PairLCBOP::init_one(int i, int j)
 
 void PairLCBOP::SR_neigh()
 {
-  int i,j,ii,jj,n,
-    allnum,           // number of atoms(both local and ghost) neighbors are stored for
-    jnum;
+  int i,j,ii,jj,n,allnum,jnum;
   double xtmp,ytmp,ztmp,delx,dely,delz,rsq,dS;
   int *ilist,*jlist,*numneigh,**firstneigh;
   int *neighptr;
@@ -270,19 +283,13 @@ void PairLCBOP::SR_neigh()
   // store all SR neighs of owned and ghost atoms
   // scan full neighbor list of I
 
-  int npage = 0;
-  int npnt = 0; // position in current page
+  ipage->reset();
 
   for (ii = 0; ii < allnum; ii++) {
     i = ilist[ii];
 
-    if (pgsize - npnt < oneatom) { // ensure at least oneatom space free at current page
-      npnt = 0;
-      npage++;
-      if (npage == maxpage) add_pages();
-    }
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     xtmp = x[i][0];
     ytmp = x[i][1];
@@ -307,14 +314,13 @@ void PairLCBOP::SR_neigh()
 
     SR_firstneigh[i] = neighptr;
     SR_numneigh[i] = n;
-    npnt += n;
-    if( npnt >= pgsize )
-      error->one(FLERR,
-                 "Neighbor list overflow, boost neigh_modify one or page");
+    ipage->vgot(n);
+    if (ipage->status())
+      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
-
   // calculate M_i
+
   for (ii = 0; ii < allnum; ii++) {
     i = ilist[ii];
 
@@ -861,21 +867,6 @@ double PairLCBOP::b(int i, int j, double rij[3],
   return bij;
 }
 
-/* ----------------------------------------------------------------------
-   add pages to SR neighbor list
-------------------------------------------------------------------------- */
-
-void PairLCBOP::add_pages(int howmany)
-{
-  int toppage = maxpage;
-  maxpage += howmany*PGDELTA;
-
-  pages = (int **)
-    memory->srealloc(pages,maxpage*sizeof(int *),"LCBOP:pages");
-  for (int i = toppage; i < maxpage; i++)
-    memory->create(pages[i],pgsize,"LCBOP:pages[i]");
-}
-
 /* ----------------------------------------------------------------------
    spline interpolation for G
 ------------------------------------------------------------------------- */
@@ -1287,11 +1278,15 @@ void PairLCBOP::spline_init() {
    memory usage of local atom-based arrays
 ------------------------------------------------------------------------- */
 
-double PairLCBOP::memory_usage() {
+double PairLCBOP::memory_usage()
+{
   double bytes = 0.0;
   bytes += maxlocal * sizeof(int);
   bytes += maxlocal * sizeof(int *);
-  bytes += maxpage * neighbor->pgsize * sizeof(int);
-  bytes += 3 * maxlocal * sizeof(double);
+
+  for (int i = 0; i < comm->nthreads; i++)
+    bytes += ipage[i].size();
+
+  bytes += 3*maxlocal * sizeof(double);
   return bytes;
 }
diff --git a/src/MANYBODY/pair_lcbop.h b/src/MANYBODY/pair_lcbop.h
index 71567cec95da7fa668391ac1e4e865706d408587..dc0ad61c7c99a5f485b99252eb508e6a5449256e 100644
--- a/src/MANYBODY/pair_lcbop.h
+++ b/src/MANYBODY/pair_lcbop.h
@@ -21,6 +21,7 @@ PairStyle(lcbop,PairLCBOP)
 #define LMP_PAIR_LCBOP_H
 
 #include "pair.h"
+#include "my_page.h"
 #include "math.h"
 #include "math_const.h"
 
@@ -42,18 +43,20 @@ class PairLCBOP : public Pair {
   int *map;                        // 0 (C) or -1 (NULL) for each type
 
   int me;
-  int maxlocal;                    // size of numneigh, firstneigh arrays
-  int maxpage;                     // # of pages currently allocated
-  int pgsize;                      // size of neighbor page
-  int oneatom;                     // max # of neighbors for one atom
 
   double cutLR;                    // LR cutoff
 
   double cutLRsq;                  // LR cutoff squared
   double cut3rebo;                 // maximum distance for 3rd SR neigh
 
+  int maxlocal;                    // size of numneigh, firstneigh arrays
+  int maxpage;                     // # of pages currently allocated
+  int pgsize;                      // size of neighbor page
+  int oneatom;                     // max # of neighbors for one atom
+  MyPage<int> *ipage;              // neighbor list pages
   int *SR_numneigh;                // # of pair neighbors for each atom
   int **SR_firstneigh;             // ptr to 1st neighbor of each atom
+
   double *N;                       // sum of cutoff fns ( f_C ) with SR neighs
   double *M;                       // sum_j f_C_ij*F(N_j - f_C_ij)
 
@@ -80,6 +83,7 @@ class PairLCBOP : public Pair {
         f_y_10,
         f_y_11;
   } F_conj_field[3][3][2];
+
   double F_conj_data[4][4][2][3]; // temporary data from file
   double gX[6];        // x coordinates for described points[# of points];
   double gC[5+1][6-1]; // coefficients for each period between described points [degree of polynomial+1][# of points-1]
@@ -98,7 +102,6 @@ class PairLCBOP : public Pair {
   void g_decompose_x( double, size_t*, double* );
   double F_conj( double, double, double, double*, double*, double* );
 
-  void add_pages(int howmany = 1);
   void read_file( char * );
 
   void spline_init();
diff --git a/src/USER-OMP/fix_shear_history_omp.cpp b/src/USER-OMP/fix_shear_history_omp.cpp
index d410f6363c57b36da87224ae28baa411aab27b0f..1fe89de3e4268505120d56a71f5b386b7720b751 100644
--- a/src/USER-OMP/fix_shear_history_omp.cpp
+++ b/src/USER-OMP/fix_shear_history_omp.cpp
@@ -31,8 +31,6 @@
 using namespace LAMMPS_NS;
 using namespace FixConst;
 
-#define MAXTOUCH 15
-
 /* ----------------------------------------------------------------------
    copy shear partner info from neighbor lists to atom arrays
    so can be exchanged with atoms
@@ -45,10 +43,10 @@ void FixShearHistoryOMP::pre_exchange()
   const int nghost = atom->nghost;
   const int nall = nlocal + nghost;
   const int nthreads = comm->nthreads;
-
-  int flag = 0;
+  maxtouch = 0;
+  
 #if defined(_OPENMP)
-#pragma omp parallel default(none) shared(flag)
+#pragma omp parallel default(none)
 #endif
   {
 
@@ -64,22 +62,23 @@ void FixShearHistoryOMP::pre_exchange()
     const int lmax = lfrom +ldelta;
     const int lto = (lmax > nlocal) ? nlocal : lmax;
 
-    const int gdelta = 1 + nghost/nthreads;
-    const int gfrom = nlocal + tid*gdelta;
-    const int gmax = gfrom + gdelta;
-    const int gto = (gmax > nall) ? nall : gmax;
-
-
-    int i,j,ii,jj,m,inum,jnum;
+    int i,j,ii,jj,m,n,inum,jnum;
     int *ilist,*jlist,*numneigh,**firstneigh;
     int *touch,**firsttouch;
     double *shear,*allshear,**firstshear;
 
-    // zero npartners for all current atoms
+    // zero npartners for all current atoms and
+    // clear page data structures for this thread
 
     for (i = lfrom; i < lto; i++) npartner[i] = 0;
 
-    // copy shear info from neighbor list atoms to atom arrays
+    MyPage <int> &ipg = ipage[tid];
+    MyPage <double[3]> &dpg = dpage[tid];
+    ipg.reset();
+    dpg.reset();
+
+    // 1st loop over neighbor list
+    // calculate nparter for each owned atom
 
     int *tag = atom->tag;
 
@@ -94,59 +93,81 @@ void FixShearHistoryOMP::pre_exchange()
     for (ii = 0; ii < inum; ii++) {
       i = ilist[ii];
       jlist = firstneigh[i];
-      allshear = firstshear[i];
       jnum = numneigh[i];
       touch = firsttouch[i];
 
       for (jj = 0; jj < jnum; jj++) {
         if (touch[jj]) {
+          if ((i >= lfrom) && (i < lto))
+            npartner[i]++;
+
           j = jlist[jj];
           j &= NEIGHMASK;
+          if ((j >= lfrom) && (j < lto))
+            npartner[j]++;
+        }
+      }
+    }
+
+    // get page chunks to store atom IDs and shear history for my atoms
+
+    for (ii = lfrom; ii < lto; ii++) {
+      i = ilist[ii];
+      n = npartner[i];
+      partner[i] = ipg.get(n);
+      shearpartner[i] = dpg.get(n);
+      if (partner[i] == NULL || shearpartner[i] == NULL)
+        error->one(FLERR,"Shear history overflow, boost neigh_modify one");
+    }
+
+    // 2nd loop over neighbor list
+    // store atom IDs and shear history for my atoms
+    // re-zero npartner to use as counter for all my atoms
+
+    for (i = lfrom; i < lto; i++) npartner[i] = 0;
+
+    for (ii = 0; ii < inum; ii++) {
+      i = ilist[ii];
+      jlist = firstneigh[i];
+      allshear = firstshear[i];
+      jnum = numneigh[i];
+      touch = firsttouch[i];
+
+      for (jj = 0; jj < jnum; jj++) {
+        if (touch[jj]) {
           shear = &allshear[3*jj];
+          j = jlist[jj];
+          j &= NEIGHMASK;
 
           if ((i >= lfrom) && (i < lto)) {
-            if (npartner[i] < MAXTOUCH) {
-              m = npartner[i];
-              partner[i][m] = tag[j];
-              shearpartner[i][m][0] = shear[0];
-              shearpartner[i][m][1] = shear[1];
-              shearpartner[i][m][2] = shear[2];
-            }
+            m = npartner[i];
+            partner[i][m] = tag[j];
+            shearpartner[i][m][0] = shear[0];
+            shearpartner[i][m][1] = shear[1];
+            shearpartner[i][m][2] = shear[2];
             npartner[i]++;
           }
 
           if ((j >= lfrom) && (j < lto)) {
-            if (npartner[j] < MAXTOUCH) {
-              m = npartner[j];
-              partner[j][m] = tag[i];
-              shearpartner[j][m][0] = -shear[0];
-              shearpartner[j][m][1] = -shear[1];
-              shearpartner[j][m][2] = -shear[2];
-            }
-            npartner[j]++;
-          }
-
-          if ((j >= gfrom) && (j < gto)) {
+            m = npartner[j];
+            partner[j][m] = tag[i];
+            shearpartner[j][m][0] = -shear[0];
+            shearpartner[j][m][1] = -shear[1];
+            shearpartner[j][m][2] = -shear[2];
             npartner[j]++;
           }
         }
       }
     }
 
-    // test for too many touching neighbors
-    int myflag = 0;
+    // set maxtouch = max # of partners of any owned atom
+    m = 0;
     for (i = lfrom; i < lto; i++)
-      if (npartner[i] >= MAXTOUCH) myflag = 1;
+      m = MAX(m,npartner[i]);
 
-    if (myflag)
 #if defined(_OPENMP)
-#pragma omp atomic
+#pragma omp critical
 #endif
-      ++flag;
+    maxtouch = MAX(m,maxtouch);
   }
-
-  int flag_all;
-  MPI_Allreduce(&flag,&flag_all,1,MPI_INT,MPI_SUM,world);
-  if (flag_all)
-    error->all(FLERR,"Too many touching neighbors - boost MAXTOUCH");
 }
diff --git a/src/USER-OMP/neigh_derive_omp.cpp b/src/USER-OMP/neigh_derive_omp.cpp
index 4356f17b62388ae1bc60500e957ea4e6c4803491..8f5fa05e81457cd78fe4c7a5c17744479f161df6 100644
--- a/src/USER-OMP/neigh_derive_omp.cpp
+++ b/src/USER-OMP/neigh_derive_omp.cpp
@@ -16,6 +16,7 @@
 #include "neigh_list.h"
 #include "atom.h"
 #include "comm.h"
+#include "my_page.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
@@ -48,27 +49,16 @@ void Neighbor::half_from_full_no_newton_omp(NeighList *list)
   int *numneigh_full = list->listfull->numneigh;
   int **firstneigh_full = list->listfull->firstneigh;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   // loop over atoms in full list
 
   for (ii = ifrom; ii < ito; ii++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      // only one thread at a time may check whether we
-      // need new neighbor list pages and then add to them.
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     // loop over parent full list
 
@@ -85,8 +75,8 @@ void Neighbor::half_from_full_no_newton_omp(NeighList *list)
     ilist[ii] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -124,27 +114,16 @@ void Neighbor::half_from_full_newton_omp(NeighList *list)
   int *numneigh_full = list->listfull->numneigh;
   int **firstneigh_full = list->listfull->firstneigh;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   // loop over parent full list
 
   for (ii = ifrom; ii < ito; ii++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      // only one thread at a time may check  whether we
-      // need new neighbor list pages and then add to them.
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     i = ilist_full[ii];
     xtmp = x[i][0];
@@ -174,8 +153,8 @@ void Neighbor::half_from_full_newton_omp(NeighList *list)
     ilist[ii] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
diff --git a/src/USER-OMP/neigh_full_omp.cpp b/src/USER-OMP/neigh_full_omp.cpp
index 8ba2ffaf2b5533abb6d95007c4f849b902186159..e609627897c6f993643e6c6703614521b42997e6 100644
--- a/src/USER-OMP/neigh_full_omp.cpp
+++ b/src/USER-OMP/neigh_full_omp.cpp
@@ -18,6 +18,7 @@
 #include "comm.h"
 #include "domain.h"
 #include "group.h"
+#include "my_page.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
@@ -57,24 +58,16 @@ void Neighbor::full_nsq_omp(NeighList *list)
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
 
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   // loop over owned atoms, storing neighbors
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -108,8 +101,8 @@ void Neighbor::full_nsq_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -152,24 +145,16 @@ void Neighbor::full_nsq_ghost_omp(NeighList *list)
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
 
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   // loop over owned & ghost atoms, storing neighbors
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -218,8 +203,8 @@ void Neighbor::full_nsq_ghost_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -266,24 +251,16 @@ void Neighbor::full_bin_omp(NeighList *list)
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
 
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   // loop over owned atoms, storing neighbors
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -322,8 +299,8 @@ void Neighbor::full_bin_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -374,24 +351,16 @@ void Neighbor::full_bin_ghost_omp(NeighList *list)
   int *stencil = list->stencil;
   int **stencilxyz = list->stencilxyz;
 
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   // loop over owned & ghost atoms, storing neighbors
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -457,8 +426,8 @@ void Neighbor::full_bin_ghost_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -510,22 +479,14 @@ void Neighbor::full_multi_omp(NeighList *list)
   int **stencil_multi = list->stencil_multi;
   double **distsq_multi = list->distsq_multi;
 
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -569,8 +530,8 @@ void Neighbor::full_multi_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
diff --git a/src/USER-OMP/neigh_gran_omp.cpp b/src/USER-OMP/neigh_gran_omp.cpp
index e8634dcf3acd49f89732d299fdcac5c1ae3e4088..5f0402f36397e230670e84b409b9a6b71c1492a5 100644
--- a/src/USER-OMP/neigh_gran_omp.cpp
+++ b/src/USER-OMP/neigh_gran_omp.cpp
@@ -40,10 +40,6 @@ void Neighbor::granular_nsq_no_newton_omp(NeighList *list)
 
   NEIGH_OMP_INIT;
 
-  if (fix_history)
-    if (nthreads > listgranhistory->maxpage)
-      listgranhistory->add_pages(nthreads - listgranhistory->maxpage);
-
 #if defined(_OPENMP)
 #pragma omp parallel default(none) shared(list,listgranhistory)
 #endif
@@ -59,6 +55,8 @@ void Neighbor::granular_nsq_no_newton_omp(NeighList *list)
   double (**shearpartner)[3];
   int **firsttouch;
   double **firstshear;
+  MyPage<int> *ipage_touch;
+  MyPage<double> *dpage_shear;
 
   double **x = atom->x;
   double *radius = atom->radius;
@@ -72,41 +70,32 @@ void Neighbor::granular_nsq_no_newton_omp(NeighList *list)
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
 
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
+
   if (fix_history) {
     npartner = fix_history->npartner;
     partner = fix_history->partner;
     shearpartner = fix_history->shearpartner;
     firsttouch = listgranhistory->firstneigh;
     firstshear = listgranhistory->firstdouble;
+    ipage_touch = listgranhistory->ipage+tid;
+    dpage_shear = listgranhistory->dpage+tid;
+    ipage_touch->reset();
+    dpage_shear->reset();
   }
 
-  int npage = tid;
-  int npnt = 0;
-
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    {
-      if (pgsize - npnt < oneatom) {
-        npnt = 0;
-        npage += nthreads;
-        if (npage >= list->maxpage) {
-          list->add_pages(nthreads);
-          if (fix_history)
-            listgranhistory->add_pages(nthreads);
-        }
-      }
-
-      n = nn = 0;
-      neighptr = &(list->pages[npage][npnt]);
-      if (fix_history) {
-        touchptr = &(listgranhistory->pages[npage][npnt]);
-        shearptr = &(listgranhistory->dpages[npage][3*npnt]);
-      }
+    n = 0;
+    neighptr = ipage.vget();
+    if (fix_history) {
+      nn = 0;
+      touchptr = ipage_touch->vget();
+      shearptr = dpage_shear->vget();
     }
-
+    
     xtmp = x[i][0];
     ytmp = x[i][1];
     ztmp = x[i][2];
@@ -158,13 +147,16 @@ void Neighbor::granular_nsq_no_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
+    ipage.vgot(n);
+    if (ipage.status())
+      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
+
     if (fix_history) {
       firsttouch[i] = touchptr;
       firstshear[i] = shearptr;
+      ipage_touch->vgot(n);
+      dpage_shear->vgot(nn);
     }
-    npnt += n;
-    if (n > oneatom)
-      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
   list->inum = nlocal;
@@ -207,22 +199,14 @@ void Neighbor::granular_nsq_newton_omp(NeighList *list)
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
 
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itag = tag[i];
     xtmp = x[i][0];
@@ -265,8 +249,8 @@ void Neighbor::granular_nsq_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -295,10 +279,6 @@ void Neighbor::granular_bin_no_newton_omp(NeighList *list)
 
   NEIGH_OMP_INIT;
 
-  if (fix_history)
-    if (nthreads > listgranhistory->maxpage)
-      listgranhistory->add_pages(nthreads - listgranhistory->maxpage);
-
 #if defined(_OPENMP)
 #pragma omp parallel default(none) shared(list,listgranhistory)
 #endif
@@ -309,6 +289,8 @@ void Neighbor::granular_bin_no_newton_omp(NeighList *list)
   double radi,radsum,cutsq;
   int *neighptr,*touchptr;
   double *shearptr;
+  MyPage<int> *ipage_touch;
+  MyPage<double> *dpage_shear;
 
   int *npartner,**partner;
   double (**shearpartner)[3];
@@ -330,39 +312,30 @@ void Neighbor::granular_bin_no_newton_omp(NeighList *list)
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
 
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
+
   if (fix_history) {
     npartner = fix_history->npartner;
     partner = fix_history->partner;
     shearpartner = fix_history->shearpartner;
     firsttouch = listgranhistory->firstneigh;
     firstshear = listgranhistory->firstdouble;
+    ipage_touch = listgranhistory->ipage+tid;
+    dpage_shear = listgranhistory->dpage+tid;
+    ipage_touch->reset();
+    dpage_shear->reset();
   }
 
-  int npage = tid;
-  int npnt = 0;
-
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    {
-      if (pgsize - npnt < oneatom) {
-        npnt = 0;
-        npage += nthreads;
-        if (npage >= list->maxpage) {
-          list->add_pages(nthreads);
-          if (fix_history)
-            listgranhistory->add_pages(nthreads);
-        }
-      }
-
-      n = nn = 0;
-      neighptr = &(list->pages[npage][npnt]);
-      if (fix_history) {
-        touchptr = &(listgranhistory->pages[npage][npnt]);
-        shearptr = &(listgranhistory->dpages[npage][3*npnt]);
-      }
+    n = 0;
+    neighptr = ipage.vget();
+    if (fix_history) {
+      nn = 0;
+      touchptr = ipage_touch->vget();
+      shearptr = dpage_shear->vget();
     }
 
     xtmp = x[i][0];
@@ -422,13 +395,16 @@ void Neighbor::granular_bin_no_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
+    ipage.vgot(n);
+    if (ipage.status())
+      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
+
     if (fix_history) {
       firsttouch[i] = touchptr;
       firstshear[i] = shearptr;
+      ipage_touch->vgot(n);
+      dpage_shear->vgot(nn);
     }
-    npnt += n;
-    if (n > oneatom)
-      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
   list->inum = nlocal;
@@ -475,22 +451,14 @@ void Neighbor::granular_bin_newton_omp(NeighList *list)
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
 
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
     n = 0;
-    neighptr = &(list->pages[npage][npnt]);
+    neighptr = ipage.vget();
 
     xtmp = x[i][0];
     ytmp = x[i][1];
@@ -543,8 +511,8 @@ void Neighbor::granular_bin_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -592,22 +560,14 @@ void Neighbor::granular_bin_newton_tri_omp(NeighList *list)
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
 
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
     n = 0;
-    neighptr = &(list->pages[npage][npnt]);
+    neighptr = ipage.vget();
 
     xtmp = x[i][0];
     ytmp = x[i][1];
@@ -648,8 +608,8 @@ void Neighbor::granular_bin_newton_tri_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
diff --git a/src/USER-OMP/neigh_half_bin_omp.cpp b/src/USER-OMP/neigh_half_bin_omp.cpp
index ba0a180410361f820ad5892cbe8057002464f3e1..f501684a19642d7a79782c644da6426f02f0e2a0 100644
--- a/src/USER-OMP/neigh_half_bin_omp.cpp
+++ b/src/USER-OMP/neigh_half_bin_omp.cpp
@@ -17,6 +17,7 @@
 #include "atom.h"
 #include "comm.h"
 #include "domain.h"
+#include "my_page.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
@@ -64,23 +65,14 @@ void Neighbor::half_bin_no_newton_omp(NeighList *list)
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -121,8 +113,8 @@ void Neighbor::half_bin_no_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -177,23 +169,14 @@ void Neighbor::half_bin_no_newton_ghost_omp(NeighList *list)
   int *stencil = list->stencil;
   int **stencilxyz = list->stencilxyz;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -263,8 +246,8 @@ void Neighbor::half_bin_no_newton_ghost_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -314,23 +297,14 @@ void Neighbor::half_bin_newton_omp(NeighList *list)
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -399,8 +373,8 @@ void Neighbor::half_bin_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -449,23 +423,14 @@ void Neighbor::half_bin_newton_tri_omp(NeighList *list)
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -513,8 +478,8 @@ void Neighbor::half_bin_newton_tri_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
diff --git a/src/USER-OMP/neigh_half_multi_omp.cpp b/src/USER-OMP/neigh_half_multi_omp.cpp
index 5cc6ba140128ea8266afdabbeaf1fa1ef8027a11..b46594870fa38100b833d3d457ee4e605e030abe 100644
--- a/src/USER-OMP/neigh_half_multi_omp.cpp
+++ b/src/USER-OMP/neigh_half_multi_omp.cpp
@@ -17,6 +17,7 @@
 #include "atom.h"
 #include "comm.h"
 #include "domain.h"
+#include "my_page.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
@@ -67,23 +68,14 @@ void Neighbor::half_multi_no_newton_omp(NeighList *list)
   int **stencil_multi = list->stencil_multi;
   double **distsq_multi = list->distsq_multi;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -129,8 +121,8 @@ void Neighbor::half_multi_no_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -182,23 +174,14 @@ void Neighbor::half_multi_newton_omp(NeighList *list)
   int **stencil_multi = list->stencil_multi;
   double **distsq_multi = list->distsq_multi;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -272,8 +255,8 @@ void Neighbor::half_multi_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -325,23 +308,14 @@ void Neighbor::half_multi_newton_tri_omp(NeighList *list)
   int **stencil_multi = list->stencil_multi;
   double **distsq_multi = list->distsq_multi;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -396,8 +370,8 @@ void Neighbor::half_multi_newton_tri_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
diff --git a/src/USER-OMP/neigh_half_nsq_omp.cpp b/src/USER-OMP/neigh_half_nsq_omp.cpp
index 6a8bee491a2b85859fbc40daca03fbb133280153..a98c4721f3d5dfe2a91be852874e191a618d5ac2 100644
--- a/src/USER-OMP/neigh_half_nsq_omp.cpp
+++ b/src/USER-OMP/neigh_half_nsq_omp.cpp
@@ -18,6 +18,7 @@
 #include "comm.h"
 #include "domain.h"
 #include "group.h"
+#include "my_page.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
@@ -58,25 +59,16 @@ void Neighbor::half_nsq_no_newton_omp(NeighList *list)
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   // loop over owned atoms, storing neighbors
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -110,8 +102,8 @@ void Neighbor::half_nsq_no_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -156,25 +148,16 @@ void Neighbor::half_nsq_no_newton_ghost_omp(NeighList *list)
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
 
-  // each thread works on its own page
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   // loop over owned & ghost atoms, storing neighbors
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -228,8 +211,8 @@ void Neighbor::half_nsq_no_newton_ghost_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
@@ -275,22 +258,14 @@ void Neighbor::half_nsq_newton_omp(NeighList *list)
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
 
-  int npage = tid;
-  int npnt = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  ipage.reset();
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-
-    neighptr = &(list->pages[npage][npnt]);
     n = 0;
+    neighptr = ipage.vget();
 
     itag = tag[i];
     itype = type[i];
@@ -341,8 +316,8 @@ void Neighbor::half_nsq_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
   NEIGH_OMP_CLOSE;
diff --git a/src/USER-OMP/neigh_respa_omp.cpp b/src/USER-OMP/neigh_respa_omp.cpp
index 28afc81f1154cd15cc55e8c41a4b7edc4a8e4e19..11c7fd9ae307a5c00ccd853461fdfbf8ec349213 100644
--- a/src/USER-OMP/neigh_respa_omp.cpp
+++ b/src/USER-OMP/neigh_respa_omp.cpp
@@ -18,6 +18,7 @@
 #include "comm.h"
 #include "domain.h"
 #include "group.h"
+#include "my_page.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
@@ -37,16 +38,8 @@ void Neighbor::respa_nsq_no_newton_omp(NeighList *list)
   NEIGH_OMP_INIT;
 
   NeighList *listinner = list->listinner;
-  if (nthreads > listinner->maxpage)
-    listinner->add_pages(nthreads - listinner->maxpage);
-
-  NeighList *listmiddle;
+  NeighList *listmiddle = list->listmiddle;
   const int respamiddle = list->respamiddle;
-  if (respamiddle) {
-    listmiddle = list->listmiddle;
-    if (nthreads > listmiddle->maxpage)
-      listmiddle->add_pages(nthreads - listmiddle->maxpage);
-  }
 
 #if defined(_OPENMP)
 #pragma omp parallel default(none) shared(list,listinner,listmiddle)
@@ -85,48 +78,30 @@ void Neighbor::respa_nsq_no_newton_omp(NeighList *list)
     firstneigh_middle = listmiddle->firstneigh;
   }
 
-  int npage = tid;
-  int npnt = 0;
-  int npage_inner = tid;
-  int npnt_inner = 0;
-  int npage_middle = tid;
-  int npnt_middle = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  MyPage<int> &ipage_inner = listinner->ipage[tid];
+  ipage.reset();
+  ipage_inner.reset();
+
+  MyPage<int> *ipage_middle;
+  if (respamiddle) {
+    ipage_middle = listmiddle->ipage + tid;
+    ipage_middle->reset();
+  }
 
   int which = 0;
   int minchange = 0;
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-   {
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-    neighptr = &(list->pages[npage][npnt]);
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner += nthreads;
-      if (npage_inner >= listinner->maxpage) listinner->add_pages(nthreads);
-    }
-    neighptr_inner = &(listinner->pages[npage_inner][npnt_inner]);
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage.vget();
+    neighptr_inner = ipage_inner.vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle += nthreads;
-        if (npage_middle >= listmiddle->maxpage) listmiddle->add_pages(nthreads);
-      }
-      neighptr_middle = &(listmiddle->pages[npage_middle][npnt_middle]);
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
-   }
 
     itype = type[i];
     xtmp = x[i][0];
@@ -172,23 +147,23 @@ void Neighbor::respa_nsq_no_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     ilist_inner[i] = i;
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage.vgot(n_inner);
+    if (ipage_inner.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       ilist_middle[i] = i;
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n_middle);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
   }
@@ -214,16 +189,8 @@ void Neighbor::respa_nsq_newton_omp(NeighList *list)
   NEIGH_OMP_INIT;
 
   NeighList *listinner = list->listinner;
-  if (nthreads > listinner->maxpage)
-    listinner->add_pages(nthreads - listinner->maxpage);
-
-  NeighList *listmiddle;
+  NeighList *listmiddle = list->listmiddle;
   const int respamiddle = list->respamiddle;
-  if (respamiddle) {
-    listmiddle = list->listmiddle;
-    if (nthreads > listmiddle->maxpage)
-      listmiddle->add_pages(nthreads - listmiddle->maxpage);
-  }
 
 #if defined(_OPENMP)
 #pragma omp parallel default(none) shared(list,listinner,listmiddle)
@@ -262,48 +229,30 @@ void Neighbor::respa_nsq_newton_omp(NeighList *list)
     firstneigh_middle = listmiddle->firstneigh;
   }
 
-  int npage = tid;
-  int npnt = 0;
-  int npage_inner = tid;
-  int npnt_inner = 0;
-  int npage_middle = tid;
-  int npnt_middle = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  MyPage<int> &ipage_inner = listinner->ipage[tid];
+  ipage.reset();
+  ipage_inner.reset();
+
+  MyPage<int> *ipage_middle;
+  if (respamiddle) {
+    ipage_middle = listmiddle->ipage + tid;
+    ipage_middle->reset();
+  }
 
   int which = 0;
   int minchange = 0;
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-   {
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-    neighptr = &(list->pages[npage][npnt]);
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner += nthreads;
-      if (npage_inner >= listinner->maxpage) listinner->add_pages(nthreads);
-    }
-    neighptr_inner = &(listinner->pages[npage_inner][npnt_inner]);
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage.vget();
+    neighptr_inner = ipage_inner.vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle += nthreads;
-        if (npage_middle >= listmiddle->maxpage) listmiddle->add_pages(nthreads);
-      }
-      neighptr_middle = &(listmiddle->pages[npage_middle][npnt_middle]);
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
-   }
 
     itag = tag[i];
     itype = type[i];
@@ -367,23 +316,23 @@ void Neighbor::respa_nsq_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     ilist_inner[i] = i;
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage.vgot(n_inner);
+    if (ipage_inner.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       ilist_middle[i] = i;
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n_middle);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
   }
@@ -412,16 +361,8 @@ void Neighbor::respa_bin_no_newton_omp(NeighList *list)
   NEIGH_OMP_INIT;
 
   NeighList *listinner = list->listinner;
-  if (nthreads > listinner->maxpage)
-    listinner->add_pages(nthreads - listinner->maxpage);
-
-  NeighList *listmiddle;
+  NeighList *listmiddle = list->listmiddle;
   const int respamiddle = list->respamiddle;
-  if (respamiddle) {
-    listmiddle = list->listmiddle;
-    if (nthreads > listmiddle->maxpage)
-      listmiddle->add_pages(nthreads - listmiddle->maxpage);
-  }
 
 #if defined(_OPENMP)
 #pragma omp parallel default(none) shared(list,listinner,listmiddle)
@@ -461,48 +402,30 @@ void Neighbor::respa_bin_no_newton_omp(NeighList *list)
     firstneigh_middle = listmiddle->firstneigh;
   }
 
-  int npage = tid;
-  int npnt = 0;
-  int npage_inner = tid;
-  int npnt_inner = 0;
-  int npage_middle = tid;
-  int npnt_middle = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  MyPage<int> &ipage_inner = listinner->ipage[tid];
+  ipage.reset();
+  ipage_inner.reset();
+
+  MyPage<int> *ipage_middle;
+  if (respamiddle) {
+    ipage_middle = listmiddle->ipage + tid;
+    ipage_middle->reset();
+  }
 
   int which = 0;
   int minchange = 0;
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-   {
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-    neighptr = &(list->pages[npage][npnt]);
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner += nthreads;
-      if (npage_inner >= listinner->maxpage) listinner->add_pages(nthreads);
-    }
-    neighptr_inner = &(listinner->pages[npage_inner][npnt_inner]);
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage.vget();
+    neighptr_inner = ipage_inner.vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle += nthreads;
-        if (npage_middle >= listmiddle->maxpage) listmiddle->add_pages(nthreads);
-      }
-      neighptr_middle = &(listmiddle->pages[npage_middle][npnt_middle]);
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
-   }
 
     itype = type[i];
     xtmp = x[i][0];
@@ -557,23 +480,23 @@ void Neighbor::respa_bin_no_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     ilist_inner[i] = i;
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage.vgot(n_inner);
+    if (ipage_inner.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       ilist_middle[i] = i;
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n_middle);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
   }
@@ -601,16 +524,8 @@ void Neighbor::respa_bin_newton_omp(NeighList *list)
   NEIGH_OMP_INIT;
 
   NeighList *listinner = list->listinner;
-  if (nthreads > listinner->maxpage)
-    listinner->add_pages(nthreads - listinner->maxpage);
-
-  NeighList *listmiddle;
+  NeighList *listmiddle = list->listmiddle;
   const int respamiddle = list->respamiddle;
-  if (respamiddle) {
-    listmiddle = list->listmiddle;
-    if (nthreads > listmiddle->maxpage)
-      listmiddle->add_pages(nthreads - listmiddle->maxpage);
-  }
 
 #if defined(_OPENMP)
 #pragma omp parallel default(none) shared(list,listinner,listmiddle)
@@ -650,48 +565,30 @@ void Neighbor::respa_bin_newton_omp(NeighList *list)
     firstneigh_middle = listmiddle->firstneigh;
   }
 
-  int npage = tid;
-  int npnt = 0;
-  int npage_inner = tid;
-  int npnt_inner = 0;
-  int npage_middle = tid;
-  int npnt_middle = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  MyPage<int> &ipage_inner = listinner->ipage[tid];
+  ipage.reset();
+  ipage_inner.reset();
+
+  MyPage<int> *ipage_middle;
+  if (respamiddle) {
+    ipage_middle = listmiddle->ipage + tid;
+    ipage_middle->reset();
+  }
 
   int which = 0;
   int minchange = 0;
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-   {
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-    neighptr = &(list->pages[npage][npnt]);
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner += nthreads;
-      if (npage_inner >= listinner->maxpage) listinner->add_pages(nthreads);
-    }
-    neighptr_inner = &(listinner->pages[npage_inner][npnt_inner]);
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage.vget();
+    neighptr_inner = ipage_inner.vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle += nthreads;
-        if (npage_middle >= listmiddle->maxpage) listmiddle->add_pages(nthreads);
-      }
-      neighptr_middle = &(listmiddle->pages[npage_middle][npnt_middle]);
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
-   }
 
     itype = type[i];
     xtmp = x[i][0];
@@ -787,23 +684,23 @@ void Neighbor::respa_bin_newton_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     ilist_inner[i] = i;
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage.vgot(n_inner);
+    if (ipage_inner.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       ilist_middle[i] = i;
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n_middle);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
   }
@@ -831,16 +728,8 @@ void Neighbor::respa_bin_newton_tri_omp(NeighList *list)
   NEIGH_OMP_INIT;
 
   NeighList *listinner = list->listinner;
-  if (nthreads > listinner->maxpage)
-    listinner->add_pages(nthreads - listinner->maxpage);
-
-  NeighList *listmiddle;
+  NeighList *listmiddle = list->listmiddle;
   const int respamiddle = list->respamiddle;
-  if (respamiddle) {
-    listmiddle = list->listmiddle;
-    if (nthreads > listmiddle->maxpage)
-      listmiddle->add_pages(nthreads - listmiddle->maxpage);
-  }
 
 #if defined(_OPENMP)
 #pragma omp parallel default(none) shared(list,listinner,listmiddle)
@@ -880,48 +769,30 @@ void Neighbor::respa_bin_newton_tri_omp(NeighList *list)
     firstneigh_middle = listmiddle->firstneigh;
   }
 
-  int npage = tid;
-  int npnt = 0;
-  int npage_inner = tid;
-  int npnt_inner = 0;
-  int npage_middle = tid;
-  int npnt_middle = 0;
+  // each thread has its own page allocator
+  MyPage<int> &ipage = list->ipage[tid];
+  MyPage<int> &ipage_inner = listinner->ipage[tid];
+  ipage.reset();
+  ipage_inner.reset();
+
+  MyPage<int> *ipage_middle;
+  if (respamiddle) {
+    ipage_middle = listmiddle->ipage + tid;
+    ipage_middle->reset();
+  }
 
   int which = 0;
   int minchange = 0;
 
   for (i = ifrom; i < ito; i++) {
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-   {
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage += nthreads;
-      if (npage >= list->maxpage) list->add_pages(nthreads);
-    }
-    neighptr = &(list->pages[npage][npnt]);
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner += nthreads;
-      if (npage_inner >= listinner->maxpage) listinner->add_pages(nthreads);
-    }
-    neighptr_inner = &(listinner->pages[npage_inner][npnt_inner]);
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage.vget();
+    neighptr_inner = ipage_inner.vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle += nthreads;
-        if (npage_middle >= listmiddle->maxpage) listmiddle->add_pages(nthreads);
-      }
-      neighptr_middle = &(listmiddle->pages[npage_middle][npnt_middle]);
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
-   }
 
     itype = type[i];
     xtmp = x[i][0];
@@ -984,23 +855,23 @@ void Neighbor::respa_bin_newton_tri_omp(NeighList *list)
     ilist[i] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage.vgot(n);
+    if (ipage.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     ilist_inner[i] = i;
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage.vgot(n_inner);
+    if (ipage_inner.status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       ilist_middle[i] = i;
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n_middle);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
   }
diff --git a/src/USER-OMP/neighbor_omp.h b/src/USER-OMP/neighbor_omp.h
index 8db58c7e20be9ede45ca12b83d0ce4364eb9bd0b..2b2ad24fe5312806e78bef1ac0487a0fa49be76d 100644
--- a/src/USER-OMP/neighbor_omp.h
+++ b/src/USER-OMP/neighbor_omp.h
@@ -25,9 +25,7 @@ namespace LAMMPS_NS {
 
 // make sure we have at least one page for each thread
 #define NEIGH_OMP_INIT                          \
-  const int nthreads = comm->nthreads;          \
-  if (nthreads > list->maxpage)                 \
-    list->add_pages(nthreads - list->maxpage)
+  const int nthreads = comm->nthreads;
 
 // get thread id and then assign each thread a fixed chunk of atoms
 #define NEIGH_OMP_SETUP(num)                    \
diff --git a/src/USER-OMP/pair_airebo_omp.cpp b/src/USER-OMP/pair_airebo_omp.cpp
index 8ab26e9bf7edf82c4b86bb9f8cb45882cd0fbda9..ec788fd6ce8e0c7c993c767aeb3b5971f14708e2 100644
--- a/src/USER-OMP/pair_airebo_omp.cpp
+++ b/src/USER-OMP/pair_airebo_omp.cpp
@@ -19,6 +19,7 @@
 #include "error.h"
 #include "force.h"
 #include "memory.h"
+#include "my_page.h"
 #include "math_special.h"
 #include "neighbor.h"
 #include "neigh_list.h"
@@ -2671,9 +2672,6 @@ void PairAIREBOOMP::REBO_neigh_thr()
     memory->create(nH,maxlocal,"AIREBO:nH");
   }
 
-  if (nthreads > maxpage)
-    add_pages(nthreads - maxpage);
-
 #if defined(_OPENMP)
 #pragma omp parallel default(none)
 #endif
@@ -2704,22 +2702,15 @@ void PairAIREBOOMP::REBO_neigh_thr()
     // store all REBO neighs of owned and ghost atoms
     // scan full neighbor list of I
 
-    int npage = tid;
-    int npnt = 0;
+    // each thread has its own page allocator
+    MyPage<int> &ipg = ipage[tid];
+    ipg.reset();
 
     for (ii = iifrom; ii < iito; ii++) {
       i = ilist[ii];
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-      if (pgsize - npnt < oneatom) {
-        npnt = 0;
-        npage += nthreads;
-        if (npage >= maxpage) add_pages(nthreads);
-      }
-      neighptr = &(pages[npage][npnt]);
       n = 0;
+      neighptr = ipg.vget();
 
       xtmp = x[i][0];
       ytmp = x[i][1];
@@ -2749,10 +2740,9 @@ void PairAIREBOOMP::REBO_neigh_thr()
 
       REBO_firstneigh[i] = neighptr;
       REBO_numneigh[i] = n;
-      npnt += n;
-
-      if (npnt >= pgsize)
-        error->one(FLERR,"REBO list overflow, boost neigh_modify one or page");
+      ipg.vgot(n);
+      if (ipg.status())
+        error->one(FLERR,"REBO list overflow, boost neigh_modify one");
     }
   }
 }
diff --git a/src/USER-OMP/pair_comb_omp.cpp b/src/USER-OMP/pair_comb_omp.cpp
index eea2c6b61027d3a55aeccf3c75b2a31c3350ec07..1a8c32f4a79c70abbc0078170674dd00d5a7c57c 100644
--- a/src/USER-OMP/pair_comb_omp.cpp
+++ b/src/USER-OMP/pair_comb_omp.cpp
@@ -19,6 +19,7 @@
 #include "group.h"
 #include "force.h"
 #include "memory.h"
+#include "my_page.h"
 #include "neighbor.h"
 #include "neigh_list.h"
 
@@ -553,15 +554,13 @@ void PairCombOMP::Short_neigh_thr()
     nmax = atom->nmax;
     memory->sfree(sht_first);
     sht_first = (int **) memory->smalloc(nmax*sizeof(int *),
-            "pair:sht_first");
+                                         "pair:sht_first");
     memory->grow(sht_num,nmax,"pair:sht_num");
     memory->grow(NCo,nmax,"pair:NCo");
-    memory->grow(bbij,nmax,nmax,"pair:bbij");
+    memory->grow(bbij,nmax,MAXNEIGH,"pair:bbij");
   }
 
   const int nthreads = comm->nthreads;
-  if (nthreads > maxpage)
-    add_pages(nthreads - maxpage);
 
 #if defined(_OPENMP)
 #pragma omp parallel default(none)
@@ -588,51 +587,44 @@ void PairCombOMP::Short_neigh_thr()
 
     const int iidelta = 1 + inum/nthreads;
     const int iifrom = tid*iidelta;
-    int iito   = iifrom + iidelta;
-    if (iito > inum) iito = inum;
+    const int iito   = ((iifrom + iidelta) > inum) ? inum : (iifrom+iidelta);
 
-    int npage = tid;
-    npntj = 0;
+    // each thread has its own page allocator
+    MyPage<int> &ipg = ipage[tid];
+    ipg.reset();
 
-    if (iifrom < inum) {
-      for (ii = iifrom; ii < iito; ii++) {
-        i = ilist[ii];
+    // create Comb neighbor list
 
-#if defined(_OPENMP)
-#pragma omp critical
-#endif
-        if(pgsize - npntj < oneatom) {
-          npntj = 0;
-          npage += nthreads;
-          if (npage >= maxpage) add_pages(nthreads);
-        }
+    for (ii = iifrom; ii < iito; ii++) {
+      i = ilist[ii];
 
-        neighptrj = &pages[npage][npntj];
-        nj = 0;
+      nj = 0;
+      neighptrj = ipg.vget();
 
-        xtmp = x[i][0];
-        ytmp = x[i][1];
-        ztmp = x[i][2];
+      xtmp = x[i][0];
+      ytmp = x[i][1];
+      ztmp = x[i][2];
 
-        jlist = firstneigh[i];
-        jnum = numneigh[i];
+      jlist = firstneigh[i];
+      jnum = numneigh[i];
 
-        for (jj = 0; jj < jnum; jj++) {
-          j = jlist[jj];
-          j &= NEIGHMASK;
+      for (jj = 0; jj < jnum; jj++) {
+        j = jlist[jj];
+        j &= NEIGHMASK;
 
-          delrj[0] = xtmp - x[j][0];
-          delrj[1] = ytmp - x[j][1];
-          delrj[2] = ztmp - x[j][2];
-          rsq = vec3_dot(delrj,delrj);
+        delrj[0] = xtmp - x[j][0];
+        delrj[1] = ytmp - x[j][1];
+        delrj[2] = ztmp - x[j][2];
+        rsq = vec3_dot(delrj,delrj);
 
-          if (rsq > cutmin) continue;
-          neighptrj[nj++] = j;
-        }
-        sht_first[i] = neighptrj;
-        sht_num[i] = nj;
-        npntj += nj;
+        if (rsq > cutmin) continue;
+        neighptrj[nj++] = j;
       }
+      sht_first[i] = neighptrj;
+      sht_num[i] = nj;
+      ipg.vgot(nj);
+      if (ipg.status())
+        error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
   }
 }
diff --git a/src/fix_shear_history.cpp b/src/fix_shear_history.cpp
index da0def908b4b8067b0c653241db8dab8c23fa713..5da158339f711f51c30a0154ef14f3b84a9ce6c1 100644
--- a/src/fix_shear_history.cpp
+++ b/src/fix_shear_history.cpp
@@ -16,6 +16,7 @@
 #include "stdio.h"
 #include "fix_shear_history.h"
 #include "atom.h"
+#include "comm.h"
 #include "neighbor.h"
 #include "neigh_list.h"
 #include "force.h"
@@ -71,8 +72,8 @@ FixShearHistory::~FixShearHistory()
   memory->destroy(npartner);
   memory->sfree(partner);
   memory->sfree(shearpartner);
-  delete ipage;
-  delete dpage;
+  delete [] ipage;
+  delete [] dpage;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -105,10 +106,18 @@ void FixShearHistory::init()
   if (oneatom != neighbor->oneatom) create = 1;
 
   if (create) {
+    delete [] ipage;
+    delete [] dpage;
+
     pgsize = neighbor->pgsize;
     oneatom = neighbor->oneatom;
-    ipage = new MyPage<int>(oneatom,pgsize);
-    dpage = new MyPage<double[3]>(oneatom,pgsize);
+    int nmypage = comm->nthreads;
+    ipage = new MyPage<int>[nmypage];
+    dpage = new MyPage<double[3]>[nmypage];
+    for (int i = 0; i < nmypage; i++) {
+      ipage[i].init(oneatom,pgsize);
+      dpage[i].init(oneatom,pgsize);
+    }
   }
 }
 
@@ -260,8 +269,13 @@ double FixShearHistory::memory_usage()
   double bytes = nmax * sizeof(int);
   bytes = nmax * sizeof(int *);
   bytes = nmax * sizeof(double *);
-  bytes += ipage->ndatum * sizeof(int);
-  bytes += dpage->ndatum * sizeof(double[3]);
+
+  int nmypage = comm->nthreads;
+  for (int i = 0; i < nmypage; i++) {
+    bytes += ipage[i].size();
+    bytes += dpage[i].size();
+  }
+
   return bytes;
 }
 
diff --git a/src/my_page.h b/src/my_page.h
index 1628a4bab6c1804a8f607b83b2d4ef9db908c890..5fdf1e13516ab3ade4f761c99c08a331b1dc4054 100644
--- a/src/my_page.h
+++ b/src/my_page.h
@@ -1,4 +1,4 @@
-/* ----------------------------------------------------------------------
+/* -*- c++ -*- ----------------------------------------------------------
    LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
    http://lammps.sandia.gov, Sandia National Laboratories
    Steve Plimpton, sjplimp@sandia.gov
@@ -24,22 +24,25 @@ usage:
 inputs:
    template T = one datum, e.g. int, double, struct, int[3]
      for int[3], access datum as ivec[i][2]
-   maxchunk = max # of datums in one chunk, default = 1
-   pagesize = # of datums in one page, default = 1024
-              should be insure big enough to store multiple chunks
-   pagedelta = # of pages to allocate at a time, default = 1
 methods:
    T *get() = return ptr to one datum
    T *get(N) = return ptr to N datums, N < maxchunk required
    T *vget() = return ptr to maxchunk datums, use as needed, then call vgot()
      all gets return NULL if error encountered
    vgot(N) = used N datums of previous vget(), N < maxchunk required
+   void init(maxchunk, pagesize, pagedelta)
+     define allocation params and allocate first page(s)
+     call right after constructor
+       can call again to reset allocation params and free previous pages
+     maxchunk = max # of datums in one chunk, default = 1
+     pagesize = # of datums in one page, default = 1024
+       should be big enough to store multiple chunks
+     pagedelta = # of pages to allocate at a time, default = 1
+     return 1 if bad params
    void reset() = clear pages w/out freeing
    int size() = return total size of allocated pages in bytes
-public variables:
-   ndatum = total # of stored datums
-   nchunk = total # of stored chunks
-   errorflag = flag for various error conditions
+   int status() = return error status
+     0 = ok, 1 = chunksize > maxchunk, 2 = allocation error
 ------------------------------------------------------------------------- */
 
 #ifndef LAMMPS_MY_PAGE_H
@@ -50,24 +53,33 @@ namespace LAMMPS_NS {
 
 template<class T>
 class MyPage {
- public:
   int ndatum;      // total # of stored datums
   int nchunk;      // total # of stored chunks
-  int errorflag;   // flag > 1 if error has occurred
-                   // 1 = invalid inputs
-                   // 2 = memory allocation error
-                   // 3 = chunk size exceeded maxchunk
 
-  MyPage(int user_maxchunk = 1, int user_pagesize = 1024, 
-         int user_pagedelta = 1) {
+ public:
+  MyPage() {
+    ndatum = nchunk = 0;
+    pages = NULL;
+    npage = 0;
+    errorflag = 0;
+  }
+
+  // (re)initialize allocation params
+  // also allocate first page(s)
+
+  int init(int user_maxchunk = 1, int user_pagesize = 1024, 
+           int user_pagedelta = 1) {
     maxchunk = user_maxchunk;
     pagesize = user_pagesize;
     pagedelta = user_pagedelta;
 
-    errorflag = 0;
-    if (maxchunk <= 0 || pagesize <= 0 || pagedelta <= 0) errorflag = 1;
-    if (maxchunk > pagesize) errorflag = 1;
-    if (errorflag) return;
+    if (maxchunk <= 0 || pagesize <= 0 || pagedelta <= 0) return 1;
+    if (maxchunk > pagesize) return 1;
+
+    // free any previously allocated pages
+
+    for (int i = 0; i < npage; i++) free(pages[i]);
+    free(pages);
 
     // initial page allocation
 
@@ -75,12 +87,13 @@ class MyPage {
     pages = NULL;
     npage = 0;
     allocate();
-    if (errorflag) return;
+    if (errorflag) return 2;
     ipage = index = 0;
     page = pages[ipage];
+    return 0;
   }
 
-  // free all pages of allocated memory
+  // free all allocated pages
 
   ~MyPage() {
     for (int i = 0; i < npage; i++) free(pages[i]);
@@ -110,7 +123,7 @@ class MyPage {
 
   T *get(int n) {
     if (n > maxchunk) {
-      errorflag = 3;
+      errorflag = 1;
       return NULL;
     }
     ndatum += n;
@@ -151,14 +164,13 @@ class MyPage {
   // error if N > maxchunk
 
   void vgot(int n) {
-    if (n > maxchunk) errorflag = 3;
+    if (n > maxchunk) errorflag = 1;
     ndatum += n;
     nchunk++;
     index += n;
   }
 
-  // reset index to beginning of first page
-  // effectively clears all pages, without freeing any memory
+  // clear all pages, without freeing any memory
 
   void reset() {
     ndatum = nchunk = 0;
@@ -166,23 +178,33 @@ class MyPage {
     page = pages[ipage];
   }
 
-  // return total size of all allocated pages
+  // return total size of allocated pages
 
-  int size() {
+  int size() const {
     return npage*pagesize*sizeof(T);
   }
 
- private:
-  int maxchunk;   // max # of datums in one requested chunk
-  int pagesize;   // # of datums in one page, default = 1024
-  int pagedelta;  // # of pages to allocate at once, default = 1
+  // return error status
 
+  int status() const {
+    return errorflag;
+  }
+
+ private:
   T **pages;      // list of allocated pages
   T *page;        // ptr to current page
   int npage;      // # of allocated pages
   int ipage;      // index of current page
   int index;      // current index on current page
   
+  int maxchunk;   // max # of datums in one requested chunk
+  int pagesize;   // # of datums in one page, default = 1024
+  int pagedelta;  // # of pages to allocate at once, default = 1
+
+  int errorflag;  // flag > 0 if error has occurred
+                  // 1 = chunk size exceeded maxchunk
+                  // 2 = memory allocation error
+
   void allocate() {
     npage += pagedelta;
     pages = (T **) realloc(pages,npage*sizeof(T *));
@@ -190,9 +212,17 @@ class MyPage {
       errorflag = 2;
       return;
     }
+
+    void *ptr;
     for (int i = npage-pagedelta; i < npage; i++) {
+#if defined(LAMMPS_MEMALIGN)
+      if (posix_memalign(&ptr, LAMMPS_MEMALIGN, pagesize*sizeof(T)))
+        errorflag = 2;
+      pages[i] = (T *) ptr;
+#else
       pages[i] = (T *) malloc(pagesize*sizeof(T));
       if (!pages[i]) errorflag = 2;
+#endif
     }
   }
 };
diff --git a/src/neigh_derive.cpp b/src/neigh_derive.cpp
index 4bda3d71d08b49080e079c5944c128979a9aa05d..343dc1dc6e494874790ec6714f07b8b742717cfd 100644
--- a/src/neigh_derive.cpp
+++ b/src/neigh_derive.cpp
@@ -14,6 +14,7 @@
 #include "neighbor.h"
 #include "neigh_list.h"
 #include "atom.h"
+#include "my_page.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
@@ -33,28 +34,21 @@ void Neighbor::half_from_full_no_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
+
   int *ilist_full = list->listfull->ilist;
   int *numneigh_full = list->listfull->numneigh;
   int **firstneigh_full = list->listfull->firstneigh;
   int inum_full = list->listfull->inum;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   // loop over atoms in full list
 
   for (ii = 0; ii < inum_full; ii++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     // loop over parent full list
 
@@ -71,8 +65,8 @@ void Neighbor::half_from_full_no_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -98,28 +92,21 @@ void Neighbor::half_from_full_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
+
   int *ilist_full = list->listfull->ilist;
   int *numneigh_full = list->listfull->numneigh;
   int **firstneigh_full = list->listfull->firstneigh;
   int inum_full = list->listfull->inum;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   // loop over parent full list
 
   for (ii = 0; ii < inum_full; ii++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     i = ilist_full[ii];
     xtmp = x[i][0];
@@ -149,8 +136,8 @@ void Neighbor::half_from_full_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -175,7 +162,8 @@ void Neighbor::skip_from(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
+
   int *ilist_skip = list->listskip->ilist;
   int *numneigh_skip = list->listskip->numneigh;
   int **firstneigh_skip = list->listskip->firstneigh;
@@ -186,8 +174,7 @@ void Neighbor::skip_from(NeighList *list)
   int **ijskip = list->ijskip;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   // loop over atoms in other list
   // skip I atom entirely if iskip is set for type[I]
@@ -198,14 +185,8 @@ void Neighbor::skip_from(NeighList *list)
     itype = type[i];
     if (iskip[itype]) continue;
 
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     // loop over parent non-skip list
 
@@ -222,8 +203,8 @@ void Neighbor::skip_from(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -255,7 +236,8 @@ void Neighbor::skip_from_granular(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
+
   int *ilist_skip = list->listskip->ilist;
   int *numneigh_skip = list->listskip->numneigh;
   int **firstneigh_skip = list->listskip->firstneigh;
@@ -269,12 +251,13 @@ void Neighbor::skip_from_granular(NeighList *list)
   NeighList *listgranhistory = list->listgranhistory;
   int **firsttouch = listgranhistory->firstneigh;
   double **firstshear = listgranhistory->firstdouble;
-  int **pages_touch = listgranhistory->pages;
-  double **pages_shear = listgranhistory->dpages;
+  MyPage<int> *ipage_touch = listgranhistory->ipage;
+  MyPage<double> *dpage_shear = listgranhistory->dpage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
+  ipage_touch->reset();
+  dpage_shear->reset();
 
   // loop over atoms in other list
   // skip I atom entirely if iskip is set for type[I]
@@ -285,21 +268,10 @@ void Neighbor::skip_from_granular(NeighList *list)
     itype = type[i];
     if (iskip[itype]) continue;
 
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) {
-        pages = list->add_pages();
-        pages_touch = listgranhistory->add_pages();
-        pages_shear = listgranhistory->dpages;
-      }
-    }
-
-    n = 0;
-    neighptr = &pages[npage][npnt];
-    nn = 0;
-    touchptr = &pages_touch[npage][npnt];
-    shearptr = &pages_shear[npage][3*npnt];
+    n = nn = 0;
+    neighptr = ipage->vget();
+    touchptr = ipage_touch->vget();
+    shearptr = dpage_shear->vget();
 
     // loop over parent non-skip granular list and its history info
 
@@ -322,11 +294,14 @@ void Neighbor::skip_from_granular(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
+    ipage->vgot(n);
+    if (ipage->status())
+      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
+
     firsttouch[i] = touchptr;
     firstshear[i] = shearptr;
-    npnt += n;
-    if (n > oneatom)
-      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
+    ipage_touch->vgot(n);
+    dpage_shear->vgot(nn);
   }
 
   list->inum = inum;
@@ -348,7 +323,8 @@ void Neighbor::skip_from_respa(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
+
   int *ilist_skip = list->listskip->ilist;
   int *numneigh_skip = list->listskip->numneigh;
   int **firstneigh_skip = list->listskip->firstneigh;
@@ -360,30 +336,29 @@ void Neighbor::skip_from_respa(NeighList *list)
   NeighList *listinner = list->listinner;
   int *numneigh_inner = listinner->numneigh;
   int **firstneigh_inner = listinner->firstneigh;
-  int **pages_inner = listinner->pages;
+  MyPage<int> *ipage_inner = listinner->ipage;
+
   int *numneigh_inner_skip = list->listskip->listinner->numneigh;
   int **firstneigh_inner_skip = list->listskip->listinner->firstneigh;
 
   NeighList *listmiddle;
-  int *numneigh_middle,**firstneigh_middle,**pages_middle;
+  int *ilist_middle,*numneigh_middle,**firstneigh_middle;
+  MyPage<int> *ipage_middle;
   int *numneigh_middle_skip,**firstneigh_middle_skip;
   int respamiddle = list->respamiddle;
   if (respamiddle) {
     listmiddle = list->listmiddle;
     numneigh_middle = listmiddle->numneigh;
     firstneigh_middle = listmiddle->firstneigh;
-    pages_middle = listmiddle->pages;
+    ipage_middle = listmiddle->ipage;
     numneigh_middle_skip = list->listskip->listmiddle->numneigh;
     firstneigh_middle_skip = list->listskip->listmiddle->firstneigh;
   }
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
-  int npage_inner = 0;
-  int npnt_inner = 0;
-  int npage_middle = 0;
-  int npnt_middle = 0;
+  ipage->reset();
+  ipage_inner->reset();
+  if (respamiddle) ipage_middle->reset();
 
   // loop over atoms in other list
   // skip I atom entirely if iskip is set for type[I]
@@ -394,32 +369,12 @@ void Neighbor::skip_from_respa(NeighList *list)
     itype = type[i];
     if (iskip[itype]) continue;
 
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-    neighptr = &pages[npage][npnt];
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner++;
-      if (npage_inner == listinner->maxpage)
-        pages_inner = listinner->add_pages();
-    }
-    neighptr_inner = &pages_inner[npage_inner][npnt_inner];
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage->vget();
+    neighptr_inner = ipage_inner->vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle++;
-        if (npage_middle == listmiddle->maxpage)
-          pages_middle = listmiddle->add_pages();
-      }
-      neighptr_middle = &pages_middle[npage_middle][npnt_middle];
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
 
     // loop over parent outer rRESPA list
@@ -460,24 +415,24 @@ void Neighbor::skip_from_respa(NeighList *list)
       }
     }
 
-    ilist[inum++] = i;
+    ilist[inum] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage_inner->vgot(n);
+    if (ipage_inner->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
   }
@@ -499,6 +454,6 @@ void Neighbor::copy_from(NeighList *list)
   list->numneigh = listcopy->numneigh;
   list->firstneigh = listcopy->firstneigh;
   list->firstdouble = listcopy->firstdouble;
-  list->pages = listcopy->pages;
-  list->dpages = listcopy->dpages;
+  list->ipage = listcopy->ipage;
+  list->dpage = listcopy->dpage;
 }
diff --git a/src/neigh_full.cpp b/src/neigh_full.cpp
index 970e3c5725c3de6cd091c22ee1bec0804dc4bd82..0c830b372026cef0f86d8793bb225806967bd393 100644
--- a/src/neigh_full.cpp
+++ b/src/neigh_full.cpp
@@ -15,6 +15,7 @@
 #include "neigh_list.h"
 #include "atom.h"
 #include "domain.h"
+#include "my_page.h"
 #include "group.h"
 #include "error.h"
 
@@ -50,24 +51,16 @@ void Neighbor::full_nsq(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   // loop over owned atoms, storing neighbors
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -101,8 +94,8 @@ void Neighbor::full_nsq(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -137,24 +130,16 @@ void Neighbor::full_nsq_ghost(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   // loop over owned & ghost atoms, storing neighbors
 
   for (i = 0; i < nall; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -203,8 +188,8 @@ void Neighbor::full_nsq_ghost(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -242,26 +227,18 @@ void Neighbor::full_bin(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   // loop over owned atoms, storing neighbors
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -300,8 +277,8 @@ void Neighbor::full_bin(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -341,27 +318,19 @@ void Neighbor::full_bin_ghost(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
   int **stencilxyz = list->stencilxyz;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   // loop over owned & ghost atoms, storing neighbors
 
   for (i = 0; i < nall; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -427,8 +396,8 @@ void Neighbor::full_bin_ghost(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -470,25 +439,17 @@ void Neighbor::full_multi(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int *nstencil_multi = list->nstencil_multi;
   int **stencil_multi = list->stencil_multi;
   double **distsq_multi = list->distsq_multi;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -532,8 +493,8 @@ void Neighbor::full_multi(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
diff --git a/src/neigh_gran.cpp b/src/neigh_gran.cpp
index b9fd604f619b28336572ef5322cab015cc260b7d..204988e23049d31aec75fa2ed91eb269688d74f0 100644
--- a/src/neigh_gran.cpp
+++ b/src/neigh_gran.cpp
@@ -41,8 +41,8 @@ void Neighbor::granular_nsq_no_newton(NeighList *list)
   double (**shearpartner)[3];
   int **firsttouch;
   double **firstshear;
-  int **pages_touch;
-  double **pages_shear;
+  MyPage<int> *ipage_touch;
+  MyPage<double> *dpage_shear;
 
   double **x = atom->x;
   double *radius = atom->radius;
@@ -60,7 +60,7 @@ void Neighbor::granular_nsq_no_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
 
   FixShearHistory *fix_history = list->fix_history;
   if (fix_history) {
@@ -70,34 +70,24 @@ void Neighbor::granular_nsq_no_newton(NeighList *list)
     listgranhistory = list->listgranhistory;
     firsttouch = listgranhistory->firstneigh;
     firstshear = listgranhistory->firstdouble;
-    pages_touch = listgranhistory->pages;
-    pages_shear = listgranhistory->dpages;
+    ipage_touch = listgranhistory->ipage;
+    dpage_shear = listgranhistory->dpage;
   }
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
+  if (fix_history) {
+    ipage_touch->reset();
+    dpage_shear->reset();
+  }
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) {
-        pages = list->add_pages();
-        if (fix_history) {
-          pages_touch = listgranhistory->add_pages();
-          pages_shear = listgranhistory->dpages;
-        }
-      }
-    }
-
     n = 0;
-    neighptr = &pages[npage][npnt];
+    neighptr = ipage->vget();
     if (fix_history) {
       nn = 0;
-      touchptr = &pages_touch[npage][npnt];
-      shearptr = &pages_shear[npage][3*npnt];
+      touchptr = ipage_touch->vget();
+      shearptr = dpage_shear->vget();
     }
 
     xtmp = x[i][0];
@@ -151,13 +141,16 @@ void Neighbor::granular_nsq_no_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
+    ipage->vgot(n);
+    if (ipage->status())
+      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
+
     if (fix_history) {
       firsttouch[i] = touchptr;
       firstshear[i] = shearptr;
+      ipage_touch->vgot(n);
+      dpage_shear->vgot(nn);
     }
-    npnt += n;
-    if (n > oneatom)
-      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
   list->inum = inum;
@@ -195,22 +188,14 @@ void Neighbor::granular_nsq_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
     n = 0;
-    neighptr = &pages[npage][npnt];
+    neighptr = ipage->vget();
 
     itag = tag[i];
     xtmp = x[i][0];
@@ -253,8 +238,8 @@ void Neighbor::granular_nsq_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -283,8 +268,8 @@ void Neighbor::granular_bin_no_newton(NeighList *list)
   double (**shearpartner)[3];
   int **firsttouch;
   double **firstshear;
-  int **pages_touch;
-  double **pages_shear;
+  MyPage<int> *ipage_touch;
+  MyPage<double> *dpage_shear;
 
   // bin local & ghost atoms
 
@@ -304,9 +289,9 @@ void Neighbor::granular_bin_no_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
+  MyPage<int> *ipage = list->ipage;
 
   FixShearHistory *fix_history = list->fix_history;
   if (fix_history) {
@@ -316,34 +301,24 @@ void Neighbor::granular_bin_no_newton(NeighList *list)
     listgranhistory = list->listgranhistory;
     firsttouch = listgranhistory->firstneigh;
     firstshear = listgranhistory->firstdouble;
-    pages_touch = listgranhistory->pages;
-    pages_shear = listgranhistory->dpages;
+    ipage_touch = listgranhistory->ipage;
+    dpage_shear = listgranhistory->dpage;
   }
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
+  if (fix_history) {
+    ipage_touch->reset();
+    dpage_shear->reset();
+  }
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) {
-        pages = list->add_pages();
-        if (fix_history) {
-          pages_touch = listgranhistory->add_pages();
-          pages_shear = listgranhistory->dpages;
-        }
-      }
-    }
-
     n = 0;
-    neighptr = &pages[npage][npnt];
+    neighptr = ipage->vget();
     if (fix_history) {
       nn = 0;
-      touchptr = &pages_touch[npage][npnt];
-      shearptr = &pages_shear[npage][3*npnt];
+      touchptr = ipage_touch->vget();
+      shearptr = dpage_shear->vget();
     }
 
     xtmp = x[i][0];
@@ -403,13 +378,16 @@ void Neighbor::granular_bin_no_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
+    ipage->vgot(n);
+    if (ipage->status())
+      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
+
     if (fix_history) {
       firsttouch[i] = touchptr;
       firstshear[i] = shearptr;
+      ipage_touch->vgot(n);
+      dpage_shear->vgot(nn);
     }
-    npnt += n;
-    if (n > oneatom)
-      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
   list->inum = inum;
@@ -447,24 +425,16 @@ void Neighbor::granular_bin_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
     n = 0;
-    neighptr = &pages[npage][npnt];
+    neighptr = ipage->vget();
 
     xtmp = x[i][0];
     ytmp = x[i][1];
@@ -517,8 +487,8 @@ void Neighbor::granular_bin_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -557,24 +527,16 @@ void Neighbor::granular_bin_newton_tri(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
     n = 0;
-    neighptr = &pages[npage][npnt];
+    neighptr = ipage->vget();
 
     xtmp = x[i][0];
     ytmp = x[i][1];
@@ -615,8 +577,8 @@ void Neighbor::granular_bin_newton_tri(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
diff --git a/src/neigh_half_bin.cpp b/src/neigh_half_bin.cpp
index 9d89274abccb9b8813a089f221dfbac5613b1575..ef55cb5ba145de4461e40f873919bee496946077 100644
--- a/src/neigh_half_bin.cpp
+++ b/src/neigh_half_bin.cpp
@@ -15,6 +15,7 @@
 #include "neigh_list.h"
 #include "atom.h"
 #include "domain.h"
+#include "my_page.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
@@ -53,24 +54,16 @@ void Neighbor::half_bin_no_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -111,8 +104,8 @@ void Neighbor::half_bin_no_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -156,25 +149,17 @@ void Neighbor::half_bin_no_newton_ghost(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
   int **stencilxyz = list->stencilxyz;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nall; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -244,8 +229,8 @@ void Neighbor::half_bin_no_newton_ghost(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -286,24 +271,16 @@ void Neighbor::half_bin_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -372,8 +349,8 @@ void Neighbor::half_bin_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -413,24 +390,16 @@ void Neighbor::half_bin_newton_tri(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -478,8 +447,8 @@ void Neighbor::half_bin_newton_tri(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
diff --git a/src/neigh_half_multi.cpp b/src/neigh_half_multi.cpp
index bee0bfffcd023aa40b3552fd3360fa870bb36046..46a76be6dbd45fa2d58b47feaa3856912b3df0d0 100644
--- a/src/neigh_half_multi.cpp
+++ b/src/neigh_half_multi.cpp
@@ -15,6 +15,7 @@
 #include "neigh_list.h"
 #include "atom.h"
 #include "domain.h"
+#include "my_page.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
@@ -55,25 +56,17 @@ void Neighbor::half_multi_no_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int *nstencil_multi = list->nstencil_multi;
   int **stencil_multi = list->stencil_multi;
   double **distsq_multi = list->distsq_multi;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -119,8 +112,8 @@ void Neighbor::half_multi_no_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -162,25 +155,17 @@ void Neighbor::half_multi_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int *nstencil_multi = list->nstencil_multi;
   int **stencil_multi = list->stencil_multi;
   double **distsq_multi = list->distsq_multi;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -254,8 +239,8 @@ void Neighbor::half_multi_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -297,25 +282,17 @@ void Neighbor::half_multi_newton_tri(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int *nstencil_multi = list->nstencil_multi;
   int **stencil_multi = list->stencil_multi;
   double **distsq_multi = list->distsq_multi;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -370,8 +347,8 @@ void Neighbor::half_multi_newton_tri(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
diff --git a/src/neigh_half_nsq.cpp b/src/neigh_half_nsq.cpp
index 7266df0ebd93de710e762aa24266e83e197b70e2..8446c9d982ff1f1f2065e25d33bb18ac85064a57 100644
--- a/src/neigh_half_nsq.cpp
+++ b/src/neigh_half_nsq.cpp
@@ -51,24 +51,16 @@ void Neighbor::half_nsq_no_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   // loop over owned atoms, storing neighbors
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -102,8 +94,8 @@ void Neighbor::half_nsq_no_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -143,24 +135,16 @@ void Neighbor::half_nsq_no_newton_ghost(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   // loop over owned & ghost atoms, storing neighbors
 
   for (i = 0; i < nall; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itype = type[i];
     xtmp = x[i][0];
@@ -214,8 +198,8 @@ void Neighbor::half_nsq_no_newton_ghost(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
@@ -256,22 +240,14 @@ void Neighbor::half_nsq_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
+  ipage->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-
-    neighptr = &pages[npage][npnt];
     n = 0;
+    neighptr = ipage->vget();
 
     itag = tag[i];
     itype = type[i];
@@ -322,8 +298,8 @@ void Neighbor::half_nsq_newton(NeighList *list)
     ilist[inum++] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
   }
 
diff --git a/src/neigh_list.cpp b/src/neigh_list.cpp
index 316b39ed7746fe3f6159a9894e7100948d36cfa3..aac201b012fe5f539b671f0777d94428a69d0098 100644
--- a/src/neigh_list.cpp
+++ b/src/neigh_list.cpp
@@ -17,6 +17,7 @@
 #include "update.h"
 #include "neighbor.h"
 #include "neigh_request.h"
+#include "my_page.h"
 #include "memory.h"
 #include "error.h"
 
@@ -28,11 +29,10 @@ enum{NSQ,BIN,MULTI};     // also in neighbor.cpp
 
 /* ---------------------------------------------------------------------- */
 
-NeighList::NeighList(LAMMPS *lmp, int size, int onesize) : Pointers(lmp)
+NeighList::NeighList(LAMMPS *lmp) :
+  Pointers(lmp)
 {
   maxatoms = 0;
-  pgsize = size;
-  oneatom = onesize;
 
   inum = gnum = 0;
   ilist = NULL;
@@ -40,11 +40,6 @@ NeighList::NeighList(LAMMPS *lmp, int size, int onesize) : Pointers(lmp)
   firstneigh = NULL;
   firstdouble = NULL;
 
-  maxpage = 0;
-  pages = NULL;
-  dpages = NULL;
-  dnum = 0;
-
   iskip = NULL;
   ijskip = NULL;
 
@@ -66,6 +61,9 @@ NeighList::NeighList(LAMMPS *lmp, int size, int onesize) : Pointers(lmp)
   nstencil_multi = NULL;
   stencil_multi = NULL;
   distsq_multi = NULL;
+
+  ipage = NULL;
+  dpage = NULL;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -78,12 +76,8 @@ NeighList::~NeighList()
     memory->sfree(firstneigh);
     memory->sfree(firstdouble);
 
-    for (int i = 0; i < maxpage; i++) memory->destroy(pages[i]);
-    memory->sfree(pages);
-    if (dnum) {
-      for (int i = 0; i < maxpage; i++) memory->destroy(dpages[i]);
-      memory->sfree(dpages);
-    }
+    delete [] ipage;
+    if (dnum) delete [] dpage;
   }
 
   delete [] iskip;
@@ -103,6 +97,28 @@ NeighList::~NeighList()
   }
 }
 
+/* ---------------------------------------------------------------------- */
+
+void NeighList::setup_pages(int pgsize_caller, int oneatom_caller, 
+                            int dnum_caller)
+{
+  pgsize = pgsize_caller;
+  oneatom = oneatom_caller;
+  dnum = dnum_caller;
+
+  int nmypage = comm->nthreads;
+  ipage = new MyPage<int>[nmypage];
+  for (int i = 0; i < nmypage; i++)
+    ipage[i].init(oneatom,pgsize,PGDELTA);
+
+  if (dnum) {
+    dpage = new MyPage<double>[nmypage];
+    for (int i = 0; i < nmypage; i++)
+      dpage[i].init(dnum*oneatom,dnum*pgsize,PGDELTA);
+  }
+  else dpage = NULL;
+}
+
 /* ----------------------------------------------------------------------
    grow atom arrays to allow for nmax atoms
    triggered by more atoms on a processor
@@ -177,30 +193,6 @@ void NeighList::stencil_allocate(int smax, int style)
   }
 }
 
-/* ----------------------------------------------------------------------
-   add PGDELTA pages to neighbor list
-------------------------------------------------------------------------- */
-
-int **NeighList::add_pages(int howmany)
-{
-  int toppage = maxpage;
-  maxpage += howmany*PGDELTA;
-
-  pages = (int **)
-    memory->srealloc(pages,maxpage*sizeof(int *),"neighlist:pages");
-  for (int i = toppage; i < maxpage; i++)
-    memory->create(pages[i],pgsize,"neighlist:pages[i]");
-
-  if (dnum) {
-    dpages = (double **)
-      memory->srealloc(dpages,maxpage*sizeof(double *),"neighlist:dpages");
-    for (int i = toppage; i < maxpage; i++)
-      memory->create(dpages[i],dnum*pgsize,"neighlist:dpages[i]");
-  }
-
-  return pages;
-}
-
 /* ----------------------------------------------------------------------
    copy skip info from request rq into list's iskip,ijskip
 ------------------------------------------------------------------------- */
@@ -272,11 +264,16 @@ bigint NeighList::memory_usage()
   bytes += memory->usage(ilist,maxatoms);
   bytes += memory->usage(numneigh,maxatoms);
   bytes += maxatoms * sizeof(int *);
-  bytes += memory->usage(pages,maxpage,pgsize);
 
+  int nmypage = comm->nthreads;
+  for (int i = 0; i < nmypage; i++)
+    bytes += ipage[i].size();
+  
   if (dnum) {
-    bytes += maxatoms * sizeof(double *);
-    bytes += memory->usage(dpages,maxpage,dnum*pgsize);
+    for (int i = 0; i < nmypage; i++) {
+      bytes += maxatoms * sizeof(double *);
+      bytes += dpage[i].size();
+    }
   }
 
   if (maxstencil) bytes += memory->usage(stencil,maxstencil);
diff --git a/src/neigh_list.h b/src/neigh_list.h
index d29f3a7db63ae41fea2b33963b5296decc9d0ebf..b5b277c218697842a73f860c54b275f07a71b6fc 100644
--- a/src/neigh_list.h
+++ b/src/neigh_list.h
@@ -41,12 +41,9 @@ class NeighList : protected Pointers {
 
   int pgsize;                      // size of each page
   int oneatom;                     // max size for one atom
-  int maxpage;                     // # of pages currently allocated
-  int **pages;                     // neighbor list pages for ints
-  double **dpages;                 // neighbor list pages for doubles
-  int dnum;                        // # of doubles for each pair (0 if none)
-
-  MyPage<int> *page;
+  int dnum;                        // # of doubles per neighbor, 0 if none
+  MyPage<int> *ipage;              // pages of neighbor indices
+  MyPage<double> *dpage;           // pages of neighbor doubles, if dnum > 0
 
   // atom types to skip when building list
   // iskip,ijskip are just ptrs to corresponding request
@@ -80,11 +77,11 @@ class NeighList : protected Pointers {
 
   class CudaNeighList *cuda_list;  // CUDA neighbor list
 
-  NeighList(class LAMMPS *, int, int);
+  NeighList(class LAMMPS *);
   ~NeighList();
+  void setup_pages(int, int, int);      // setup page data structures
   void grow(int);                       // grow maxlocal
   void stencil_allocate(int, int);      // allocate stencil arrays
-  int **add_pages(int howmany=1);       // add pages to neigh list
   void copy_skip_info(int *, int **);   // copy skip info from a neigh request
   void print_attributes();              // debug routine
   int get_maxlocal() {return maxatoms;}
diff --git a/src/neigh_respa.cpp b/src/neigh_respa.cpp
index c152eba2b8d8cd52ec0642457a5915f27e4e3241..8d1ea14ae5e5d88d2f374082ee3a762659ea4969 100644
--- a/src/neigh_respa.cpp
+++ b/src/neigh_respa.cpp
@@ -16,6 +16,7 @@
 #include "atom.h"
 #include "domain.h"
 #include "group.h"
+#include "my_page.h"
 #include "error.h"
 
 using namespace LAMMPS_NS;
@@ -54,64 +55,40 @@ void Neighbor::respa_nsq_no_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
 
   NeighList *listinner = list->listinner;
   int *ilist_inner = listinner->ilist;
   int *numneigh_inner = listinner->numneigh;
   int **firstneigh_inner = listinner->firstneigh;
-  int **pages_inner = listinner->pages;
+  MyPage<int> *ipage_inner = listinner->ipage;
 
   NeighList *listmiddle;
-  int *ilist_middle,*numneigh_middle,**firstneigh_middle,**pages_middle;
+  int *ilist_middle,*numneigh_middle,**firstneigh_middle;
+  MyPage<int> *ipage_middle;
   int respamiddle = list->respamiddle;
   if (respamiddle) {
     listmiddle = list->listmiddle;
     ilist_middle = listmiddle->ilist;
     numneigh_middle = listmiddle->numneigh;
     firstneigh_middle = listmiddle->firstneigh;
-    pages_middle = listmiddle->pages;
+    ipage_middle = listmiddle->ipage;
   }
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
-  int npage_inner = 0;
-  int npnt_inner = 0;
-  int npage_middle = 0;
-  int npnt_middle = 0;
-
   int which = 0;
   int minchange = 0;
+  ipage->reset();
+  ipage_inner->reset();
+  if (respamiddle) ipage_middle->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-    neighptr = &pages[npage][npnt];
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner++;
-      if (npage_inner == listinner->maxpage)
-        pages_inner = listinner->add_pages();
-    }
-    neighptr_inner = &pages_inner[npage_inner][npnt_inner];
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage->vget();
+    neighptr_inner = ipage_inner->vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle++;
-        if (npage_middle == listmiddle->maxpage)
-          pages_middle = listmiddle->add_pages();
-      }
-      neighptr_middle = &pages_middle[npage_middle][npnt_middle];
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
 
     itype = type[i];
@@ -158,23 +135,23 @@ void Neighbor::respa_nsq_no_newton(NeighList *list)
     ilist[inum] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     ilist_inner[inum] = i;
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage_inner->vgot(n_inner);
+    if (ipage_inner->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       ilist_middle[inum] = i;
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n_middle);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
 
@@ -221,64 +198,40 @@ void Neighbor::respa_nsq_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
+  MyPage<int> *ipage = list->ipage;
 
   NeighList *listinner = list->listinner;
   int *ilist_inner = listinner->ilist;
   int *numneigh_inner = listinner->numneigh;
   int **firstneigh_inner = listinner->firstneigh;
-  int **pages_inner = listinner->pages;
+  MyPage<int> *ipage_inner = listinner->ipage;
 
   NeighList *listmiddle;
-  int *ilist_middle,*numneigh_middle,**firstneigh_middle,**pages_middle;
+  int *ilist_middle,*numneigh_middle,**firstneigh_middle;
+  MyPage<int> *ipage_middle;
   int respamiddle = list->respamiddle;
   if (respamiddle) {
     listmiddle = list->listmiddle;
     ilist_middle = listmiddle->ilist;
     numneigh_middle = listmiddle->numneigh;
     firstneigh_middle = listmiddle->firstneigh;
-    pages_middle = listmiddle->pages;
+    ipage_middle = listmiddle->ipage;
   }
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
-  int npage_inner = 0;
-  int npnt_inner = 0;
-  int npage_middle = 0;
-  int npnt_middle = 0;
-
   int which = 0;
   int minchange = 0;
+  ipage->reset();
+  ipage_inner->reset();
+  if (respamiddle) ipage_middle->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-    neighptr = &pages[npage][npnt];
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner++;
-      if (npage_inner == listinner->maxpage)
-        pages_inner = listinner->add_pages();
-    }
-    neighptr_inner = &pages_inner[npage_inner][npnt_inner];
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage->vget();
+    neighptr_inner = ipage_inner->vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle++;
-        if (npage_middle == listmiddle->maxpage)
-          pages_middle = listmiddle->add_pages();
-      }
-      neighptr_middle = &pages_middle[npage_middle][npnt_middle];
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
 
     itag = tag[i];
@@ -343,23 +296,23 @@ void Neighbor::respa_nsq_newton(NeighList *list)
     ilist[inum] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     ilist_inner[inum] = i;
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage_inner->vgot(n_inner);
+    if (ipage_inner->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       ilist_middle[inum] = i;
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n_middle);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
 
@@ -406,66 +359,42 @@ void Neighbor::respa_bin_no_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
+  MyPage<int> *ipage = list->ipage;
 
   NeighList *listinner = list->listinner;
   int *ilist_inner = listinner->ilist;
   int *numneigh_inner = listinner->numneigh;
   int **firstneigh_inner = listinner->firstneigh;
-  int **pages_inner = listinner->pages;
+  MyPage<int> *ipage_inner = listinner->ipage;
 
   NeighList *listmiddle;
-  int *ilist_middle,*numneigh_middle,**firstneigh_middle,**pages_middle;
+  int *ilist_middle,*numneigh_middle,**firstneigh_middle;
+  MyPage<int> *ipage_middle;
   int respamiddle = list->respamiddle;
   if (respamiddle) {
     listmiddle = list->listmiddle;
     ilist_middle = listmiddle->ilist;
     numneigh_middle = listmiddle->numneigh;
     firstneigh_middle = listmiddle->firstneigh;
-    pages_middle = listmiddle->pages;
+    ipage_middle = listmiddle->ipage;
   }
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
-  int npage_inner = 0;
-  int npnt_inner = 0;
-  int npage_middle = 0;
-  int npnt_middle = 0;
-
   int which = 0;
   int minchange = 0;
+  ipage->reset();
+  ipage_inner->reset();
+  if (respamiddle) ipage_middle->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-    neighptr = &pages[npage][npnt];
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner++;
-      if (npage_inner == listinner->maxpage)
-        pages_inner = listinner->add_pages();
-    }
-    neighptr_inner = &pages_inner[npage_inner][npnt_inner];
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage->vget();
+    neighptr_inner = ipage_inner->vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle++;
-        if (npage_middle == listmiddle->maxpage)
-          pages_middle = listmiddle->add_pages();
-      }
-      neighptr_middle = &pages_middle[npage_middle][npnt_middle];
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
 
     itype = type[i];
@@ -521,23 +450,23 @@ void Neighbor::respa_bin_no_newton(NeighList *list)
     ilist[inum] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     ilist_inner[inum] = i;
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage_inner->vgot(n_inner);
+    if (ipage_inner->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       ilist_middle[inum] = i;
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n_middle);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
 
@@ -583,66 +512,42 @@ void Neighbor::respa_bin_newton(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
+  MyPage<int> *ipage = list->ipage;
 
   NeighList *listinner = list->listinner;
   int *ilist_inner = listinner->ilist;
   int *numneigh_inner = listinner->numneigh;
   int **firstneigh_inner = listinner->firstneigh;
-  int **pages_inner = listinner->pages;
+  MyPage<int> *ipage_inner = listinner->ipage;
 
   NeighList *listmiddle;
-  int *ilist_middle,*numneigh_middle,**firstneigh_middle,**pages_middle;
+  int *ilist_middle,*numneigh_middle,**firstneigh_middle;
+  MyPage<int> *ipage_middle;
   int respamiddle = list->respamiddle;
   if (respamiddle) {
     listmiddle = list->listmiddle;
     ilist_middle = listmiddle->ilist;
     numneigh_middle = listmiddle->numneigh;
     firstneigh_middle = listmiddle->firstneigh;
-    pages_middle = listmiddle->pages;
+    ipage_middle = listmiddle->ipage;
   }
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
-  int npage_inner = 0;
-  int npnt_inner = 0;
-  int npage_middle = 0;
-  int npnt_middle = 0;
-
   int which = 0;
   int minchange = 0;
+  ipage->reset();
+  ipage_inner->reset();
+  if (respamiddle) ipage_middle->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-    neighptr = &pages[npage][npnt];
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner++;
-      if (npage_inner == listinner->maxpage)
-        pages_inner = listinner->add_pages();
-    }
-    neighptr_inner = &pages_inner[npage_inner][npnt_inner];
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage->vget();
+    neighptr_inner = ipage_inner->vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle++;
-        if (npage_middle == listmiddle->maxpage)
-          pages_middle = listmiddle->add_pages();
-      }
-      neighptr_middle = &pages_middle[npage_middle][npnt_middle];
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
 
     itype = type[i];
@@ -739,23 +644,23 @@ void Neighbor::respa_bin_newton(NeighList *list)
     ilist[inum] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     ilist_inner[inum] = i;
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage_inner->vgot(n_inner);
+    if (ipage_inner->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       ilist_middle[inum] = i;
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n_middle);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
 
@@ -801,66 +706,42 @@ void Neighbor::respa_bin_newton_tri(NeighList *list)
   int *ilist = list->ilist;
   int *numneigh = list->numneigh;
   int **firstneigh = list->firstneigh;
-  int **pages = list->pages;
   int nstencil = list->nstencil;
   int *stencil = list->stencil;
+  MyPage<int> *ipage = list->ipage;
 
   NeighList *listinner = list->listinner;
   int *ilist_inner = listinner->ilist;
   int *numneigh_inner = listinner->numneigh;
   int **firstneigh_inner = listinner->firstneigh;
-  int **pages_inner = listinner->pages;
+  MyPage<int> *ipage_inner = listinner->ipage;
 
   NeighList *listmiddle;
-  int *ilist_middle,*numneigh_middle,**firstneigh_middle,**pages_middle;
+  int *ilist_middle,*numneigh_middle,**firstneigh_middle;
+  MyPage<int> *ipage_middle;
   int respamiddle = list->respamiddle;
   if (respamiddle) {
     listmiddle = list->listmiddle;
     ilist_middle = listmiddle->ilist;
     numneigh_middle = listmiddle->numneigh;
     firstneigh_middle = listmiddle->firstneigh;
-    pages_middle = listmiddle->pages;
+    ipage_middle = listmiddle->ipage;
   }
 
   int inum = 0;
-  int npage = 0;
-  int npnt = 0;
-  int npage_inner = 0;
-  int npnt_inner = 0;
-  int npage_middle = 0;
-  int npnt_middle = 0;
-
   int which = 0;
   int minchange = 0;
+  ipage->reset();
+  ipage_inner->reset();
+  if (respamiddle) ipage_middle->reset();
 
   for (i = 0; i < nlocal; i++) {
-
-    if (pgsize - npnt < oneatom) {
-      npnt = 0;
-      npage++;
-      if (npage == list->maxpage) pages = list->add_pages();
-    }
-    neighptr = &pages[npage][npnt];
-    n = 0;
-
-    if (pgsize - npnt_inner < oneatom) {
-      npnt_inner = 0;
-      npage_inner++;
-      if (npage_inner == listinner->maxpage)
-        pages_inner = listinner->add_pages();
-    }
-    neighptr_inner = &pages_inner[npage_inner][npnt_inner];
-    n_inner = 0;
-
+    n = n_inner = 0;
+    neighptr = ipage->vget();
+    neighptr_inner = ipage_inner->vget();
     if (respamiddle) {
-      if (pgsize - npnt_middle < oneatom) {
-        npnt_middle = 0;
-        npage_middle++;
-        if (npage_middle == listmiddle->maxpage)
-          pages_middle = listmiddle->add_pages();
-      }
-      neighptr_middle = &pages_middle[npage_middle][npnt_middle];
       n_middle = 0;
+      neighptr_middle = ipage_middle->vget();
     }
 
     itype = type[i];
@@ -924,23 +805,23 @@ void Neighbor::respa_bin_newton_tri(NeighList *list)
     ilist[inum] = i;
     firstneigh[i] = neighptr;
     numneigh[i] = n;
-    npnt += n;
-    if (n > oneatom)
+    ipage->vgot(n);
+    if (ipage->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     ilist_inner[inum] = i;
     firstneigh_inner[i] = neighptr_inner;
     numneigh_inner[i] = n_inner;
-    npnt_inner += n_inner;
-    if (n_inner > oneatom)
+    ipage_inner->vgot(n_inner);
+    if (ipage_inner->status())
       error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
 
     if (respamiddle) {
       ilist_middle[inum] = i;
       firstneigh_middle[i] = neighptr_middle;
       numneigh_middle[i] = n_middle;
-      npnt_middle += n_middle;
-      if (n_middle > oneatom)
+      ipage_middle->vgot(n_middle);
+      if (ipage_middle->status())
         error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
     }
 
diff --git a/src/neighbor.cpp b/src/neighbor.cpp
index 0f597b64d773d1b6d02f5ada9ba71b1eebaa69cf..b9ecb4c47136318fbeb9902a40f5548d69268161 100644
--- a/src/neighbor.cpp
+++ b/src/neighbor.cpp
@@ -456,9 +456,9 @@ void Neighbor::init()
     // pass list ptr back to requestor (except for Command class)
 
     for (i = 0; i < nlist; i++) {
-      lists[i] = new NeighList(lmp,pgsize,oneatom);
+      lists[i] = new NeighList(lmp);
+      lists[i]->setup_pages(pgsize,oneatom,requests[i]->dnum);
       lists[i]->index = i;
-      lists[i]->dnum = requests[i]->dnum;
 
       if (requests[i]->pair) {
         Pair *pair = (Pair *) requests[i]->requestor;
@@ -620,14 +620,11 @@ void Neighbor::init()
     for (i = 0; i < nlist; i++) lists[i]->print_attributes();
 #endif
 
-    // allocate atom arrays and 1st pages of lists that store them
+    // allocate atom arrays for neighbor lists that store them
 
     maxatom = atom->nmax;
     for (i = 0; i < nlist; i++)
-      if (lists[i]->growflag) {
-        lists[i]->grow(maxatom);
-        lists[i]->add_pages();
-      }
+      if (lists[i]->growflag) lists[i]->grow(maxatom);
 
     // setup 3 vectors of pairwise neighbor lists
     // blist = lists whose pair_build() is invoked every reneighbor