Simplifying stencils

src/nstencil.cpp

@@ -82,7 +82,6 @@ NStencil::NStencil(LAMMPS *lmp) : Pointers(lmp)
   stencil_half = nullptr;
   stencil_skip = nullptr;
   stencil_bin_type = nullptr;
-  stencil_cutsq = nullptr;
 
   dimension = domain->dimension;
 }
@@ -122,7 +121,6 @@ NStencil::~NStencil()
     memory->destroy(stencil_half);
     memory->destroy(stencil_skip);
     memory->destroy(stencil_bin_type);
-    memory->destroy(stencil_cutsq);
     
     memory->destroy(stencil_sx_multi2);
     memory->destroy(stencil_sy_multi2);
@@ -280,8 +278,6 @@ void NStencil::create_setup()
                    "neighstencil:stencil_skip");
     memory->create(stencil_bin_type, n+1, n+1, 
                    "neighstencil:stencil_bin_type");
-    memory->create(stencil_cutsq, n+1, n+1, 
-                   "neighstencil:stencil_cutsq");
 
     memory->create(stencil_sx_multi2, n+1, n+1, 
                    "neighstencil:stencil_sx_multi2");               
@@ -346,7 +342,7 @@ void NStencil::create_setup()
         stencil_mbiny_multi2[i][j] = mbiny_multi2[bin_type];
         stencil_mbinz_multi2[i][j] = mbinz_multi2[bin_type];
         
-        stencil_range = sqrt(stencil_cutsq[i][j]);
+        stencil_range = sqrt(cutneighsq[i][j]);
         
         sx = static_cast<int> (stencil_range*bininvx_multi2[bin_type]);
         if (sx*binsizex < stencil_range) sx++;

src/nstencil.h

@@ -49,7 +49,6 @@ class NStencil : protected Pointers {
   bool **stencil_half;             // flag creation of a half stencil for itype-jtype
   bool **stencil_skip;             // skip creation of itype-jtype stencils (for newton on)  
   int **stencil_bin_type;          // what type to use for bin information
-  double **stencil_cutsq;          // cutoff used for stencil size
 
   NStencil(class LAMMPS *);
   virtual ~NStencil();

src/nstencil_full_multi2_2d.cpp

@@ -32,32 +32,13 @@ void NStencilFullMulti22d::set_stencil_properties()
   int n = atom->ntypes;
   int i, j;
   
-  // like -> like => use half stencil
+  // Always look up neighbor using full stencil and neighbor's bin
-  for (i = 1; i <= n; i++) {
-    stencil_half[i][i] = 0;
-    stencil_skip[i][i] = 0;
-    stencil_bin_type[i][i] = i;
-    stencil_cutsq[i][i] = cutneighsq[i][i];
-  }
-
-  // smaller -> larger => use existing newtoff stencil in larger bin
-  // larger -> smaller => use multi-like stencil for small-large in smaller bin
-  // If types are same cutoff, use existing like-like stencil.
 
   for (i = 1; i <= n; i++) {
     for (j = 1; j <= n; j++) {
-      if(i == j) continue;
-
       stencil_half[i][j] = 0;
       stencil_skip[i][j] = 0;
-      
-      if(cutneighsq[i][i] <= cutneighsq[j][j]){
-        stencil_cutsq[i][j] = cutneighsq[j][j];   
       stencil_bin_type[i][j] = j;
-      } else {
-        stencil_cutsq[i][j] = cutneighsq[i][j];   
-        stencil_bin_type[i][j] = j;
-      }
     }
   }
 }
@@ -87,7 +68,7 @@ void NStencilFullMulti22d::create()
       
       bin_type = stencil_bin_type[itype][jtype];
       
-      cutsq = stencil_cutsq[itype][jtype];
+      cutsq = cutneighsq[itype][jtype];
       
       for (j = -sy; j <= sy; j++)
         for (i = -sx; i <= sx; i++)

src/nstencil_full_multi2_3d.cpp

@@ -32,32 +32,14 @@ void NStencilFullMulti23d::set_stencil_properties()
   int n = atom->ntypes;
   int i, j;
   
