diff --git a/src/OPENMP/npair_multi_omp.cpp b/src/OPENMP/npair_multi_omp.cpp index bbdce96ffc..e5db08e157 100644 --- a/src/OPENMP/npair_multi_omp.cpp +++ b/src/OPENMP/npair_multi_omp.cpp @@ -118,15 +118,16 @@ void NPairMultiOmp::build(NeighList *list) if (icollection == jcollection) jbin = ibin; else jbin = coord2bin(x[i], jcollection); - // loop over all atoms in surrounding bins in stencil including self - // skip i = j - // use full stencil for all collection combinations - s = stencil_multi[icollection][jcollection]; ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { js = binhead_multi[jcollection][jbin + s[k]]; + + // own-bin for half stencil + if (HALF) + if (flag_half_multi[icollection][jcollection] && s[k] == 0) js = bins[i]; + for (j = js; j >= 0; j = bins[j]) { if (!HALF) { // Full neighbor list @@ -144,40 +145,52 @@ void NPairMultiOmp::build(NeighList *list) // below = lower z or (equal z and lower y) or (equal zy and lower x) // (equal zyx and j <= i) // latter excludes self-self interaction but allows superposed atoms - if (x[j][2] < ztmp) continue; - if (x[j][2] == ztmp) { - if (x[j][1] < ytmp) continue; - if (x[j][1] == ytmp) { - if (x[j][0] < xtmp) continue; - if (x[j][0] == xtmp && j <= i) continue; + + // if same size (same collection), use half stencil + if (flag_half_multi[icollection][jcollection]) { + if (x[j][2] < ztmp) continue; + if (x[j][2] == ztmp) { + if (x[j][1] < ytmp) continue; + if (x[j][1] == ytmp) { + if (x[j][0] < xtmp) continue; + if (x[j][0] == xtmp && j <= i) continue; + } } } } else { // Half neighbor list, newton on, orthonormal - // store every pair for every bin in stencil,except for i's bin + // if same size: uses half stencil so includes a check of the central bin + if (flag_half_multi[icollection][jcollection]){ + if (s[k] == 0) { + // if same collection, + // if j is owned atom, store it, since j is beyond i in linked list + // if j is ghost, only store if j coords are "above and to the right" of i - if (stencil[k] == 0) { - // if j is owned atom, store it, since j is beyond i in linked list - // if j is ghost, only store if j coords are "above and to the "right" of i - if (j >= nlocal) { - if (x[j][2] < ztmp) continue; - if (x[j][2] == ztmp) { - if (x[j][1] < ytmp) continue; - if (x[j][1] == ytmp && x[j][0] < xtmp) continue; + // if different collections, + // if j is owned atom, store it if j > i + // if j is ghost, only store if j coords are "above and to the right" of i + + if ((icollection != jcollection) && (j < i)) continue; + + if (j >= nlocal) { + if (x[j][2] < ztmp) continue; + if (x[j][2] == ztmp) { + if (x[j][1] < ytmp) continue; + if (x[j][1] == ytmp && x[j][0] < xtmp) continue; + } } } } } jtype = type[j]; - if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue; + if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue; delx = xtmp - x[j][0]; dely = ytmp - x[j][1]; delz = ztmp - x[j][2]; rsq = delx*delx + dely*dely + delz*delz; - if (SIZE) { radsum = radius[i] + radius[j]; cut = radsum + skin; diff --git a/src/npair.cpp b/src/npair.cpp index c1615411c0..6fddcfef54 100644 --- a/src/npair.cpp +++ b/src/npair.cpp @@ -174,6 +174,8 @@ void NPair::copy_stencil_info() nstencil_multi = ns->nstencil_multi; stencil_multi = ns->stencil_multi; + + flag_half_multi = ns->flag_half_multi; } /* ---------------------------------------------------------------------- diff --git a/src/npair.h b/src/npair.h index 3eeb1d48f4..661196e586 100644 --- a/src/npair.h +++ b/src/npair.h @@ -96,6 +96,7 @@ class NPair : protected Pointers { int *nstencil_multi_old; int **stencil_multi_old; double **distsq_multi_old; + bool **flag_half_multi; int **nstencil_multi; int ***stencil_multi; diff --git a/src/npair_multi.cpp b/src/npair_multi.cpp index 9f719f403a..76f6d90f6e 100644 --- a/src/npair_multi.cpp +++ b/src/npair_multi.cpp @@ -112,6 +112,11 @@ void NPairMulti::build(NeighList *list) for (k = 0; k < ns; k++) { js = binhead_multi[jcollection][jbin + s[k]]; + + // own-bin for half stencil + if (HALF) + if (flag_half_multi[icollection][jcollection] && s[k] == 0) js = bins[i]; + for (j = js; j >= 0; j = bins[j]) { if (!HALF) { // Full neighbor list @@ -129,26 +134,39 @@ void NPairMulti::build(NeighList *list) // below = lower z or (equal z and lower y) or (equal zy and lower x) // (equal zyx and j <= i) // latter excludes self-self interaction but allows superposed atoms - if (x[j][2] < ztmp) continue; - if (x[j][2] == ztmp) { - if (x[j][1] < ytmp) continue; - if (x[j][1] == ytmp) { - if (x[j][0] < xtmp) continue; - if (x[j][0] == xtmp && j <= i) continue; + + // if same size (same collection), use half stencil + if (flag_half_multi[icollection][jcollection]) { + if (x[j][2] < ztmp) continue; + if (x[j][2] == ztmp) { + if (x[j][1] < ytmp) continue; + if (x[j][1] == ytmp) { + if (x[j][0] < xtmp) continue; + if (x[j][0] == xtmp && j <= i) continue; + } } } } else { // Half neighbor list, newton on, orthonormal - // store every pair for every bin in stencil,except for i's bin + // if same size: uses half stencil so includes a check of the central bin + if (flag_half_multi[icollection][jcollection]) { + if (s[k] == 0) { + // if same collection, + // if j is owned atom, store it, since j is beyond i in linked list + // if j is ghost, only store if j coords are "above and to the right" of i - if (stencil[k] == 0) { - // if j is owned atom, store it, since j is beyond i in linked list - // if j is ghost, only store if j coords are "above and to the "right" of i - if (j >= nlocal) { - if (x[j][2] < ztmp) continue; - if (x[j][2] == ztmp) { - if (x[j][1] < ytmp) continue; - if (x[j][1] == ytmp && x[j][0] < xtmp) continue; + // if different collections, + // if j is owned atom, store it if j > i + // if j is ghost, only store if j coords are "above and to the right" of i + + if ((icollection != jcollection) && (j < i)) continue; + + if (j >= nlocal) { + if (x[j][2] < ztmp) continue; + if (x[j][2] == ztmp) { + if (x[j][1] < ytmp) continue; + if (x[j][1] == ytmp && x[j][0] < xtmp) continue; + } } } } diff --git a/src/nstencil_multi.cpp b/src/nstencil_multi.cpp index 6e9fefea1e..202dace29d 100644 --- a/src/nstencil_multi.cpp +++ b/src/nstencil_multi.cpp @@ -87,13 +87,13 @@ void NStencilMulti::create() bin_collection = bin_collection_multi[icollection][jcollection]; cutsq = cutcollectionsq[icollection][jcollection]; - half_flag = flag_half_multi[icollection][jcollection]; // Half and ortho stencils include central bin first // This preserves the historical order of the neighbor list // as the old npair classes used to separately parse the central bin first - if (half_flag && (!TRI)) stencil[nstencil++] = 0; + if (HALF && (!TRI)) + if (half_flag) stencil_multi[icollection][jcollection][ns++] = 0; // For half stencils, only the upper plane is needed int sy_min = sy; @@ -106,15 +106,16 @@ void NStencilMulti::create() for (k = -sz_min; k <= sz; k++) { for (j = -sy_min; j <= sy; j++) { for (i = -sx; i <= sx; i++) { - // Now only include "upper right" bins for half and ortho stencils - if (HALF) { - if (half_flag && (!DIM_3D) && (!TRI)) - if (! (j > 0 || (j == 0 && i > 0))) continue; - if (half_flag && DIM_3D && (!TRI)) - if (! (k > 0 || j > 0 || (j == 0 && i > 0))) continue; + if (HALF && (!TRI)) { + if (half_flag) { + if (DIM_3D) { + if (! (k > 0 || j > 0 || (j == 0 && i > 0))) continue; + } else { + if (! (j > 0 || (j == 0 && i > 0))) continue; + } + } } - if (bin_distance_multi(i,j,k,bin_collection) < cutsq) stencil_multi[icollection][jcollection][ns++] = k * mbiny * mbinx + j * mbinx + i; }