From 2fa8f3d114b7c0364ba458dd3021aea8d723fe0d Mon Sep 17 00:00:00 2001
From: sjplimp <sjplimp@f3b2605a-c512-4ea7-a41b-209d697bcdaa>
Date: Thu, 8 Mar 2007 22:07:17 +0000
Subject: [PATCH] git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@378
 f3b2605a-c512-4ea7-a41b-209d697bcdaa

---
 src/comm.cpp         |  26 +++++----
 src/comm.h           |   1 +
 src/create_atoms.cpp |  18 +++---
 src/domain.cpp       |  82 +++++++++++++++++++-------
 src/domain.h         |   5 +-
 src/fix_recenter.cpp |   2 +-
 src/neighbor.cpp     | 133 ++++++++++++++++++++++++-------------------
 src/neighbor.h       |   2 +-
 8 files changed, 164 insertions(+), 105 deletions(-)

diff --git a/src/comm.cpp b/src/comm.cpp
index be2116f190..42f267ec08 100644
--- a/src/comm.cpp
+++ b/src/comm.cpp
@@ -171,17 +171,19 @@ void Comm::set_procs()
 
 void Comm::setup()
 {
-  // for triclinic, tilted planes are neigh cutoff apart
-  // but cutneigh is distance between those planes along xyz (non-tilted) axes
+  // cutcomm = distance at which ghost atoms need to be acquired
+  // for orthogonal, cutcomm = box coords = cutneigh in all 3 dims
+  // for triclinic, cutneigh = distance between tilted planes in box coords
+  //   cutcomm = lamda coords = distance between those same planes
+  //                            can be different for each dim
 
-  double cutneigh[3];
   double *prd,*prd_border,*sublo,*subhi;
 
   if (triclinic == 0) {
     prd = prd_border = domain->prd;
     sublo = domain->sublo;
     subhi = domain->subhi;
-    cutneigh[0] = cutneigh[1] = cutneigh[2] = neighbor->cutneigh;
+    cutcomm[0] = cutcomm[1] = cutcomm[2] = neighbor->cutneigh;
   } else {
     prd = domain->prd;
     prd_border = domain->prd_lamda;
@@ -190,19 +192,19 @@ void Comm::setup()
     double *h_inv = domain->h_inv;
     double length;
     length = sqrt(h_inv[0]*h_inv[0] + h_inv[5]*h_inv[5] + h_inv[4]*h_inv[4]);
-    cutneigh[0] = neighbor->cutneigh * length;
+    cutcomm[0] = neighbor->cutneigh * length;
     length = sqrt(h_inv[1]*h_inv[1] + h_inv[3]*h_inv[3]);
-    cutneigh[1] = neighbor->cutneigh * length;
+    cutcomm[1] = neighbor->cutneigh * length;
     length = h_inv[2];
-    cutneigh[2] = neighbor->cutneigh * length;
+    cutcomm[2] = neighbor->cutneigh * length;
   }
 
   // need = # of procs I need atoms from in each dim
   // for 2d, don't communicate in z
 
-  need[0] = static_cast<int> (cutneigh[0] * procgrid[0] / prd_border[0]) + 1;
-  need[1] = static_cast<int> (cutneigh[1] * procgrid[1] / prd_border[1]) + 1;
-  need[2] = static_cast<int> (cutneigh[2] * procgrid[2] / prd_border[2]) + 1;
+  need[0] = static_cast<int> (cutcomm[0] * procgrid[0] / prd_border[0]) + 1;
+  need[1] = static_cast<int> (cutcomm[1] * procgrid[1] / prd_border[1]) + 1;
+  need[2] = static_cast<int> (cutcomm[2] * procgrid[2] / prd_border[2]) + 1;
   if (force->dimension == 2) need[2] = 0;
 
   // if non-periodic, do not communicate further than procgrid-1 away
@@ -247,7 +249,7 @@ void Comm::setup()
 	recvproc[iswap] = procneigh[dim][1];
 	if (ineed < 2) slablo[iswap] = -BIG;
 	else slablo[iswap] = 0.5 * (sublo[dim] + subhi[dim]);
-	slabhi[iswap] = sublo[dim] + cutneigh[dim];
+	slabhi[iswap] = sublo[dim] + cutcomm[dim];
 	if (myloc[dim] == 0) {
 	  if (periodicity[dim] == 0)
 	    slabhi[iswap] = slablo[iswap] - 1.0;
@@ -267,7 +269,7 @@ void Comm::setup()
       } else {
 	sendproc[iswap] = procneigh[dim][1];
 	recvproc[iswap] = procneigh[dim][0];
-	slablo[iswap] = subhi[dim] - cutneigh[dim];
+	slablo[iswap] = subhi[dim] - cutcomm[dim];
 	if (ineed < 2) slabhi[iswap] = BIG;
 	else slabhi[iswap] = 0.5 * (sublo[dim] + subhi[dim]);
 	if (myloc[dim] == procgrid[dim]-1) {
diff --git a/src/comm.h b/src/comm.h
index a6eaf8ce04..3725ebb975 100644
--- a/src/comm.h
+++ b/src/comm.h
@@ -29,6 +29,7 @@ class Comm : protected Pointers {
   int need[3];                      // procs I need atoms from in each dim
   int maxforward_fix,maxreverse_fix;   // comm sizes called from Fix,Pair
   int maxforward_pair,maxreverse_pair;
+  double cutcomm[3];                // cutoffs used for acquiring ghost atoms
   
   Comm(class LAMMPS *);
   ~Comm();
diff --git a/src/create_atoms.cpp b/src/create_atoms.cpp
index e422894525..33bf007b5f 100644
--- a/src/create_atoms.cpp
+++ b/src/create_atoms.cpp
@@ -90,19 +90,19 @@ void CreateAtoms::command(int narg, char **arg)
     } else error->all("Illegal create_atoms command");
   }
 
-  // convert 8 corners of my sub-box from box coords to lattice coords
-  // min to max = bounding box around the pts in lattice space
+  // convert 8 corners of my subdomain from box coords to lattice coords
+  // for orthogonal, use corner pts of my subbox
+  // for triclinic, use bounding box of my subbox
+  // xyz min to max = bounding box around the domain corners in lattice space
 
   int triclinic = domain->triclinic;
-  double *bboxlo,*bboxhi;
+  double bboxlo[3],bboxhi[3];
 
   if (triclinic == 0) {
-    bboxlo = domain->sublo;
-    bboxhi = domain->subhi;
-  } else {
-    bboxlo = domain->sublo_bound;
-    bboxhi = domain->subhi_bound;
-  }
+    bboxlo[0] = domain->sublo[0]; bboxhi[0] = domain->subhi[0];
+    bboxlo[1] = domain->sublo[1]; bboxhi[1] = domain->subhi[1];
+    bboxlo[2] = domain->sublo[2]; bboxhi[2] = domain->subhi[2];
+  } else domain->bbox(domain->sublo_lamda,domain->subhi_lamda,bboxlo,bboxhi);
 
   double xmin,ymin,zmin,xmax,ymax,zmax;
   xmin = ymin = zmin = BIG;
diff --git a/src/domain.cpp b/src/domain.cpp
index 699f612c68..99ae3a4d1f 100644
--- a/src/domain.cpp
+++ b/src/domain.cpp
@@ -219,7 +219,6 @@ void Domain::set_lamda_box()
    set local subbox params
    assumes global box is defined and proc assignment has been made
    for uppermost proc, insure subhi = boxhi (in case round-off occurs)
-   for triclinic, lo/hi_bound is a bbox around all 8 tilted pts of sub-box
 ------------------------------------------------------------------------- */
 
 void Domain::set_local_box()
@@ -242,26 +241,6 @@ void Domain::set_local_box()
     if (myloc[2] < procgrid[2]-1)
       subhi[2] = boxlo[2] + (myloc[2]+1) * zprd / procgrid[2];
     else subhi[2] = boxhi[2];
-
-  } else {
-    sublo[0] = h[0]*sublo_lamda[0] + h[5]*sublo_lamda[1] + 
-      h[4]*sublo_lamda[2] + boxlo[0];
-    sublo[1] = h[1]*sublo_lamda[1] + h[3]*sublo_lamda[2] + boxlo[1];
-    sublo[2] = h[2]*sublo_lamda[2] + boxlo[2];
-
-    subhi[0] = sublo[0] + xprd/procgrid[0];
-    subhi[1] = sublo[1] + yprd/procgrid[1];
-    subhi[2] = sublo[2] + zprd/procgrid[2];
-
-    sublo_bound[0] = MIN(sublo[0],sublo[0] + xy/procgrid[1]);
-    sublo_bound[0] = MIN(sublo_bound[0],sublo_bound[0] + xz/procgrid[2]);
-    sublo_bound[1] = MIN(sublo[1],sublo[1] + yz/procgrid[2]);
-    sublo_bound[2] = sublo[2];
-
-    subhi_bound[0] = MAX(subhi[0],subhi[0] + xy/procgrid[1]);
-    subhi_bound[0] = MAX(subhi_bound[0],subhi_bound[0] + xz/procgrid[2]);
-    subhi_bound[1] = MAX(subhi[1],subhi[1] + yz/procgrid[2]);
-    subhi_bound[2] = subhi[2];
   }
 }
 