-  // like -> like => use half stencil
+  // Always look up neighbor using full stencil and neighbor's bin
-  for (i = 1; i <= n; i++) {
+  // Stencil cutoff set by i-j cutoff
-    stencil_half[i][i] = 0;
-    stencil_skip[i][i] = 0;
-    stencil_bin_type[i][i] = i;
-    stencil_cutsq[i][i] = cutneighsq[i][i];
-  }
-
-  // smaller -> larger => use existing newtoff stencil in larger bin
-  // larger -> smaller => use multi-like stencil for small-large in smaller bin
-  // If types are same cutoff, use existing like-like stencil.
 
   for (i = 1; i <= n; i++) {
     for (j = 1; j <= n; j++) {
-      if(i == j) continue;
-
       stencil_half[i][j] = 0;
       stencil_skip[i][j] = 0;
-      
-      if(cutneighsq[i][i] <= cutneighsq[j][j]){
-        stencil_cutsq[i][j] = cutneighsq[j][j];   
       stencil_bin_type[i][j] = j;
-      } else {
-        stencil_cutsq[i][j] = cutneighsq[i][j];   
-        stencil_bin_type[i][j] = j;
-      }
     }
   }
 }
@@ -89,7 +71,7 @@ void NStencilFullMulti23d::create()
       
       bin_type = stencil_bin_type[itype][jtype];
       
-      cutsq = stencil_cutsq[itype][jtype];
+      cutsq = cutneighsq[itype][jtype];
       
       for (k = -sz; k <= sz; k++)
         for (j = -sy; j <= sy; j++)

src/nstencil_half_multi2_2d.cpp

@@ -33,26 +33,16 @@ void NStencilHalfMulti22d::set_stencil_properties()
   int n = atom->ntypes;
   int i, j;
 
-  // like -> like => use half stencil
-  for (i = 1; i <= n; i++) {
-    stencil_half[i][i] = 1;
-    stencil_skip[i][i] = 0;
-    stencil_bin_type[i][i] = i;
-    stencil_cutsq[i][i] = cutneighsq[i][i];
-  }
-
   // Cross types: use full stencil, looking one way through hierarchy
   // smaller -> larger => use full stencil in larger bin
-  // larger -> smaller => no nstecil required
+  // larger -> smaller => no nstencil required
-  // If cut offs are same, use existing type-type stencil
+  // If cut offs are same, use half stencil
 
   for (i = 1; i <= n; i++) {
     for (j = 1; j <= n; j++) {
-      if(i == j) continue;
       if(cutneighsq[i][i] > cutneighsq[j][j]) continue;
 
       stencil_skip[i][j] = 0;
-      stencil_cutsq[i][j] = cutneighsq[i][j];          
       
       if(cutneighsq[i][i] == cutneighsq[j][j]){
         stencil_half[i][j] = 1;
@@ -90,7 +80,7 @@ void NStencilHalfMulti22d::create()
       
       bin_type = stencil_bin_type[itype][jtype];      
       
-      cutsq = stencil_cutsq[itype][jtype];
+      cutsq = cutneighsq[itype][jtype];
       
       if (stencil_half[itype][jtype]) {
         for (j = 0; j <= sy; j++)

src/nstencil_half_multi2_2d_tri.cpp

@@ -33,26 +33,16 @@ void NStencilHalfMulti22dTri::set_stencil_properties()
   int n = atom->ntypes;
   int i, j;
 
-  // like -> like => use half stencil
-  for (i = 1; i <= n; i++) {
-    stencil_half[i][i] = 1;
-    stencil_skip[i][i] = 0;
-    stencil_bin_type[i][i] = i;
-    stencil_cutsq[i][i] = cutneighsq[i][i];
-  }
-
   // Cross types: use full stencil, looking one way through hierarchy
   // smaller -> larger => use full stencil in larger bin
-  // larger -> smaller => no nstecil required
+  // larger -> smaller => no nstencil required
-  // If cut offs are same, use existing type-type stencil
+  // If cut offs are same, use half stencil
 
   for (i = 1; i <= n; i++) {
     for (j = 1; j <= n; j++) {
-      if(i == j) continue;
       if(cutneighsq[i][i] > cutneighsq[j][j]) continue;
 
       stencil_skip[i][j] = 0;
-      stencil_cutsq[i][j] = cutneighsq[i][j];          
       
       if(cutneighsq[i][i] == cutneighsq[j][j]){
         stencil_half[i][j] = 1;
@@ -90,7 +80,7 @@ void NStencilHalfMulti22dTri::create()
       
       bin_type = stencil_bin_type[itype][jtype];
       
-      cutsq = stencil_cutsq[itype][jtype];
+      cutsq = cutneighsq[itype][jtype];
       
       if (stencil_half[itype][jtype]) {
         for (j = 0; j <= sy; j++)

src/nstencil_half_multi2_3d.cpp

@@ -33,26 +33,16 @@ void NStencilHalfMulti23d::set_stencil_properties()
   int n = atom->ntypes;
   int i, j;
 
-  // like -> like => use half stencil
-  for (i = 1; i <= n; i++) {
-    stencil_half[i][i] = 1;
-    stencil_skip[i][i] = 0;
-    stencil_bin_type[i][i] = i;
-    stencil_cutsq[i][i] = cutneighsq[i][i];
-  }
-
   // Cross types: use full stencil, looking one way through hierarchy
   // smaller -> larger => use full stencil in larger bin
-  // larger -> smaller => no nstecil required
+  // larger -> smaller => no nstencil required
-  // If cut offs are same, use existing type-type stencil
+  // If cut offs are same, use half stencil
 
   for (i = 1; i <= n; i++) {
     for (j = 1; j <= n; j++) {
-      if(i == j) continue;
       if(cutneighsq[i][i] > cutneighsq[j][j]) continue;
 
       stencil_skip[i][j] = 0;
-      stencil_cutsq[i][j] = cutneighsq[i][j];          
       
       if(cutneighsq[i][i] == cutneighsq[j][j]){
         stencil_half[i][j] = 1;
@@ -92,7 +82,7 @@ void NStencilHalfMulti23d::create()
       
       bin_type = stencil_bin_type[itype][jtype];
       
-      cutsq = stencil_cutsq[itype][jtype];
+      cutsq = cutneighsq[itype][jtype];
       
       if (stencil_half[itype][jtype]) {
         for (k = 0; k <= sz; k++)

src/nstencil_half_multi2_3d_tri.cpp

@@ -33,26 +33,16 @@ void NStencilHalfMulti23dTri::set_stencil_properties()
   int n = atom->ntypes;
   int i, j;
   
-  // like -> like => use half stencil
-  for (i = 1; i <= n; i++) {
-    stencil_half[i][i] = 1;
-    stencil_skip[i][i] = 0;
-    stencil_bin_type[i][i] = i;
-    stencil_cutsq[i][i] = cutneighsq[i][i];
-  }
-
   // Cross types: use full stencil, looking one way through hierarchy
   // smaller -> larger => use full stencil in larger bin
-  // larger -> smaller => no nstecil required
+  // larger -> smaller => no nstencil required
-  // If cut offs are same, use existing type-type stencil
+  // If cut offs are same, use half stencil
 
   for (i = 1; i <= n; i++) {
     for (j = 1; j <= n; j++) {
-      if(i == j) continue;
       if(cutneighsq[i][i] > cutneighsq[j][j]) continue;
 
       stencil_skip[i][j] = 0;
-      stencil_cutsq[i][j] = cutneighsq[i][j];          
       
       if(cutneighsq[i][i] == cutneighsq[j][j]){
         stencil_half[i][j] = 1;
@@ -92,7 +82,7 @@ void NStencilHalfMulti23dTri::create()
       
       bin_type = stencil_bin_type[itype][jtype];
       
-      cutsq = stencil_cutsq[itype][jtype];
+      cutsq = cutneighsq[itype][jtype];
       
       if (stencil_half[itype][jtype]) {
         for (k = 0; k <= sz; k++)