@@ -861,3 +840,64 @@ void Domain::x2lamda(double *x, double *lamda)
   lamda[1] = h_inv[1]*delta[1] + h_inv[3]*delta[2];
   lamda[2] = h_inv[2]*delta[2];
 }
+
+/* ----------------------------------------------------------------------
+   convert 8 lamda corner pts of lo/hi box to box coords
+   return bboxlo/hi = bounding box around 8 corner pts in box coords
+------------------------------------------------------------------------- */
+
+void Domain::bbox(double *lo, double *hi, double *bboxlo, double *bboxhi)
+{
+  double x[3];
+
+  bboxlo[0] = bboxlo[1] = bboxlo[2] = BIG;
+  bboxhi[0] = bboxhi[1] = bboxhi[2] = -BIG;
+
+  x[0] = lo[0]; x[1] = lo[1]; x[2] = lo[2];
+  lamda2x(x,x);
+  bboxlo[0] = MIN(bboxlo[0],x[0]); bboxhi[0] = MAX(bboxhi[0],x[0]);
+  bboxlo[1] = MIN(bboxlo[1],x[1]); bboxhi[1] = MAX(bboxhi[1],x[1]);
+  bboxlo[2] = MIN(bboxlo[2],x[2]); bboxhi[2] = MAX(bboxhi[2],x[2]);
+
+  x[0] = hi[0]; x[1] = lo[1]; x[2] = lo[2];
+  lamda2x(x,x);
+  bboxlo[0] = MIN(bboxlo[0],x[0]); bboxhi[0] = MAX(bboxhi[0],x[0]);
+  bboxlo[1] = MIN(bboxlo[1],x[1]); bboxhi[1] = MAX(bboxhi[1],x[1]);
+  bboxlo[2] = MIN(bboxlo[2],x[2]); bboxhi[2] = MAX(bboxhi[2],x[2]);
+
+  x[0] = lo[0]; x[1] = hi[1]; x[2] = lo[2];
+  lamda2x(x,x);
+  bboxlo[0] = MIN(bboxlo[0],x[0]); bboxhi[0] = MAX(bboxhi[0],x[0]);
+  bboxlo[1] = MIN(bboxlo[1],x[1]); bboxhi[1] = MAX(bboxhi[1],x[1]);
+  bboxlo[2] = MIN(bboxlo[2],x[2]); bboxhi[2] = MAX(bboxhi[2],x[2]);
+
+  x[0] = hi[0]; x[1] = hi[1]; x[2] = lo[2];
+  lamda2x(x,x);
+  bboxlo[0] = MIN(bboxlo[0],x[0]); bboxhi[0] = MAX(bboxhi[0],x[0]);
+  bboxlo[1] = MIN(bboxlo[1],x[1]); bboxhi[1] = MAX(bboxhi[1],x[1]);
+  bboxlo[2] = MIN(bboxlo[2],x[2]); bboxhi[2] = MAX(bboxhi[2],x[2]);
+
+  x[0] = lo[0]; x[1] = lo[1]; x[2] = hi[2];
+  lamda2x(x,x);
+  bboxlo[0] = MIN(bboxlo[0],x[0]); bboxhi[0] = MAX(bboxhi[0],x[0]);
+  bboxlo[1] = MIN(bboxlo[1],x[1]); bboxhi[1] = MAX(bboxhi[1],x[1]);
+  bboxlo[2] = MIN(bboxlo[2],x[2]); bboxhi[2] = MAX(bboxhi[2],x[2]);
+
+  x[0] = hi[0]; x[1] = lo[1]; x[2] = hi[2];
+  lamda2x(x,x);
+  bboxlo[0] = MIN(bboxlo[0],x[0]); bboxhi[0] = MAX(bboxhi[0],x[0]);
+  bboxlo[1] = MIN(bboxlo[1],x[1]); bboxhi[1] = MAX(bboxhi[1],x[1]);
+  bboxlo[2] = MIN(bboxlo[2],x[2]); bboxhi[2] = MAX(bboxhi[2],x[2]);
+
+  x[0] = lo[0]; x[1] = hi[1]; x[2] = hi[2];
+  lamda2x(x,x);
+  bboxlo[0] = MIN(bboxlo[0],x[0]); bboxhi[0] = MAX(bboxhi[0],x[0]);
+  bboxlo[1] = MIN(bboxlo[1],x[1]); bboxhi[1] = MAX(bboxhi[1],x[1]);
+  bboxlo[2] = MIN(bboxlo[2],x[2]); bboxhi[2] = MAX(bboxhi[2],x[2]);
+
+  x[0] = hi[0]; x[1] = hi[1]; x[2] = hi[2];
+  lamda2x(x,x);
+  bboxlo[0] = MIN(bboxlo[0],x[0]); bboxhi[0] = MAX(bboxhi[0],x[0]);
+  bboxlo[1] = MIN(bboxlo[1],x[1]); bboxhi[1] = MAX(bboxhi[1],x[1]);
+  bboxlo[2] = MIN(bboxlo[2],x[2]); bboxhi[2] = MAX(bboxhi[2],x[2]);
+}
diff --git a/src/domain.h b/src/domain.h
index 16a6949271..ebde41d841 100644
--- a/src/domain.h
+++ b/src/domain.h
@@ -66,10 +66,8 @@ class Domain : protected Pointers {
   double sublo[3],subhi[3];              // sub-box bounds on this proc
 
                                          // triclinic box
-                                         // sublo = lower-left corner with tilt
-                                         // subhi = upper-right as if untilted
+                                         // sublo/hi = undefined
   double sublo_lamda[3],subhi_lamda[3];  // bounds of subbox in lamda
-  double sublo_bound[3],subhi_bound[3];  // bounding box of tilted subdomain
 
                                          // triclinic box
   double xy,xz,yz;                       // 3 tilt factors
@@ -105,6 +103,7 @@ class Domain : protected Pointers {
   void x2lamda(int);
   void lamda2x(double *, double *);
   void x2lamda(double *, double *);
+  void bbox(double *, double *, double *, double *);
 };
 
 }
diff --git a/src/fix_recenter.cpp b/src/fix_recenter.cpp
index f870a00eab..88f8f2327e 100644
--- a/src/fix_recenter.cpp
+++ b/src/fix_recenter.cpp
@@ -135,7 +135,7 @@ void FixRecenter::init()
 void FixRecenter::initial_integrate()
 {
   // target COM
-  // bounding box around domain works for both orthog and triclinic
+  // bounding box around domain works for both orthogonal and triclinic
 
   double xtarget,ytarget,ztarget;
   double *bboxlo,*bboxhi;
diff --git a/src/neighbor.cpp b/src/neighbor.cpp
index dea05f5244..41f6efd494 100644
--- a/src/neighbor.cpp
+++ b/src/neighbor.cpp
@@ -19,6 +19,7 @@
 #include "neighbor.h"
 #include "atom.h"
 #include "atom_vec.h"
+#include "comm.h"
 #include "force.h"
 #include "pair.h"
 #include "domain.h"
@@ -868,15 +869,16 @@ void Neighbor::build_full()
 /* ----------------------------------------------------------------------
    setup neighbor binning parameters
    bin numbering is global: 0 = 0.0 to binsize, 1 = binsize to 2*binsize
-			    nbin-1 = prd-binsize to binsize
-			    nbin = prd to prd+binsize
-                            -1 = -binsize to 0.0
+			    nbin-1 = bbox-binsize to binsize
+			    nbin = bbox to bbox+binsize
+                            -1,-2,etc = -binsize to 0.0, -2*size to -size, etc
    code will work for any binsize
      since next(xyz) and stencil extend as far as necessary
      binsize = 1/2 of cutoff is roughly optimal
-   prd must be filled exactly by integer # of bins
+   for orthogonal boxes, prd must be filled exactly by integer # of bins
      so procs on both sides of PBC see same bin boundary
-     triclinic is an exception, since stencil & neigh list built differently
+   for triclinic, tilted simulation box cannot contain integer # of bins
+     stencil & neigh list built differently to account for this
    mbinlo = lowest global bin any of my ghost atoms could fall into
    mbinhi = highest global bin any of my ghost atoms could fall into
    mbin = number of bins I need in a dimension
@@ -886,28 +888,41 @@ void Neighbor::build_full()
 
 void Neighbor::setup_bins()
 {
-  // bbox = bounding box of entire domain
-  // bboxlo,bboxhi = bounding box of proc's sub-domain
-  // for triclinic, bounding box surrounds all 8 corner pts
+  // bbox lo/hi = bounding box of entire domain
+  // bbox = size of bbox of entire domain
+  // for triclinic, bbox bounds all 8 corners of tilted box
+  // bsubbox lo/hi = bounding box of my subdomain extended by cutneigh
+  // for triclinic, subdomain with cutneigh extension is first computed
+  //   in lamda coords, including tilt factor via cutcomm,
+  //   domain->bbox then computes bbox of corner pts transformed to box coords
 
-  double bbox[3];
-  double *bboxlo,*bboxhi;
+  double bbox[3],bsubboxlo[3],bsubboxhi[3];
 
   if (triclinic == 0) {
-    boxlo = domain->boxlo;
-    boxhi = domain->boxhi;
-    bboxlo = domain->sublo;
-    bboxhi = domain->subhi;
+    bboxlo = domain->boxlo;
+    bboxhi = domain->boxhi;
+    bsubboxlo[0] = domain->sublo[0] - cutneigh;
+    bsubboxlo[1] = domain->sublo[1] - cutneigh;
+    bsubboxlo[2] = domain->sublo[2] - cutneigh;
+    bsubboxhi[0] = domain->subhi[0] + cutneigh;
+    bsubboxhi[1] = domain->subhi[1] + cutneigh;
+    bsubboxhi[2] = domain->subhi[2] + cutneigh;
   } else {
-    boxlo = domain->boxlo_bound;
-    boxhi = domain->boxhi_bound;
-    bboxlo = domain->sublo_bound;
-    bboxhi = domain->subhi_bound;
+    bboxlo = domain->boxlo_bound;
+    bboxhi = domain->boxhi_bound;
+    double lo[3],hi[3];
+    lo[0] = domain->sublo_lamda[0] - comm->cutcomm[0];
+    lo[1] = domain->sublo_lamda[1] - comm->cutcomm[1];
+    lo[2] = domain->sublo_lamda[2] - comm->cutcomm[2];
+    hi[0] = domain->subhi_lamda[0] + comm->cutcomm[0];
+    hi[1] = domain->subhi_lamda[1] + comm->cutcomm[1];
+    hi[2] = domain->subhi_lamda[2] + comm->cutcomm[2];
+    domain->bbox(lo,hi,bsubboxlo,bsubboxhi);
   }
 
-  bbox[0] = boxhi[0] - boxlo[0];
-  bbox[1] = boxhi[1] - boxlo[1];
-  bbox[2] = boxhi[2] - boxlo[2];
+  bbox[0] = bboxhi[0] - bboxlo[0];
+  bbox[1] = bboxhi[1] - bboxlo[1];
+  bbox[2] = bboxhi[2] - bboxlo[2];
 
   // test for too many global bins in any dimension due to huge domain
 
@@ -919,6 +934,7 @@ void Neighbor::setup_bins()
 
   // divide box into bins
   // optimal size is roughly 1/2 the cutoff
+  // use one bin even if cutoff >> bbox
 
   nbinx = static_cast<int> (2.0*bbox[0]*cutneighinv);
   nbiny = static_cast<int> (2.0*bbox[1]*cutneighinv);
@@ -940,28 +956,29 @@ void Neighbor::setup_bins()
 
   // mbinlo/hi = lowest and highest global bins my ghost atoms could be in
   // coord = lowest and highest values of coords for my ghost atoms
-  //         add in SMALL for round-off safety
+  // static_cast(-1.5) = -1, so subract additional -1
+  // add in SMALL for round-off safety
 
-  double coord;
   int mbinxhi,mbinyhi,mbinzhi;
+  double coord;
 
-  coord = bboxlo[0] - cutneigh - SMALL*bbox[0];
-  mbinxlo = static_cast<int> ((coord-boxlo[0])*bininvx);
-  if (coord < 0.0) mbinxlo = mbinxlo - 1;
-  coord = bboxhi[0] + cutneigh + SMALL*bbox[0];
-  mbinxhi = static_cast<int> ((coord-boxlo[0])*bininvx);
+  coord = bsubboxlo[0] - SMALL*bbox[0];
+  mbinxlo = static_cast<int> ((coord-bboxlo[0])*bininvx);
+  if (coord < bboxlo[0]) mbinxlo = mbinxlo - 1;
+  coord = bsubboxhi[0] + SMALL*bbox[0];
+  mbinxhi = static_cast<int> ((coord-bboxlo[0])*bininvx);
 
-  coord = bboxlo[1] - cutneigh - SMALL*bbox[1];
-  mbinylo = static_cast<int> ((coord-boxlo[1])*bininvy);
-  if (coord < 0.0) mbinylo = mbinylo - 1;
-  coord = bboxhi[1] + cutneigh + SMALL*bbox[1];
-  mbinyhi = static_cast<int> ((coord-boxlo[1])*bininvy);
+  coord = bsubboxlo[1] - SMALL*bbox[1];
+  mbinylo = static_cast<int> ((coord-bboxlo[1])*bininvy);
+  if (coord < bboxlo[1]) mbinylo = mbinylo - 1;
+  coord = bsubboxhi[1] + SMALL*bbox[1];
+  mbinyhi = static_cast<int> ((coord-bboxlo[1])*bininvy);
 
-  coord = bboxlo[2] - cutneigh - SMALL*bbox[2];
-  mbinzlo = static_cast<int> ((coord-boxlo[2])*bininvz);
-  if (coord < 0.0) mbinzlo = mbinzlo - 1;
-  coord = bboxhi[2] + cutneigh + SMALL*bbox[2];
-  mbinzhi = static_cast<int> ((coord-boxlo[2])*bininvz);
+  coord = bsubboxlo[2] - SMALL*bbox[2];
+  mbinzlo = static_cast<int> ((coord-bboxlo[2])*bininvz);
+  if (coord < bboxlo[2]) mbinzlo = mbinzlo - 1;
+  coord = bsubboxhi[2] + SMALL*bbox[2];
+  mbinzhi = static_cast<int> ((coord-bboxlo[2])*bininvz);
 
   // extend bins by 1 to insure stencil extent is included
 
@@ -1415,37 +1432,37 @@ void Neighbor::bin_atoms()
 
 /* ----------------------------------------------------------------------
    convert atom coords into local bin #
-   only ghost atoms will have coord >= boxhi or coord < boxlo
-   take special care to insure ghosts are put in correct bins
-   this is necessary so that both procs on either side of PBC
-     treat a pair of atoms straddling the PBC in a consistent way
-     triclinic is an exception, since stencil & neigh list built differently
+   for orthogonal, only ghost atoms will have coord >= bboxhi or coord < bboxlo
+     take special care to insure ghosts are put in correct bins
+     this is necessary so that both procs on either side of PBC
+       treat a pair of atoms straddling the PBC in a consistent way
+   for triclinic, doesn't matter since stencil & neigh list built differently
 ------------------------------------------------------------------------- */
 
 int Neighbor::coord2bin(double *x)
 {
   int ix,iy,iz;
 
-  if (x[0] >= boxhi[0])
-    ix = static_cast<int> ((x[0]-boxhi[0])*bininvx) + nbinx - mbinxlo;
-  else if (x[0] >= boxlo[0])
-    ix = static_cast<int> ((x[0]-boxlo[0])*bininvx) - mbinxlo;
+  if (x[0] >= bboxhi[0])
+    ix = static_cast<int> ((x[0]-bboxhi[0])*bininvx) + nbinx - mbinxlo;
+  else if (x[0] >= bboxlo[0])
+    ix = static_cast<int> ((x[0]-bboxlo[0])*bininvx) - mbinxlo;
   else
-    ix = static_cast<int> ((x[0]-boxlo[0])*bininvx) - mbinxlo - 1;
+    ix = static_cast<int> ((x[0]-bboxlo[0])*bininvx) - mbinxlo - 1;
   
-  if (x[1] >= boxhi[1])
-    iy = static_cast<int> ((x[1]-boxhi[1])*bininvy) + nbiny - mbinylo;
-  else if (x[1] >= boxlo[1])
-    iy = static_cast<int> ((x[1]-boxlo[1])*bininvy) - mbinylo;
+  if (x[1] >= bboxhi[1])
+    iy = static_cast<int> ((x[1]-bboxhi[1])*bininvy) + nbiny - mbinylo;
+  else if (x[1] >= bboxlo[1])
+    iy = static_cast<int> ((x[1]-bboxlo[1])*bininvy) - mbinylo;
   else
-    iy = static_cast<int> ((x[1]-boxlo[1])*bininvy) - mbinylo - 1;
+    iy = static_cast<int> ((x[1]-bboxlo[1])*bininvy) - mbinylo - 1;
   
-  if (x[2] >= boxhi[2])
-    iz = static_cast<int> ((x[2]-boxhi[2])*bininvz) + nbinz - mbinzlo;
-  else if (x[2] >= boxlo[2])
-    iz = static_cast<int> ((x[2]-boxlo[2])*bininvz) - mbinzlo;
+  if (x[2] >= bboxhi[2])
+    iz = static_cast<int> ((x[2]-bboxhi[2])*bininvz) + nbinz - mbinzlo;
+  else if (x[2] >= bboxlo[2])
+    iz = static_cast<int> ((x[2]-bboxlo[2])*bininvz) - mbinzlo;
   else
-    iz = static_cast<int> ((x[2]-boxlo[2])*bininvz) - mbinzlo - 1;
+    iz = static_cast<int> ((x[2]-bboxlo[2])*bininvz) - mbinzlo - 1;
 
   return (iz*mbiny*mbinx + iy*mbinx + ix + 1);
 }
diff --git a/src/neighbor.h b/src/neighbor.h
index a4354b49c0..73c9afcf3c 100644
--- a/src/neighbor.h
+++ b/src/neighbor.h
@@ -118,7 +118,7 @@ class Neighbor : protected Pointers {
   double bininvx,bininvy,bininvz;
 
   int triclinic;                   // 0 if domain is orthog, 1 if triclinic
-  double *boxlo,*boxhi;            // copies of domain bounds
+  double *bboxlo,*bboxhi;          // copy of full domain bounding box
 
   int **pages;                     // half neighbor list pages
   int maxpage;                     // # of half pages currently allocated