diff --git a/src/GRANULAR/pair_gran_hooke_history.cpp b/src/GRANULAR/pair_gran_hooke_history.cpp index 7b7586d355..46b41126aa 100644 --- a/src/GRANULAR/pair_gran_hooke_history.cpp +++ b/src/GRANULAR/pair_gran_hooke_history.cpp @@ -43,6 +43,7 @@ PairGranHookeHistory::PairGranHookeHistory(LAMMPS *lmp) : Pair(lmp) single_enable = 1; no_virial_fdotr_compute = 1; centroidstressflag = CENTROID_NOTAVAIL; + finitecutflag = 1; history = 1; size_history = 3; @@ -796,3 +797,23 @@ double PairGranHookeHistory::memory_usage() double bytes = nmax * sizeof(double); return bytes; } + +/* ---------------------------------------------------------------------- + self-interaction range of particle +------------------------------------------------------------------------- */ + +double PairGranHookeHistory::atom2cut(int i) +{ + double cut = atom->radius[i]*2; + return cut; +} + +/* ---------------------------------------------------------------------- + maximum interaction range for two finite particles +------------------------------------------------------------------------- */ + +double PairGranHookeHistory::radii2cut(double r1, double r2) +{ + double cut = r1+r2; + return cut; +} \ No newline at end of file diff --git a/src/GRANULAR/pair_gran_hooke_history.h b/src/GRANULAR/pair_gran_hooke_history.h index c27ce8a9af..9575c49d0c 100644 --- a/src/GRANULAR/pair_gran_hooke_history.h +++ b/src/GRANULAR/pair_gran_hooke_history.h @@ -42,6 +42,8 @@ class PairGranHookeHistory : public Pair { int pack_forward_comm(int, int *, double *, int, int *); void unpack_forward_comm(int, int, double *); double memory_usage(); + double atom2cut(int); + double radii2cut(double,double); protected: double kn,kt,gamman,gammat,xmu; diff --git a/src/GRANULAR/pair_granular.cpp b/src/GRANULAR/pair_granular.cpp index 215926e23e..0f1fbd49f3 100644 --- a/src/GRANULAR/pair_granular.cpp +++ b/src/GRANULAR/pair_granular.cpp @@ -67,6 +67,7 @@ PairGranular::PairGranular(LAMMPS *lmp) : Pair(lmp) single_enable = 1; no_virial_fdotr_compute = 1; centroidstressflag = CENTROID_NOTAVAIL; + finitecutflag = 1; single_extra = 12; svector = new double[single_extra]; @@ -1819,3 +1820,50 @@ void PairGranular::transfer_history(double* source, double* target) for (int i = 0; i < size_history; i++) target[i] = history_transfer_factors[i]*source[i]; } + +/* ---------------------------------------------------------------------- + self-interaction range of particle +------------------------------------------------------------------------- */ + +double PairGranular::atom2cut(int i) +{ + double cut; + + cut = atom->radius[i]*2; + if(beyond_contact) { + int itype = atom->type[i] + if(normal_model[itype][itype] == JKR) { + cut += pulloffdistance(cut, cut, itype, itype); + } + } + + return cut; +} + +/* ---------------------------------------------------------------------- + maximum interaction range for two finite particles +------------------------------------------------------------------------- */ + +double PairGranular::radii2cut(double r1, double r2) +{ + double cut = 0.0; + + if(beyond_contact) { + int n = atom->ntypes; + double temp; + + // Check all combinations of i and j to find theoretical maximum pull off distance + for(int i = 0; i < n; i++){ + for(int j = 0; j < n; j++){ + if(normal_model[i][j] == JKR) { + temp = pulloffdistance(r1, r2, i, j); + if(temp > cut) cut = temp; + } + } + } + } + + cut += r1 + r2; + + return cut; +} diff --git a/src/GRANULAR/pair_granular.h b/src/GRANULAR/pair_granular.h index 7ef4f2af98..7fe47b3275 100644 --- a/src/GRANULAR/pair_granular.h +++ b/src/GRANULAR/pair_granular.h @@ -40,6 +40,8 @@ class PairGranular : public Pair { int pack_forward_comm(int, int *, double *, int, int *); void unpack_forward_comm(int, int, double *); double memory_usage(); + double atom2cut(int); + double radii2cut(double,double); protected: double dt; diff --git a/src/USER-DPD/nstencil_half_bin_2d_newton_ssa.cpp b/src/USER-DPD/nstencil_half_bin_2d_ssa.cpp similarity index 95% rename from src/USER-DPD/nstencil_half_bin_2d_newton_ssa.cpp rename to src/USER-DPD/nstencil_half_bin_2d_ssa.cpp index 803a80f3e3..aca7b9df09 100644 --- a/src/USER-DPD/nstencil_half_bin_2d_newton_ssa.cpp +++ b/src/USER-DPD/nstencil_half_bin_2d_ssa.cpp @@ -16,13 +16,13 @@ James Larentzos and Timothy I. Mattox (Engility Corporation) ------------------------------------------------------------------------- */ -#include "nstencil_half_bin_2d_newton_ssa.h" +#include "nstencil_half_bin_2d_ssa.h" using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -NStencilHalfBin2dNewtonSSA::NStencilHalfBin2dNewtonSSA(LAMMPS *lmp) : +NStencilHalfBin2dSSA::NStencilHalfBin2dSSA(LAMMPS *lmp) : NStencilSSA(lmp) {} /* ---------------------------------------------------------------------- @@ -37,7 +37,7 @@ NStencilHalfBin2dNewtonSSA::NStencilHalfBin2dNewtonSSA(LAMMPS *lmp) : to locate all the Active Interaction Region (AIR) ghosts for SSA ------------------------------------------------------------------------- */ -void NStencilHalfBin2dNewtonSSA::create() +void NStencilHalfBin2dSSA::create() { int i,j,pos = 0; nstencil_ssa[0] = 0; // redundant info, but saves a conditional diff --git a/src/USER-DPD/nstencil_half_bin_2d_newton_ssa.h b/src/USER-DPD/nstencil_half_bin_2d_ssa.h similarity index 66% rename from src/USER-DPD/nstencil_half_bin_2d_newton_ssa.h rename to src/USER-DPD/nstencil_half_bin_2d_ssa.h index 1d5cc3f6b2..a50c852670 100644 --- a/src/USER-DPD/nstencil_half_bin_2d_newton_ssa.h +++ b/src/USER-DPD/nstencil_half_bin_2d_ssa.h @@ -13,23 +13,23 @@ #ifdef NSTENCIL_CLASS -NStencilStyle(half/bin/2d/newton/ssa, - NStencilHalfBin2dNewtonSSA, - NS_HALF | NS_BIN | NS_2D | NS_NEWTON | NS_SSA | NS_ORTHO | NS_GHOST) +NStencilStyle(half/bin/2d/ssa, + NStencilHalfBin2dSSA, + NS_HALF | NS_BIN | NS_2D | NS_SSA | NS_ORTHO | NS_GHOST) #else -#ifndef LMP_NSTENCIL_HALF_BIN_2D_NEWTON_SSA_H -#define LMP_NSTENCIL_HALF_BIN_2D_NEWTON_SSA_H +#ifndef LMP_NSTENCIL_HALF_BIN_2D_SSA_H +#define LMP_NSTENCIL_HALF_BIN_2D_SSA_H #include "nstencil_ssa.h" namespace LAMMPS_NS { -class NStencilHalfBin2dNewtonSSA : public NStencilSSA { +class NStencilHalfBin2dSSA : public NStencilSSA { public: - NStencilHalfBin2dNewtonSSA(class LAMMPS *); - ~NStencilHalfBin2dNewtonSSA() {} + NStencilHalfBin2dSSA(class LAMMPS *); + ~NStencilHalfBin2dSSA() {} void create(); }; diff --git a/src/USER-DPD/nstencil_half_bin_3d_newton_ssa.cpp b/src/USER-DPD/nstencil_half_bin_3d_ssa.cpp similarity index 97% rename from src/USER-DPD/nstencil_half_bin_3d_newton_ssa.cpp rename to src/USER-DPD/nstencil_half_bin_3d_ssa.cpp index c20dabe28f..8e33f5414f 100644 --- a/src/USER-DPD/nstencil_half_bin_3d_newton_ssa.cpp +++ b/src/USER-DPD/nstencil_half_bin_3d_ssa.cpp @@ -16,13 +16,13 @@ James Larentzos and Timothy I. Mattox (Engility Corporation) ------------------------------------------------------------------------- */ -#include "nstencil_half_bin_3d_newton_ssa.h" +#include "nstencil_half_bin_3d_ssa.h" using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -NStencilHalfBin3dNewtonSSA::NStencilHalfBin3dNewtonSSA(LAMMPS *lmp) : +NStencilHalfBin3dSSA::NStencilHalfBin3dSSA(LAMMPS *lmp) : NStencilSSA(lmp) {} /* ---------------------------------------------------------------------- @@ -37,7 +37,7 @@ NStencilHalfBin3dNewtonSSA::NStencilHalfBin3dNewtonSSA(LAMMPS *lmp) : to locate all the Active Interaction Region (AIR) ghosts for SSA ------------------------------------------------------------------------- */ -void NStencilHalfBin3dNewtonSSA::create() +void NStencilHalfBin3dSSA::create() { int i,j,k,pos = 0; nstencil_ssa[0] = 0; // redundant info, but saves a conditional diff --git a/src/USER-DPD/nstencil_half_bin_3d_newton_ssa.h b/src/USER-DPD/nstencil_half_bin_3d_ssa.h similarity index 66% rename from src/USER-DPD/nstencil_half_bin_3d_newton_ssa.h rename to src/USER-DPD/nstencil_half_bin_3d_ssa.h index 450a696e46..94f8fb4142 100644 --- a/src/USER-DPD/nstencil_half_bin_3d_newton_ssa.h +++ b/src/USER-DPD/nstencil_half_bin_3d_ssa.h @@ -13,23 +13,23 @@ #ifdef NSTENCIL_CLASS -NStencilStyle(half/bin/3d/newton/ssa, - NStencilHalfBin3dNewtonSSA, - NS_HALF | NS_BIN | NS_3D | NS_NEWTON | NS_SSA | NS_ORTHO | NS_GHOST) +NStencilStyle(half/bin/3d/ssa, + NStencilHalfBin3dSSA, + NS_HALF | NS_BIN | NS_3D | NS_SSA | NS_ORTHO | NS_GHOST) #else -#ifndef LMP_NSTENCIL_HALF_BIN_3D_NEWTON_SSA_H -#define LMP_NSTENCIL_HALF_BIN_3D_NEWTON_SSA_H +#ifndef LMP_NSTENCIL_HALF_BIN_3D_SSA_H +#define LMP_NSTENCIL_HALF_BIN_3D_SSA_H #include "nstencil_ssa.h" namespace LAMMPS_NS { -class NStencilHalfBin3dNewtonSSA : public NStencilSSA { +class NStencilHalfBin3dSSA : public NStencilSSA { public: - NStencilHalfBin3dNewtonSSA(class LAMMPS *); - ~NStencilHalfBin3dNewtonSSA() {} + NStencilHalfBin3dSSA(class LAMMPS *); + ~NStencilHalfBin3dSSA() {} void create(); }; diff --git a/src/USER-OMP/npair_full_multi_omp.cpp b/src/USER-OMP/npair_full_multi_omp.cpp index 4b38e660d4..4fc6fb4a68 100755 --- a/src/USER-OMP/npair_full_multi_omp.cpp +++ b/src/USER-OMP/npair_full_multi_omp.cpp @@ -14,6 +14,7 @@ #include "omp_compat.h" #include "npair_full_multi_omp.h" #include "npair_omp.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -30,7 +31,7 @@ NPairFullMultiOmp::NPairFullMultiOmp(LAMMPS *lmp) : NPair(lmp) {} /* ---------------------------------------------------------------------- binned neighbor list construction for all neighbors - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent every neighbor pair appears in list of both atoms i and j ------------------------------------------------------------------------- */ @@ -46,7 +47,7 @@ void NPairFullMultiOmp::build(NeighList *list) #endif NPAIR_OMP_SETUP(nlocal); - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,which,ns,imol,iatom; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,which,ns,imol,iatom; tagint tagprev; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; int *neighptr,*s; @@ -54,6 +55,7 @@ void NPairFullMultiOmp::build(NeighList *list) // loop over each atom, storing neighbors + int *collection = neighbor->collection; double **x = atom->x; int *type = atom->type; int *mask = atom->mask; @@ -80,7 +82,7 @@ void NPairFullMultiOmp::build(NeighList *list) neighptr = ipage.vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -92,22 +94,22 @@ void NPairFullMultiOmp::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + 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 group combinations + // use full stencil for all collection combinations - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >= 0; j = bins[j]) { if (i == j) continue; diff --git a/src/USER-OMP/npair_half_multi_newtoff_omp.cpp b/src/USER-OMP/npair_half_multi_newtoff_omp.cpp index 93d150f445..ed0e770b2f 100755 --- a/src/USER-OMP/npair_half_multi_newtoff_omp.cpp +++ b/src/USER-OMP/npair_half_multi_newtoff_omp.cpp @@ -14,6 +14,7 @@ #include "omp_compat.h" #include "npair_half_multi_newtoff_omp.h" #include "npair_omp.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -30,7 +31,7 @@ NPairHalfMultiNewtoffOmp::NPairHalfMultiNewtoffOmp(LAMMPS *lmp) : NPair(lmp) {} /* ---------------------------------------------------------------------- binned neighbor list construction with partial Newton's 3rd law - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks own bin and other bins in stencil pair stored once if i,j are both owned and i < j pair stored by me if j is ghost (also stored by proc owning j) @@ -48,7 +49,7 @@ void NPairHalfMultiNewtoffOmp::build(NeighList *list) #endif NPAIR_OMP_SETUP(nlocal); - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,which,ns,imol,iatom; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,which,ns,imol,iatom; tagint tagprev; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; int *neighptr,*s; @@ -56,6 +57,7 @@ void NPairHalfMultiNewtoffOmp::build(NeighList *list) // loop over each atom, storing neighbors + int *collection = neighbor->collection; double **x = atom->x; int *type = atom->type; int *mask = atom->mask; @@ -82,7 +84,7 @@ void NPairHalfMultiNewtoffOmp::build(NeighList *list) neighptr = ipage.vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -94,24 +96,24 @@ void NPairHalfMultiNewtoffOmp::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // loop over all atoms in other bins in stencil including self // only store pair if i < j // stores own/own pairs only once // stores own/ghost pairs on both procs - // use full stencil for all group combinations + // use full stencil for all collection combinations - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >=0; j = bins[j]) { if (j <= i) continue; diff --git a/src/USER-OMP/npair_half_multi_newton_omp.cpp b/src/USER-OMP/npair_half_multi_newton_omp.cpp index 11a84a3040..13d286ab5b 100755 --- a/src/USER-OMP/npair_half_multi_newton_omp.cpp +++ b/src/USER-OMP/npair_half_multi_newton_omp.cpp @@ -14,6 +14,7 @@ #include "omp_compat.h" #include "npair_half_multi_newton_omp.h" #include "npair_omp.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -30,7 +31,7 @@ NPairHalfMultiNewtonOmp::NPairHalfMultiNewtonOmp(LAMMPS *lmp) : NPair(lmp) {} /* ---------------------------------------------------------------------- binned neighbor list construction with full Newton's 3rd law - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks its own bin and other bins in Newton stencil every pair stored exactly once by some processor ------------------------------------------------------------------------- */ @@ -47,7 +48,7 @@ void NPairHalfMultiNewtonOmp::build(NeighList *list) #endif NPAIR_OMP_SETUP(nlocal); - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,which,ns,imol,iatom; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,which,ns,imol,iatom; tagint tagprev; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; int *neighptr,*s; @@ -55,6 +56,7 @@ void NPairHalfMultiNewtonOmp::build(NeighList *list) // loop over each atom, storing neighbors + int *collection = neighbor->collection; double **x = atom->x; int *type = atom->type; int *mask = atom->mask; @@ -81,7 +83,7 @@ void NPairHalfMultiNewtonOmp::build(NeighList *list) neighptr = ipage.vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -93,29 +95,29 @@ void NPairHalfMultiNewtonOmp::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // if same size: uses half stencil so check central bin - if(cutmultisq[igroup][igroup] == cutmultisq[jgroup][jgroup]){ + if(cutcollectionsq[icollection][icollection] == cutcollectionsq[jcollection][jcollection]){ - if(igroup == jgroup) js = bins[i]; - else js = binhead_multi[jgroup][jbin]; + if(icollection == jcollection) js = bins[i]; + else js = binhead_multi[jcollection][jbin]; - // if same group, + // 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 different groups, + // 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 for (j = js; j >= 0; j = bins[j]) { - if(igroup != jgroup and j < i) continue; + if(icollection != jcollection and j < i) continue; if (j >= nlocal) { if (x[j][2] < ztmp) continue; @@ -151,16 +153,16 @@ void NPairHalfMultiNewtonOmp::build(NeighList *list) } } - // for all groups, loop over all atoms in other bins in stencil, store every pair + // for all collections, loop over all atoms in other bins in stencil, store every pair // stencil is empty if i larger than j // stencil is half if i same size as j // stencil is full if i smaller than j - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >= 0; j = bins[j]) { jtype = type[j]; diff --git a/src/USER-OMP/npair_half_multi_newton_tri_omp.cpp b/src/USER-OMP/npair_half_multi_newton_tri_omp.cpp index b5882580a8..cb758673d4 100755 --- a/src/USER-OMP/npair_half_multi_newton_tri_omp.cpp +++ b/src/USER-OMP/npair_half_multi_newton_tri_omp.cpp @@ -14,6 +14,7 @@ #include "omp_compat.h" #include "npair_half_multi_newton_tri_omp.h" #include "npair_omp.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -31,7 +32,7 @@ NPairHalfMultiNewtonTriOmp::NPairHalfMultiNewtonTriOmp(LAMMPS *lmp) : /* ---------------------------------------------------------------------- binned neighbor list construction with Newton's 3rd law for triclinic - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks its own bin and other bins in triclinic stencil every pair stored exactly once by some processor ------------------------------------------------------------------------- */ @@ -48,7 +49,7 @@ void NPairHalfMultiNewtonTriOmp::build(NeighList *list) #endif NPAIR_OMP_SETUP(nlocal); - int i,j,k,n,itype,jtype,ibin,jbin,igroup,jgroup,which,ns,imol,iatom; + int i,j,k,n,itype,jtype,ibin,jbin,icollection,jcollection,which,ns,imol,iatom; tagint tagprev; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; int *neighptr,*s; @@ -56,6 +57,7 @@ void NPairHalfMultiNewtonTriOmp::build(NeighList *list) // loop over each atom, storing neighbors + int *collection = neighbor->collection; double **x = atom->x; int *type = atom->type; int *mask = atom->mask; @@ -82,7 +84,7 @@ void NPairHalfMultiNewtonTriOmp::build(NeighList *list) neighptr = ipage.vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -94,12 +96,12 @@ void NPairHalfMultiNewtonTriOmp::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // loop over all atoms in bins in stencil // stencil is empty if i larger than j @@ -110,15 +112,15 @@ void NPairHalfMultiNewtonTriOmp::build(NeighList *list) // (equal zyx and j <= i) // latter excludes self-self interaction but allows superposed atoms - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >= 0; j = bins[j]) { - // if same size (same group), use half stencil - if(cutmultisq[igroup][igroup] == cutmultisq[jgroup][jgroup]){ + // if same size (same collection), use half stencil + if(cutcollectionsq[icollection][icollection] == cutcollectionsq[jcollection][jcollection]){ if (x[j][2] < ztmp) continue; if (x[j][2] == ztmp) { if (x[j][1] < ytmp) continue; diff --git a/src/USER-OMP/npair_half_size_multi_newtoff_omp.cpp b/src/USER-OMP/npair_half_size_multi_newtoff_omp.cpp index 0d588b2f24..0efb4f067a 100755 --- a/src/USER-OMP/npair_half_size_multi_newtoff_omp.cpp +++ b/src/USER-OMP/npair_half_size_multi_newtoff_omp.cpp @@ -14,6 +14,7 @@ #include "omp_compat.h" #include "npair_half_size_multi_newtoff_omp.h" #include "npair_omp.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -29,7 +30,7 @@ NPairHalfSizeMultiNewtoffOmp::NPairHalfSizeMultiNewtoffOmp(LAMMPS *lmp) : NPair( /* ---------------------------------------------------------------------- size particles binned neighbor list construction with partial Newton's 3rd law - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks own bin and other bins in stencil pair stored once if i,j are both owned and i < j pair stored by me if j is ghost (also stored by proc owning j) @@ -47,7 +48,7 @@ void NPairHalfSizeMultiNewtoffOmp::build(NeighList *list) #endif NPAIR_OMP_SETUP(nlocal); - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,ns; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,ns; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; double radi,radsum,cutdistsq; int *neighptr,*s; @@ -55,6 +56,7 @@ void NPairHalfSizeMultiNewtoffOmp::build(NeighList *list) // loop over each atom, storing neighbors + int *collection = neighbor->collection; double **x = atom->x; double *radius = atom->radius; int *type = atom->type; @@ -75,7 +77,7 @@ void NPairHalfSizeMultiNewtoffOmp::build(NeighList *list) neighptr = ipage.vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -83,24 +85,24 @@ void NPairHalfSizeMultiNewtoffOmp::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // loop over all atoms in other bins in stencil including self // only store pair if i < j // stores own/own pairs only once // stores own/ghost pairs on both procs - // use full stencil for all group combinations + // use full stencil for all collection combinations - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >=0; j = bins[j]) { if (j <= i) continue; diff --git a/src/USER-OMP/npair_half_size_multi_newton_omp.cpp b/src/USER-OMP/npair_half_size_multi_newton_omp.cpp index cec3edd02b..68885fb453 100755 --- a/src/USER-OMP/npair_half_size_multi_newton_omp.cpp +++ b/src/USER-OMP/npair_half_size_multi_newton_omp.cpp @@ -14,6 +14,7 @@ #include "omp_compat.h" #include "npair_half_size_multi_newton_omp.h" #include "npair_omp.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -29,7 +30,7 @@ NPairHalfSizeMultiNewtonOmp::NPairHalfSizeMultiNewtonOmp(LAMMPS *lmp) : NPair(lm /* ---------------------------------------------------------------------- size particles binned neighbor list construction with full Newton's 3rd law - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks its own bin and other bins in Newton stencil every pair stored exactly once by some processor ------------------------------------------------------------------------- */ @@ -46,7 +47,7 @@ void NPairHalfSizeMultiNewtonOmp::build(NeighList *list) #endif NPAIR_OMP_SETUP(nlocal); - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,ns; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,ns; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; double radi,radsum,cutdistsq; int *neighptr,*s; @@ -54,6 +55,7 @@ void NPairHalfSizeMultiNewtonOmp::build(NeighList *list) // loop over each atom, storing neighbors + int *collection = neighbor->collection; double **x = atom->x; double *radius = atom->radius; int *type = atom->type; @@ -74,7 +76,7 @@ void NPairHalfSizeMultiNewtonOmp::build(NeighList *list) neighptr = ipage.vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -82,29 +84,29 @@ void NPairHalfSizeMultiNewtonOmp::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // if same size: uses half stencil so check central bin - if(cutmultisq[igroup][igroup] == cutmultisq[jgroup][jgroup]){ + if(cutcollectionsq[icollection][icollection] == cutcollectionsq[jcollection][jcollection]){ - if(igroup == jgroup) js = bins[i]; - else js = binhead_multi[jgroup][jbin]; + if(icollection == jcollection) js = bins[i]; + else js = binhead_multi[jcollection][jbin]; - // if same group, + // 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 different groups, + // 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 for (j = js; j >= 0; j = bins[j]) { - if(igroup != jgroup and j < i) continue; + if(icollection != jcollection and j < i) continue; if (j >= nlocal) { if (x[j][2] < ztmp) continue; @@ -133,16 +135,16 @@ void NPairHalfSizeMultiNewtonOmp::build(NeighList *list) } } - // for all groups, loop over all atoms in other bins in stencil, store every pair + // for all collections, loop over all atoms in other bins in stencil, store every pair // stencil is empty if i larger than j // stencil is half if i same size as j // stencil is full if i smaller than j - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >= 0; j = bins[j]) { jtype = type[j]; diff --git a/src/USER-OMP/npair_half_size_multi_newton_tri_omp.cpp b/src/USER-OMP/npair_half_size_multi_newton_tri_omp.cpp index 396e870846..76f61b8792 100755 --- a/src/USER-OMP/npair_half_size_multi_newton_tri_omp.cpp +++ b/src/USER-OMP/npair_half_size_multi_newton_tri_omp.cpp @@ -14,6 +14,7 @@ #include "omp_compat.h" #include "npair_half_size_multi_newton_tri_omp.h" #include "npair_omp.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -30,7 +31,7 @@ NPairHalfSizeMultiNewtonTriOmp::NPairHalfSizeMultiNewtonTriOmp(LAMMPS *lmp) : /* ---------------------------------------------------------------------- size particles binned neighbor list construction with Newton's 3rd law for triclinic - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks its own bin and other bins in triclinic stencil every pair stored exactly once by some processor ------------------------------------------------------------------------- */ @@ -47,7 +48,7 @@ void NPairHalfSizeMultiNewtonTriOmp::build(NeighList *list) #endif NPAIR_OMP_SETUP(nlocal); - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,ns; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,ns; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; double radi,radsum,cutdistsq; int *neighptr,*s; @@ -55,6 +56,7 @@ void NPairHalfSizeMultiNewtonTriOmp::build(NeighList *list) // loop over each atom, storing neighbors + int *collection = neighbor->collection; double **x = atom->x; double *radius = atom->radius; int *type = atom->type; @@ -75,7 +77,7 @@ void NPairHalfSizeMultiNewtonTriOmp::build(NeighList *list) neighptr = ipage.vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -83,12 +85,12 @@ void NPairHalfSizeMultiNewtonTriOmp::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // loop over all atoms in bins in stencil @@ -100,15 +102,15 @@ void NPairHalfSizeMultiNewtonTriOmp::build(NeighList *list) // (equal zyx and j <= i) // latter excludes self-self interaction but allows superposed atoms - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >= 0; j = bins[j]) { - // if same size (same group), use half stencil - if(cutmultisq[igroup][igroup] == cutmultisq[jgroup][jgroup]){ + // if same size (same collection), use half stencil + if(cutcollectionsq[icollection][icollection] == cutcollectionsq[jcollection][jcollection]){ if (x[j][2] < ztmp) continue; if (x[j][2] == ztmp) { if (x[j][1] < ytmp) continue; diff --git a/src/comm.cpp b/src/comm.cpp index 867eab6256..3f2d4ecde6 100644 --- a/src/comm.cpp +++ b/src/comm.cpp @@ -58,6 +58,7 @@ Comm::Comm(LAMMPS *lmp) : Pointers(lmp) bordergroup = 0; cutghostuser = 0.0; cutusermulti = nullptr; + cutusermultiold = nullptr; ghost_velocity = 0; user_procgrid[0] = user_procgrid[1] = user_procgrid[2] = 0; @@ -118,6 +119,7 @@ Comm::~Comm() memory->destroy(ysplit); memory->destroy(zsplit); memory->destroy(cutusermulti); + memory->destroy(cutusermultiold); delete [] customfile; delete [] outfile; } @@ -146,8 +148,13 @@ void Comm::copy_arrays(Comm *oldcomm) } if (oldcomm->cutusermulti) { - memory->create(cutusermulti,atom->ntypes+1,"comm:cutusermulti"); - memcpy(cutusermulti,oldcomm->cutusermulti,atom->ntypes+1); + memory->create(cutusermulti,ncollections,"comm:cutusermulti"); + memcpy(cutusermulti,oldcomm->cutusermulti,ncollections); + } + + if (oldcomm->cutusermultiold) { + memory->create(cutusermultiold,atom->ntypes+1,"comm:cutusermultiold"); + memcpy(cutusermultiold,oldcomm->cutusermultiold,atom->ntypes+1); } if (customfile) { @@ -242,13 +249,16 @@ void Comm::init() for (int i = 0; i < nfix; i++) if (fix[i]->maxexchange_dynamic) maxexchange_fix_dynamic = 1; + if(mode == Comm::MULTI and neighbor->style != Neighbor::MULTI) + error->all(FLERR,"Cannot use comm mode multi without multi-style neighbor lists"); + if(multi_reduce){ if (force->newton == 0) error->all(FLERR,"Cannot use multi/reduce communication with Newton off"); if (neighbor->any_full()) error->all(FLERR,"Cannot use multi/reduce communication with a full neighbor list"); - if (neighbor->style != Neighbor::MULTI) - error->all(FLERR,"Cannot use multi/reduce communication without multi-style neighbor lists"); + if(mode != Comm::MULTI) + error->all(FLERR,"Cannot use multi/reduce communication without mode multi"); } } @@ -285,13 +295,20 @@ void Comm::modify_params(int narg, char **arg) if (strcmp(arg[iarg+1],"single") == 0) { // need to reset cutghostuser when switching comm mode if (mode == Comm::MULTI) cutghostuser = 0.0; + if (mode == Comm::MULTIOLD) cutghostuser = 0.0; memory->destroy(cutusermulti); cutusermulti = nullptr; mode = Comm::SINGLE; } else if (strcmp(arg[iarg+1],"multi") == 0) { // need to reset cutghostuser when switching comm mode if (mode == Comm::SINGLE) cutghostuser = 0.0; + if (mode == Comm::MULTIOLD) cutghostuser = 0.0; mode = Comm::MULTI; + } else if (strcmp(arg[iarg+1],"multi/old") == 0) { + // need to reset cutghostuser when switching comm mode + if (mode == Comm::SINGLE) cutghostuser = 0.0; + if (mode == Comm::MULTI) cutghostuser = 0.0; + mode = Comm::MULTIOLD; } else error->all(FLERR,"Illegal comm_modify command"); iarg += 2; } else if (strcmp(arg[iarg],"group") == 0) { @@ -308,6 +325,9 @@ void Comm::modify_params(int narg, char **arg) if (mode == Comm::MULTI) error->all(FLERR, "Use cutoff/multi keyword to set cutoff in multi mode"); + if (mode == Comm::MULTIOLD) + error->all(FLERR, + "Use cutoff/multi/old keyword to set cutoff in multi mode"); cutghostuser = utils::numeric(FLERR,arg[iarg+1],false,lmp); if (cutghostuser < 0.0) error->all(FLERR,"Invalid cutoff in comm_modify command"); @@ -317,18 +337,20 @@ void Comm::modify_params(int narg, char **arg) double cut; if (mode == Comm::SINGLE) error->all(FLERR,"Use cutoff keyword to set cutoff in single mode"); - if (domain->box_exist == 0) + if (mode == Comm::MULTIOLD) + error->all(FLERR,"Use cutoff/multi/old keyword to set cutoff in multi/old mode"); + if (domain->box_exist == 0) error->all(FLERR, "Cannot set cutoff/multi before simulation box is defined"); - const int ntypes = atom->ntypes; + ncollections = neighbor->ncollections; if (iarg+3 > narg) error->all(FLERR,"Illegal comm_modify command"); if (cutusermulti == nullptr) { - memory->create(cutusermulti,ntypes+1,"comm:cutusermulti"); - for (i=0; i < ntypes+1; ++i) + memory->create(cutusermulti,ncollections,"comm:cutusermulti"); + for (i=0; i < ncollections; ++i) cutusermulti[i] = -1.0; } - utils::bounds(FLERR,arg[iarg+1],1,ntypes,nlo,nhi,error); + utils::bounds(FLERR,arg[iarg+1],1,ncollections,nlo,nhi,error); cut = utils::numeric(FLERR,arg[iarg+2],false,lmp); cutghostuser = MAX(cutghostuser,cut); if (cut < 0.0) @@ -336,9 +358,35 @@ void Comm::modify_params(int narg, char **arg) for (i=nlo; i<=nhi; ++i) cutusermulti[i] = cut; iarg += 3; - } else if (strcmp(arg[iarg],"multi/reduce") == 0) { + } else if (strcmp(arg[iarg],"cutoff/multi/old") == 0) { + int i,nlo,nhi; + double cut; if (mode == Comm::SINGLE) - error->all(FLERR,"Use multi/reduce in mode multi only"); + error->all(FLERR,"Use cutoff keyword to set cutoff in single mode"); + if (mode == Comm::MULTI) + error->all(FLERR,"Use cutoff/multi keyword to set cutoff in multi mode"); + if (domain->box_exist == 0) + error->all(FLERR, + "Cannot set cutoff/multi before simulation box is defined"); + const int ntypes = atom->ntypes; + if (iarg+3 > narg) + error->all(FLERR,"Illegal comm_modify command"); + if (cutusermultiold == nullptr) { + memory->create(cutusermultiold,ntypes+1,"comm:cutusermultiold"); + for (i=0; i < ntypes+1; ++i) + cutusermultiold[i] = -1.0; + } + utils::bounds(FLERR,arg[iarg+1],1,ntypes,nlo,nhi,error); + cut = utils::numeric(FLERR,arg[iarg+2],false,lmp); + cutghostuser = MAX(cutghostuser,cut); + if (cut < 0.0) + error->all(FLERR,"Invalid cutoff in comm_modify command"); + for (i=nlo; i<=nhi; ++i) + cutusermultiold[i] = cut; + iarg += 3; + } else if (strcmp(arg[iarg],"reduce/multi") == 0) { + if (mode == Comm::SINGLE) + error->all(FLERR,"Use reduce/multi in mode multi only"); multi_reduce = 1; iarg += 1; } else if (strcmp(arg[iarg],"vel") == 0) { diff --git a/src/comm.h b/src/comm.h index 4ae4faf7fd..877941523c 100644 --- a/src/comm.h +++ b/src/comm.h @@ -25,14 +25,15 @@ class Comm : protected Pointers { // LAYOUT_NONUNIFORM = logical bricks, but diff sizes via LB // LAYOUT_TILED = general tiling, due to RCB LB enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; - int mode; // 0 = single cutoff, 1 = multi-type cutoff - enum{SINGLE,MULTI}; + int mode; // 0 = single cutoff, 1 = multi-collection cutoff, 2 = multiold-type cutoff + enum{SINGLE,MULTI,MULTIOLD}; int me,nprocs; // proc info int ghost_velocity; // 1 if ghost atoms have velocity, 0 if not double cutghost[3]; // cutoffs used for acquiring ghost atoms double cutghostuser; // user-specified ghost cutoff (mode == 0) - double *cutusermulti; // per type user ghost cutoff (mode == 1) + double *cutusermulti; // per collection user ghost cutoff (mode == 1) + double *cutusermultiold; // per type user ghost cutoff (mode == 2) int recv_from_partition; // recv proc layout from this partition int send_to_partition; // send my proc layout to this partition // -1 if no recv or send @@ -152,7 +153,9 @@ class Comm : protected Pointers { int ncores; // # of cores per node int coregrid[3]; // 3d grid of cores within a node int user_coregrid[3]; // user request for cores in each dim - int multi_reduce; // 1 if multi cutoff is intra-type cutoff + int multi_reduce; // 1 if multi cutoff is intra-collection cutoff + int ncollections; // number of collection cutoffs defined for multi + int ncollections_prior; // value of ncollections at last setup void init_exchange(); int rendezvous_irregular(int, char *, int, int, int *, diff --git a/src/comm_brick.cpp b/src/comm_brick.cpp index 939669fb75..ca0fea8757 100644 --- a/src/comm_brick.cpp +++ b/src/comm_brick.cpp @@ -46,6 +46,7 @@ CommBrick::CommBrick(LAMMPS *lmp) : size_reverse_send(nullptr), size_reverse_recv(nullptr), slablo(nullptr), slabhi(nullptr), multilo(nullptr), multihi(nullptr), cutghostmulti(nullptr), pbc_flag(nullptr), pbc(nullptr), firstrecv(nullptr), + multioldlo(nullptr), multioldhi(nullptr), cutghostmultiold(nullptr), sendlist(nullptr), localsendlist(nullptr), maxsendlist(nullptr), buf_send(nullptr), buf_recv(nullptr) { @@ -65,6 +66,11 @@ CommBrick::~CommBrick() memory->destroy(cutghostmulti); } + if (mode == Comm::MULTIOLD) { + free_multiold(); + memory->destroy(cutghostmultiold); + } + if (sendlist) for (int i = 0; i < maxswap; i++) memory->destroy(sendlist[i]); if (localsendlist) memory->destroy(localsendlist); memory->sfree(sendlist); @@ -101,6 +107,9 @@ void CommBrick::init_buffers() multilo = multihi = nullptr; cutghostmulti = nullptr; + multioldlo = multioldhi = nullptr; + cutghostmultiold = nullptr; + buf_send = buf_recv = nullptr; maxsend = maxrecv = BUFMIN; grow_send(maxsend,2); @@ -128,23 +137,37 @@ void CommBrick::init() init_exchange(); if (bufextra > bufextra_old) grow_send(maxsend+bufextra,2); - // memory for multi-style communication + // memory for multi style communication + // allocate in setup if (mode == Comm::MULTI && multilo == nullptr) { allocate_multi(maxswap); - memory->create(cutghostmulti,atom->ntypes+1,3,"comm:cutghostmulti"); - } - if (mode == Comm::SINGLE && multilo) { + memory->create(cutghostmulti,ncollections,3,"comm:cutghostmulti"); + ncollections_prior = ncollections; + } + if ((mode == Comm::SINGLE or mode == Comm::MULTIOLD) && multilo) { free_multi(); memory->destroy(cutghostmulti); } + + // memory for multi/old-style communication + + if (mode == Comm::MULTIOLD && multioldlo == nullptr) { + allocate_multiold(maxswap); + memory->create(cutghostmultiold,atom->ntypes+1,3,"comm:cutghostmultiold"); + } + if ((mode == Comm::SINGLE or mode == Comm::MULTI) && multioldlo) { + free_multiold(); + memory->destroy(cutghostmultiold); + } } /* ---------------------------------------------------------------------- setup spatial-decomposition communication patterns function of neighbor cutoff(s) & cutghostuser & current box size single mode sets slab boundaries (slablo,slabhi) based on max cutoff - multi mode sets type-dependent slab boundaries (multilo,multihi) + multi mode sets collection-dependent slab boundaries (multilo,multihi) + multi/old mode sets type-dependent slab boundaries (multioldlo,multioldhi) ------------------------------------------------------------------------- */ void CommBrick::setup() @@ -156,8 +179,10 @@ void CommBrick::setup() // neigh->cutghost = distance between tilted planes in box coords // cutghost is in lamda coords = distance between those planes // for multi: - // cutghostmulti = same as cutghost, only for each atom type - + // cutghostmulti = same as cutghost, only for each atom collection + // for multi/old: + // cutghostmultiold = same as cutghost, only for each atom type + int i,j; int ntypes = atom->ntypes; double *prd,*sublo,*subhi; @@ -167,49 +192,58 @@ void CommBrick::setup() error->warning(FLERR,"Communication cutoff is 0.0. No ghost atoms " "will be generated. Atoms may get lost."); - if (mode == Comm::MULTI) { - if (multi_reduce) { - // If using multi/reduce, communicate itype particles a distance equal - // to the max of itype-jtype group interaction - // only consider smaller jtype groups - int igroup, jgroup; - double **cutmultisq = neighbor->cutmultisq; - int *map_type_multi = neighbor->map_type_multi; - for (i = 1; i <= ntypes; i++) { - - if (cutusermulti) { - cutghostmulti[i][0] = cutusermulti[i]; - cutghostmulti[i][1] = cutusermulti[i]; - cutghostmulti[i][2] = cutusermulti[i]; - } else { - cutghostmulti[i][0] = 0.0; - cutghostmulti[i][1] = 0.0; - cutghostmulti[i][2] = 0.0; - } - - igroup = map_type_multi[i]; - for (j = 1; j <= ntypes; j++){ - jgroup = map_type_multi[j]; - - if(cutmultisq[jgroup][jgroup] > cutmultisq[igroup][igroup]) continue; - cutghostmulti[i][0] = MAX(cutghostmulti[i][0],sqrt(cutmultisq[igroup][jgroup])); - cutghostmulti[i][1] = MAX(cutghostmulti[i][1],sqrt(cutmultisq[igroup][jgroup])); - cutghostmulti[i][2] = MAX(cutghostmulti[i][2],sqrt(cutmultisq[igroup][jgroup])); - } + // build initial collection array + neighbor->build_collection(0); + + if(cutusermulti and ncollections != neighbor->ncollections) + error->all(FLERR, "Cannot change number of collections after defining comm_modify multi/cutoff"); + else ncollections = neighbor->ncollections; + + // reallocate memory for multi-style communication at setup if ncollections change + if(ncollections_prior != ncollections){ + if(multilo) free_multi(); + if(cutghostmulti) memory->destroy(cutghostmulti); + + allocate_multi(maxswap); + memory->create(cutghostmulti,ncollections,3,"comm:cutghostmulti"); + ncollections_prior = ncollections; + } + + double **cutcollectionsq = neighbor->cutcollectionsq; + // If using multi/reduce, communicate particles a distance equal + // to the max cutoff with equally sized or smaller collections + // If not, communicate the maximum cutoff of the entire collection + for (i = 0; i < ncollections; i++) { + if (cutusermulti) { + cutghostmulti[i][0] = cutusermulti[i]; + cutghostmulti[i][1] = cutusermulti[i]; + cutghostmulti[i][2] = cutusermulti[i]; + } else { + cutghostmulti[i][0] = 0.0; + cutghostmulti[i][1] = 0.0; + cutghostmulti[i][2] = 0.0; } - } else { - // otherwise, communicate a distance equal to the maximum interaction distance for each type - double *cuttype = neighbor->cuttype; - for (i = 1; i <= ntypes; i++) { - double tmp = 0.0; - if (cutusermulti) tmp = cutusermulti[i]; - cutghostmulti[i][0] = MAX(tmp,cuttype[i]); - cutghostmulti[i][1] = MAX(tmp,cuttype[i]); - cutghostmulti[i][2] = MAX(tmp,cuttype[i]); + + for (j = 0; j < ncollections; j++){ + if(multi_reduce and cutcollectionsq[j][j] > cutcollectionsq[i][i]) continue; + cutghostmulti[i][0] = MAX(cutghostmulti[i][0],sqrt(cutcollectionsq[i][j])); + cutghostmulti[i][1] = MAX(cutghostmulti[i][1],sqrt(cutcollectionsq[i][j])); + cutghostmulti[i][2] = MAX(cutghostmulti[i][2],sqrt(cutcollectionsq[i][j])); } } } + + if (mode == Comm::MULTIOLD) { + double *cuttype = neighbor->cuttype; + for (i = 1; i <= ntypes; i++) { + double tmp = 0.0; + if (cutusermultiold) tmp = cutusermultiold[i]; + cutghostmultiold[i][0] = MAX(tmp,cuttype[i]); + cutghostmultiold[i][1] = MAX(tmp,cuttype[i]); + cutghostmultiold[i][2] = MAX(tmp,cuttype[i]); + } + } if (triclinic == 0) { prd = domain->prd; @@ -228,13 +262,21 @@ void CommBrick::setup() cutghost[1] = cut * length1; length2 = h_inv[2]; cutghost[2] = cut * length2; - if (mode == Comm::MULTI){ - for (i = 1; i <= ntypes; i++) { + if (mode == Comm::MULTI) { + for (i = 0; i < ncollections; i++) { cutghostmulti[i][0] *= length0; cutghostmulti[i][1] *= length1; cutghostmulti[i][2] *= length2; } } + + if (mode == Comm::MULTIOLD) { + for (i = 1; i <= ntypes; i++) { + cutghostmultiold[i][0] *= length0; + cutghostmultiold[i][1] *= length1; + cutghostmultiold[i][2] *= length2; + } + } } // recvneed[idim][0/1] = # of procs away I recv atoms from, within cutghost @@ -378,12 +420,18 @@ void CommBrick::setup() if (ineed < 2) slablo[iswap] = -BIG; else slablo[iswap] = 0.5 * (sublo[dim] + subhi[dim]); slabhi[iswap] = sublo[dim] + cutghost[dim]; - } else { - for (i = 1; i <= ntypes; i++) { + } else if (mode == Comm::MULTI) { + for (i = 0; i < ncollections; i++) { if (ineed < 2) multilo[iswap][i] = -BIG; else multilo[iswap][i] = 0.5 * (sublo[dim] + subhi[dim]); multihi[iswap][i] = sublo[dim] + cutghostmulti[i][dim]; } + } else { + for (i = 1; i <= ntypes; i++) { + if (ineed < 2) multioldlo[iswap][i] = -BIG; + else multioldlo[iswap][i] = 0.5 * (sublo[dim] + subhi[dim]); + multioldhi[iswap][i] = sublo[dim] + cutghostmultiold[i][dim]; + } } if (myloc[dim] == 0) { pbc_flag[iswap] = 1; @@ -401,12 +449,18 @@ void CommBrick::setup() slablo[iswap] = subhi[dim] - cutghost[dim]; if (ineed < 2) slabhi[iswap] = BIG; else slabhi[iswap] = 0.5 * (sublo[dim] + subhi[dim]); - } else { - for (i = 1; i <= ntypes; i++) { + } else if (mode == Comm::MULTI) { + for (i = 0; i < ncollections; i++) { multilo[iswap][i] = subhi[dim] - cutghostmulti[i][dim]; if (ineed < 2) multihi[iswap][i] = BIG; else multihi[iswap][i] = 0.5 * (sublo[dim] + subhi[dim]); } + } else { + for (i = 1; i <= ntypes; i++) { + multioldlo[iswap][i] = subhi[dim] - cutghostmultiold[i][dim]; + if (ineed < 2) multioldhi[iswap][i] = BIG; + else multioldhi[iswap][i] = 0.5 * (sublo[dim] + subhi[dim]); + } } if (myloc[dim] == procgrid[dim]-1) { pbc_flag[iswap] = 1; @@ -727,15 +781,19 @@ void CommBrick::exchange() void CommBrick::borders() { - int i,n,itype,iswap,dim,ineed,twoneed; - int nsend,nrecv,sendflag,nfirst,nlast,ngroup; + int i,n,itype,icollection,iswap,dim,ineed,twoneed; + int nsend,nrecv,sendflag,nfirst,nlast,ngroup,nprior; double lo,hi; int *type; + int *collection; double **x; double *buf,*mlo,*mhi; MPI_Request request; AtomVec *avec = atom->avec; + // After exchanging/sorting, need to reconstruct collection array for border communication + if(mode == Comm::MULTI) neighbor->build_collection(0); + // do swaps over all 3 dimensions iswap = 0; @@ -756,10 +814,14 @@ void CommBrick::borders() if (mode == Comm::SINGLE) { lo = slablo[iswap]; hi = slabhi[iswap]; - } else { - type = atom->type; + } else if (mode == Comm::MULTI) { + collection = neighbor->collection; mlo = multilo[iswap]; mhi = multihi[iswap]; + } else { + type = atom->type; + mlo = multioldlo[iswap]; + mhi = multioldhi[iswap]; } if (ineed % 2 == 0) { nfirst = nlast; @@ -788,6 +850,14 @@ void CommBrick::borders() if (nsend == maxsendlist[iswap]) grow_list(iswap,nsend); sendlist[iswap][nsend++] = i; } + } else if (mode == Comm::MULTI) { + for (i = nfirst; i < nlast; i++) { + icollection = collection[i]; + if (x[i][dim] >= mlo[icollection] && x[i][dim] <= mhi[icollection]) { + if (nsend == maxsendlist[iswap]) grow_list(iswap,nsend); + sendlist[iswap][nsend++] = i; + } + } } else { for (i = nfirst; i < nlast; i++) { itype = type[i]; @@ -878,7 +948,10 @@ void CommBrick::borders() size_reverse_send[iswap] = nrecv*size_reverse; size_reverse_recv[iswap] = nsend*size_reverse; firstrecv[iswap] = atom->nlocal + atom->nghost; + nprior = atom->nlocal + atom->nghost; atom->nghost += nrecv; + if(neighbor->style == Neighbor::MULTI) neighbor->build_collection(nprior); + iswap++; } } @@ -1431,6 +1504,12 @@ void CommBrick::grow_swap(int n) allocate_multi(n); } + if (mode == Comm::MULTIOLD) { + free_multiold(); + allocate_multiold(n); + } + + sendlist = (int **) memory->srealloc(sendlist,n*sizeof(int *),"comm:sendlist"); memory->grow(maxsendlist,n,"comm:maxsendlist"); @@ -1462,15 +1541,26 @@ void CommBrick::allocate_swap(int n) } /* ---------------------------------------------------------------------- - allocation of multi-type swap info + allocation of multi-collection swap info ------------------------------------------------------------------------- */ void CommBrick::allocate_multi(int n) { - multilo = memory->create(multilo,n,atom->ntypes+1,"comm:multilo"); - multihi = memory->create(multihi,n,atom->ntypes+1,"comm:multihi"); + multilo = memory->create(multilo,n,ncollections,"comm:multilo"); + multihi = memory->create(multihi,n,ncollections,"comm:multihi"); } +/* ---------------------------------------------------------------------- + allocation of multi/old-type swap info +------------------------------------------------------------------------- */ + +void CommBrick::allocate_multiold(int n) +{ + multioldlo = memory->create(multioldlo,n,atom->ntypes+1,"comm:multioldlo"); + multioldhi = memory->create(multioldhi,n,atom->ntypes+1,"comm:multioldhi"); +} + + /* ---------------------------------------------------------------------- free memory for swaps ------------------------------------------------------------------------- */ @@ -1492,7 +1582,7 @@ void CommBrick::free_swap() } /* ---------------------------------------------------------------------- - free memory for multi-type swaps + free memory for multi-collection swaps ------------------------------------------------------------------------- */ void CommBrick::free_multi() @@ -1502,6 +1592,17 @@ void CommBrick::free_multi() multilo = multihi = nullptr; } +/* ---------------------------------------------------------------------- + free memory for multi/old-type swaps +------------------------------------------------------------------------- */ + +void CommBrick::free_multiold() +{ + memory->destroy(multioldlo); + memory->destroy(multioldhi); + multioldlo = multioldhi = nullptr; +} + /* ---------------------------------------------------------------------- extract data potentially useful to other classes ------------------------------------------------------------------------- */ diff --git a/src/comm_brick.h b/src/comm_brick.h index 6165f54de5..a1794fdb7b 100644 --- a/src/comm_brick.h +++ b/src/comm_brick.h @@ -61,8 +61,10 @@ class CommBrick : public Comm { int *size_reverse_send; // # to send in each reverse comm int *size_reverse_recv; // # to recv in each reverse comm double *slablo,*slabhi; // bounds of slab to send at each swap - double **multilo,**multihi; // bounds of slabs for multi-type swap - double **cutghostmulti; // cutghost on a per-type basis + double **multilo,**multihi; // bounds of slabs for multi-collection swap + double **multioldlo,**multioldhi; // bounds of slabs for multi-type swap + double **cutghostmulti; // cutghost on a per-collection basis + double **cutghostmultiold; // cutghost on a per-type basis int *pbc_flag; // general flag for sending atoms thru PBC int **pbc; // dimension flags for PBC adjustments @@ -84,11 +86,13 @@ class CommBrick : public Comm { virtual void grow_send(int, int); // reallocate send buffer virtual void grow_recv(int); // free/allocate recv buffer virtual void grow_list(int, int); // reallocate one sendlist - virtual void grow_swap(int); // grow swap and multi arrays + virtual void grow_swap(int); // grow swap, multi, and multi/old arrays virtual void allocate_swap(int); // allocate swap arrays virtual void allocate_multi(int); // allocate multi arrays + virtual void allocate_multiold(int); // allocate multi/old arrays virtual void free_swap(); // free swap arrays virtual void free_multi(); // free multi arrays + virtual void free_multiold(); // free multi/old arrays }; } diff --git a/src/comm_tiled.cpp b/src/comm_tiled.cpp index 2e916e7e1e..accdbf7b04 100644 --- a/src/comm_tiled.cpp +++ b/src/comm_tiled.cpp @@ -53,6 +53,7 @@ CommTiled::CommTiled(LAMMPS *lmp) : Comm(lmp) overlap = nullptr; rcbinfo = nullptr; cutghostmulti = nullptr; + cutghostmultiold = nullptr; init_buffers(); } @@ -81,6 +82,7 @@ CommTiled::~CommTiled() deallocate_swap(maxswap); memory->sfree(rcbinfo); memory->destroy(cutghostmulti); + memory->destroy(cutghostmultiold); } /* ---------------------------------------------------------------------- @@ -98,7 +100,9 @@ void CommTiled::init_buffers() overlap = nullptr; rcbinfo = nullptr; cutghostmulti = nullptr; + cutghostmultiold = nullptr; sendbox_multi = nullptr; + sendbox_multiold = nullptr; maxswap = 6; allocate_swap(maxswap); @@ -115,9 +119,10 @@ void CommTiled::init() nswap = 2*domain->dimension; - memory->destroy(cutghostmulti); - if (mode == Comm::MULTI) - memory->create(cutghostmulti,atom->ntypes+1,3,"comm:cutghostmulti"); + memory->destroy(cutghostmultiold); + if (mode == Comm::MULTIOLD) + memory->create(cutghostmultiold,atom->ntypes+1,3,"comm:cutghostmultiold"); + int bufextra_old = bufextra; init_exchange(); @@ -173,49 +178,62 @@ void CommTiled::setup() // set cutoff for comm forward and comm reverse // check that cutoff < any periodic box length - + if (mode == Comm::MULTI) { - if (multi_reduce) { - // If using multi/reduce, communicate itype particles a distance equal - // to the max of itype-jtype group interaction - // only consider smaller jtype groups - int igroup, jgroup; - double **cutmultisq = neighbor->cutmultisq; - int *map_type_multi = neighbor->map_type_multi; - for (i = 1; i <= ntypes; i++) { - - if (cutusermulti) { - cutghostmulti[i][0] = cutusermulti[i]; - cutghostmulti[i][1] = cutusermulti[i]; - cutghostmulti[i][2] = cutusermulti[i]; - } else { - cutghostmulti[i][0] = 0.0; - cutghostmulti[i][1] = 0.0; - cutghostmulti[i][2] = 0.0; - } + // build collection from scratch as it is needed for atom exchange + neighbor->build_collection(0); + + if(cutusermulti and ncollections != neighbor->ncollections) + error->all(FLERR, "Cannot change number of collections after defining comm_modify multi/cutoff"); + + ncollections = neighbor->ncollections; + + // allocate memory for multi-style communication at setup as ncollections can change + if(ncollections_prior != ncollections){ + memory->destroy(cutghostmulti); + if (mode == Comm::MULTI) + memory->create(cutghostmulti,ncollections,3,"comm:cutghostmulti"); - igroup = map_type_multi[i]; - for (j = 1; j <= ntypes; j++){ - jgroup = map_type_multi[j]; - - if(cutmultisq[jgroup][jgroup] > cutmultisq[igroup][igroup]) continue; - cutghostmulti[i][0] = MAX(cutghostmulti[i][0],sqrt(cutmultisq[igroup][jgroup])); - cutghostmulti[i][1] = MAX(cutghostmulti[i][1],sqrt(cutmultisq[igroup][jgroup])); - cutghostmulti[i][2] = MAX(cutghostmulti[i][2],sqrt(cutmultisq[igroup][jgroup])); - } + for(i = 0; i < maxswap; i ++) + grow_swap_send_multi(i,DELTA_PROCS); + + ncollections_prior = ncollections; + } + + double **cutcollectionsq = neighbor->cutcollectionsq; + // If using multi/reduce, communicate particles a distance equal + // to the max cutoff with equally sized or smaller collections + // If not, communicate the maximum cutoff of the entire collection + for (i = 0; i < ncollections; i++) { + if (cutusermulti) { + cutghostmulti[i][0] = cutusermulti[i]; + cutghostmulti[i][1] = cutusermulti[i]; + cutghostmulti[i][2] = cutusermulti[i]; + } else { + cutghostmulti[i][0] = 0.0; + cutghostmulti[i][1] = 0.0; + cutghostmulti[i][2] = 0.0; } - } else { - // otherwise, communicate a distance equal to the maximum interaction distance for each type - double *cuttype = neighbor->cuttype; - for (i = 1; i <= ntypes; i++) { - double tmp = 0.0; - if (cutusermulti) tmp = cutusermulti[i]; - cutghostmulti[i][0] = MAX(tmp,cuttype[i]); - cutghostmulti[i][1] = MAX(tmp,cuttype[i]); - cutghostmulti[i][2] = MAX(tmp,cuttype[i]); + + for (j = 0; j < ncollections; j++){ + if(multi_reduce and cutcollectionsq[j][j] > cutcollectionsq[i][i]) continue; + cutghostmulti[i][0] = MAX(cutghostmulti[i][0],sqrt(cutcollectionsq[i][j])); + cutghostmulti[i][1] = MAX(cutghostmulti[i][1],sqrt(cutcollectionsq[i][j])); + cutghostmulti[i][2] = MAX(cutghostmulti[i][2],sqrt(cutcollectionsq[i][j])); } } } + + if (mode == Comm::MULTIOLD) { + double *cuttype = neighbor->cuttype; + for (i = 1; i <= ntypes; i++) { + double tmp = 0.0; + if (cutusermultiold) tmp = cutusermultiold[i]; + cutghostmultiold[i][0] = MAX(tmp,cuttype[i]); + cutghostmultiold[i][1] = MAX(tmp,cuttype[i]); + cutghostmultiold[i][2] = MAX(tmp,cuttype[i]); + } + } double cut = get_comm_cutoff(); if ((cut == 0.0) && (me == 0)) @@ -233,12 +251,20 @@ void CommTiled::setup() length2 = h_inv[2]; cutghost[2] = cut * length2; if (mode == Comm::MULTI) { - for (i = 1; i <= ntypes; i++) { + for (i = 0; i < ncollections; i++) { cutghostmulti[i][0] *= length0; cutghostmulti[i][1] *= length1; cutghostmulti[i][2] *= length2; } } + + if (mode == Comm::MULTIOLD) { + for (i = 1; i <= ntypes; i++) { + cutghostmultiold[i][0] *= length0; + cutghostmultiold[i][1] *= length1; + cutghostmultiold[i][2] *= length2; + } + } } if ((periodicity[0] && cutghost[0] > prd[0]) || @@ -376,7 +402,7 @@ void CommTiled::setup() // extend sbox in those lower dims to include ghost atoms // single mode and multi mode - double oboxlo[3],oboxhi[3],sbox[6],sbox_multi[6]; + double oboxlo[3],oboxhi[3],sbox[6],sbox_multi[6],sbox_multiold[6]; if (mode == Comm::SINGLE) { for (i = 0; i < noverlap; i++) { @@ -465,7 +491,7 @@ void CommTiled::setup() sbox[5] = MIN(oboxhi[2],hi2[2]); } - for (int itype = 1; itype <= atom->ntypes; itype++) { + for (int icollection = 0; icollection < ncollections; icollection++) { sbox_multi[0] = sbox[0]; sbox_multi[1] = sbox[1]; sbox_multi[2] = sbox[2]; @@ -476,36 +502,112 @@ void CommTiled::setup() sbox_multi[idim] = sublo[idim]; if (i < noverlap1) sbox_multi[3+idim] = - MIN(sbox_multi[3+idim]+cutghostmulti[itype][idim],subhi[idim]); + MIN(sbox_multi[3+idim]+cutghostmulti[icollection][idim],subhi[idim]); else sbox_multi[3+idim] = - MIN(sbox_multi[3+idim]-prd[idim]+cutghostmulti[itype][idim], + MIN(sbox_multi[3+idim]-prd[idim]+cutghostmulti[icollection][idim], subhi[idim]); } else { if (i < noverlap1) sbox_multi[idim] = - MAX(sbox_multi[idim]-cutghostmulti[itype][idim],sublo[idim]); + MAX(sbox_multi[idim]-cutghostmulti[icollection][idim],sublo[idim]); else sbox_multi[idim] = - MAX(sbox_multi[idim]+prd[idim]-cutghostmulti[itype][idim], + MAX(sbox_multi[idim]+prd[idim]-cutghostmulti[icollection][idim], sublo[idim]); sbox_multi[3+idim] = subhi[idim]; } if (idim >= 1) { if (sbox_multi[0] == oboxlo[0]) - sbox_multi[0] -= cutghostmulti[itype][idim]; + sbox_multi[0] -= cutghostmulti[icollection][idim]; if (sbox_multi[3] == oboxhi[0]) - sbox_multi[3] += cutghostmulti[itype][idim]; + sbox_multi[3] += cutghostmulti[icollection][idim]; } if (idim == 2) { if (sbox_multi[1] == oboxlo[1]) - sbox_multi[1] -= cutghostmulti[itype][idim]; + sbox_multi[1] -= cutghostmulti[icollection][idim]; if (sbox_multi[4] == oboxhi[1]) - sbox_multi[4] += cutghostmulti[itype][idim]; + sbox_multi[4] += cutghostmulti[icollection][idim]; } - memcpy(sendbox_multi[iswap][i][itype],sbox_multi,6*sizeof(double)); + memcpy(sendbox_multi[iswap][i][icollection],sbox_multi,6*sizeof(double)); + } + } + } + + if (mode == Comm::MULTIOLD) { + for (i = 0; i < noverlap; i++) { + pbc_flag[iswap][i] = 0; + pbc[iswap][i][0] = pbc[iswap][i][1] = pbc[iswap][i][2] = + pbc[iswap][i][3] = pbc[iswap][i][4] = pbc[iswap][i][5] = 0; + + (this->*box_other)(idim,idir,overlap[i],oboxlo,oboxhi); + + if (i < noverlap1) { + sbox[0] = MAX(oboxlo[0],lo1[0]); + sbox[1] = MAX(oboxlo[1],lo1[1]); + sbox[2] = MAX(oboxlo[2],lo1[2]); + sbox[3] = MIN(oboxhi[0],hi1[0]); + sbox[4] = MIN(oboxhi[1],hi1[1]); + sbox[5] = MIN(oboxhi[2],hi1[2]); + } else { + pbc_flag[iswap][i] = 1; + if (idir == 0) pbc[iswap][i][idim] = 1; + else pbc[iswap][i][idim] = -1; + if (triclinic) { + if (idim == 1) pbc[iswap][i][5] = pbc[iswap][i][idim]; + if (idim == 2) pbc[iswap][i][4] = pbc[iswap][i][3] = pbc[iswap][i][idim]; + } + sbox[0] = MAX(oboxlo[0],lo2[0]); + sbox[1] = MAX(oboxlo[1],lo2[1]); + sbox[2] = MAX(oboxlo[2],lo2[2]); + sbox[3] = MIN(oboxhi[0],hi2[0]); + sbox[4] = MIN(oboxhi[1],hi2[1]); + sbox[5] = MIN(oboxhi[2],hi2[2]); + } + + for (int itype = 1; itype <= atom->ntypes; itype++) { + sbox_multiold[0] = sbox[0]; + sbox_multiold[1] = sbox[1]; + sbox_multiold[2] = sbox[2]; + sbox_multiold[3] = sbox[3]; + sbox_multiold[4] = sbox[4]; + sbox_multiold[5] = sbox[5]; + if (idir == 0) { + sbox_multiold[idim] = sublo[idim]; + if (i < noverlap1) + sbox_multiold[3+idim] = + MIN(sbox_multiold[3+idim]+cutghostmultiold[itype][idim],subhi[idim]); + else + sbox_multiold[3+idim] = + MIN(sbox_multiold[3+idim]-prd[idim]+cutghostmultiold[itype][idim], + subhi[idim]); + } else { + if (i < noverlap1) + sbox_multiold[idim] = + MAX(sbox_multiold[idim]-cutghostmultiold[itype][idim],sublo[idim]); + else + sbox_multiold[idim] = + MAX(sbox_multiold[idim]+prd[idim]-cutghostmultiold[itype][idim], + sublo[idim]); + sbox_multiold[3+idim] = subhi[idim]; + } + + if (idim >= 1) { + if (sbox_multiold[0] == oboxlo[0]) + sbox_multiold[0] -= cutghostmultiold[itype][idim]; + if (sbox_multiold[3] == oboxhi[0]) + sbox_multiold[3] += cutghostmultiold[itype][idim]; + } + if (idim == 2) { + if (sbox_multiold[1] == oboxlo[1]) + sbox_multiold[1] -= cutghostmultiold[itype][idim]; + if (sbox_multiold[4] == oboxhi[1]) + sbox_multiold[4] += cutghostmultiold[itype][idim]; + } + + memcpy(sendbox_multiold[iswap][i][itype],sbox_multiold,6*sizeof(double)); } } } @@ -937,11 +1039,14 @@ void CommTiled::exchange() void CommTiled::borders() { - int i,m,n,nlast,nsend,nrecv,ngroup,ncount,ncountall; + int i,m,n,nlast,nsend,nrecv,ngroup,nprior,ncount,ncountall; double xlo,xhi,ylo,yhi,zlo,zhi; double *bbox; double **x; AtomVec *avec = atom->avec; + + // After exchanging, need to reconstruct collection array for border communication + if(mode == Comm::MULTI) neighbor->build_collection(0); // send/recv max one = max # of atoms in single send/recv for any swap // send/recv max all = max # of atoms in all sends/recvs within any swap @@ -1008,6 +1113,56 @@ void CommTiled::borders() smaxone = MAX(smaxone,ncount); ncountall += ncount; + } else if (mode == Comm::MULTI) { + + int* collection=neighbor->collection; + int icollection; + ncount = 0; + + if (!bordergroup) { + for (i = 0; i < nlast; i++) { + icollection=collection[i]; + bbox = sendbox_multi[iswap][m][icollection]; + xlo = bbox[0]; ylo = bbox[1]; zlo = bbox[2]; + xhi = bbox[3]; yhi = bbox[4]; zhi = bbox[5]; + if (x[i][0] >= xlo && x[i][0] < xhi && + x[i][1] >= ylo && x[i][1] < yhi && + x[i][2] >= zlo && x[i][2] < zhi) { + if (ncount == maxsendlist[iswap][m]) grow_list(iswap,m,ncount); + sendlist[iswap][m][ncount++] = i; + } + } + } else { + ngroup = atom->nfirst; + for (i = 0; i < ngroup; i++) { + icollection=collection[i]; + bbox = sendbox_multi[iswap][m][icollection]; + xlo = bbox[0]; ylo = bbox[1]; zlo = bbox[2]; + xhi = bbox[3]; yhi = bbox[4]; zhi = bbox[5]; + if (x[i][0] >= xlo && x[i][0] < xhi && + x[i][1] >= ylo && x[i][1] < yhi && + x[i][2] >= zlo && x[i][2] < zhi) { + if (ncount == maxsendlist[iswap][m]) grow_list(iswap,m,ncount); + sendlist[iswap][m][ncount++] = i; + } + } + for (i = atom->nlocal; i < nlast; i++) { + icollection=collection[i]; + bbox = sendbox_multi[iswap][m][icollection]; + xlo = bbox[0]; ylo = bbox[1]; zlo = bbox[2]; + xhi = bbox[3]; yhi = bbox[4]; zhi = bbox[5]; + if (x[i][0] >= xlo && x[i][0] < xhi && + x[i][1] >= ylo && x[i][1] < yhi && + x[i][2] >= zlo && x[i][2] < zhi) { + if (ncount == maxsendlist[iswap][m]) grow_list(iswap,m,ncount); + sendlist[iswap][m][ncount++] = i; + } + } + } + + sendnum[iswap][m] = ncount; + smaxone = MAX(smaxone,ncount); + ncountall += ncount; } else { int* type=atom->type; @@ -1017,7 +1172,7 @@ void CommTiled::borders() if (!bordergroup) { for (i = 0; i < nlast; i++) { itype=type[i]; - bbox = sendbox_multi[iswap][m][itype]; + bbox = sendbox_multiold[iswap][m][itype]; xlo = bbox[0]; ylo = bbox[1]; zlo = bbox[2]; xhi = bbox[3]; yhi = bbox[4]; zhi = bbox[5]; if (x[i][0] >= xlo && x[i][0] < xhi && @@ -1031,7 +1186,7 @@ void CommTiled::borders() ngroup = atom->nfirst; for (i = 0; i < ngroup; i++) { itype=type[i]; - bbox = sendbox_multi[iswap][m][itype]; + bbox = sendbox_multiold[iswap][m][itype]; xlo = bbox[0]; ylo = bbox[1]; zlo = bbox[2]; xhi = bbox[3]; yhi = bbox[4]; zhi = bbox[5]; if (x[i][0] >= xlo && x[i][0] < xhi && @@ -1043,7 +1198,7 @@ void CommTiled::borders() } for (i = atom->nlocal; i < nlast; i++) { itype=type[i]; - bbox = sendbox_multi[iswap][m][itype]; + bbox = sendbox_multiold[iswap][m][itype]; xlo = bbox[0]; ylo = bbox[1]; zlo = bbox[2]; xhi = bbox[3]; yhi = bbox[4]; zhi = bbox[5]; if (x[i][0] >= xlo && x[i][0] < xhi && @@ -1179,8 +1334,11 @@ void CommTiled::borders() // increment ghost atoms n = nrecvproc[iswap]; - if (n) + if (n) { + nprior = atom->nghost + atom->nlocal; atom->nghost += forward_recv_offset[iswap][n-1] + recvnum[iswap][n-1]; + if(neighbor->style == Neighbor::MULTI) neighbor->build_collection(nprior); + } } // For molecular systems we lose some bits for local atom indices due @@ -2120,6 +2278,7 @@ void CommTiled::allocate_swap(int n) pbc = new int**[n]; sendbox = new double**[n]; sendbox_multi = new double***[n]; + sendbox_multiold = new double***[n]; maxsendlist = new int*[n]; sendlist = new int**[n]; @@ -2134,6 +2293,7 @@ void CommTiled::allocate_swap(int n) pbc[i] = nullptr; sendbox[i] = nullptr; sendbox_multi[i] = nullptr; + sendbox_multiold[i] = nullptr; maxsendlist[i] = nullptr; sendlist[i] = nullptr; } @@ -2182,8 +2342,9 @@ void CommTiled::grow_swap_send(int i, int n, int nold) memory->create(pbc[i],n,6,"comm:pbc_flag"); memory->destroy(sendbox[i]); memory->create(sendbox[i],n,6,"comm:sendbox"); - memory->destroy(sendbox_multi[i]); - memory->create(sendbox_multi[i],n,atom->ntypes+1,6,"comm:sendbox_multi"); + grow_swap_send_multi(i,n); + memory->destroy(sendbox_multiold[i]); + memory->create(sendbox_multiold[i],n,atom->ntypes+1,6,"comm:sendbox_multiold"); delete [] maxsendlist[i]; maxsendlist[i] = new int[n]; @@ -2215,6 +2376,19 @@ void CommTiled::grow_swap_recv(int i, int n) size_reverse_send[i] = new int[n]; } + +/* ---------------------------------------------------------------------- + grow info for swap I for multi as ncollections can change +------------------------------------------------------------------------- */ + +void CommTiled::grow_swap_send_multi(int i, int n) +{ + memory->destroy(sendbox_multi[i]); + + if(ncollections > 0) + memory->create(sendbox_multi[i],n,ncollections,6,"comm:sendbox_multi"); +} + /* ---------------------------------------------------------------------- deallocate swap info ------------------------------------------------------------------------- */ @@ -2243,6 +2417,7 @@ void CommTiled::deallocate_swap(int n) memory->destroy(pbc[i]); memory->destroy(sendbox[i]); memory->destroy(sendbox_multi[i]); + memory->destroy(sendbox_multiold[i]); delete [] maxsendlist[i]; @@ -2265,6 +2440,7 @@ void CommTiled::deallocate_swap(int n) delete [] pbc; delete [] sendbox; delete [] sendbox_multi; + delete [] sendbox_multiold; delete [] maxsendlist; delete [] sendlist; diff --git a/src/comm_tiled.h b/src/comm_tiled.h index fa2c76e6e9..fdb175117e 100644 --- a/src/comm_tiled.h +++ b/src/comm_tiled.h @@ -77,9 +77,12 @@ class CommTiled : public Comm { double ***sendbox; // bounding box of atoms to send per swap/proc - double **cutghostmulti; // cutghost on a per-type basis + double **cutghostmulti; // cutghost on a per-collection basis + double **cutghostmultiold; // cutghost on a per-type basis double ****sendbox_multi; // bounding box of atoms to send // per swap/proc for multi comm + double ****sendbox_multiold; // bounding box of atoms to send + // per swap/proc for multi/old comm // exchange comm info, proc lists do not include self @@ -148,6 +151,7 @@ class CommTiled : public Comm { void grow_list(int, int, int); // reallocate sendlist for one swap/proc void allocate_swap(int); // allocate swap arrays void grow_swap_send(int, int, int); // grow swap arrays for send and recv + void grow_swap_send_multi(int, int); // grow multi swap arrays for send and recv void grow_swap_recv(int, int); void deallocate_swap(int); // deallocate swap arrays @@ -163,10 +167,6 @@ E: Cannot yet use comm_style tiled with triclinic box Self-explanatory. -E: Cannot yet use comm_style tiled with multi-mode comm - -Self-explanatory. - E: Communication cutoff for comm_style tiled cannot exceed periodic box length Self-explanatory. diff --git a/src/init b/src/init new file mode 100644 index 0000000000..617c9effa3 --- /dev/null +++ b/src/init @@ -0,0 +1,789 @@ +angle_charmm.cpp: if (comm->me == 0) { +angle_cosine.cpp: if (comm->me == 0) utils::sfread(FLERR,&k[1],sizeof(double),atom->nangletypes,fp,nullptr,error); +angle_cosine_periodic.cpp: if (comm->me == 0) { +angle_cosine_squared.cpp: if (comm->me == 0) { +angle.cpp: memory->create(eatom,comm->nthreads*maxeatom,"angle:eatom"); +angle.cpp: memory->create(vatom,comm->nthreads*maxvatom,6,"angle:vatom"); +angle.cpp: memory->create(cvatom,comm->nthreads*maxcvatom,9,"angle:cvatom"); +angle.cpp: double bytes = comm->nthreads*maxeatom * sizeof(double); +angle.cpp: bytes += comm->nthreads*maxvatom*6 * sizeof(double); +angle.cpp: bytes += comm->nthreads*maxcvatom*9 * sizeof(double); +angle_deprecated.cpp: if (lmp->comm->me == 0) +angle_harmonic.cpp: if (comm->me == 0) { +angle_hybrid.cpp: const int nthreads = comm->nthreads; +angle_hybrid.cpp: if (comm->nthreads > 1) { +angle_hybrid.cpp: int me = comm->me; +angle_table.cpp: if (comm->me == 0) { +angle_zero.cpp: if (comm->me == 0) { +atom.cpp: if (comm->layout != Comm::LAYOUT_TILED) { +atom.cpp: if (comm->myloc[0] == 0) sublo[0] -= epsilon[0]; +atom.cpp: if (comm->myloc[0] == comm->procgrid[0]-1) subhi[0] += epsilon[0]; +atom.cpp: if (comm->myloc[1] == 0) sublo[1] -= epsilon[1]; +atom.cpp: if (comm->myloc[1] == comm->procgrid[1]-1) subhi[1] += epsilon[1]; +atom.cpp: if (comm->myloc[2] == 0) sublo[2] -= epsilon[2]; +atom.cpp: if (comm->myloc[2] == comm->procgrid[2]-1) subhi[2] += epsilon[2]; +atom.cpp: if (comm->mysplit[0][0] == 0.0) sublo[0] -= epsilon[0]; +atom.cpp: if (comm->mysplit[0][1] == 1.0) subhi[0] += epsilon[0]; +atom.cpp: if (comm->mysplit[1][0] == 0.0) sublo[1] -= epsilon[1]; +atom.cpp: if (comm->mysplit[1][1] == 1.0) subhi[1] += epsilon[1]; +atom.cpp: if (comm->mysplit[2][0] == 0.0) sublo[2] -= epsilon[2]; +atom.cpp: if (comm->mysplit[2][1] == 1.0) subhi[2] += epsilon[2]; +atom.cpp: if (comm->me == 0) { +atom.cpp: called by comm->exchange() if atom_modify first group is set +atom.cpp: always called between comm->exchange() and comm->borders() +atom.cpp: if (comm->me == 0) +atom_map.cpp: int nper = static_cast (natoms/comm->nprocs); +atom_vec_body.cpp: printf("Proc %d, step %ld, flag %d\n",comm->me,update->ntimestep,flag); +atom_vec_body.cpp: printf("Proc %d, step %ld, flag %d\n",comm->me,update->ntimestep,flag); +atom_vec_body.cpp: printf("Proc %d, step %ld, flag %d\n",comm->me,update->ntimestep,flag); +atom_vec_body.cpp: printf("Proc %d, step %ld, flag %d\n",comm->me,update->ntimestep,flag); +atom_vec_body.cpp: printf("Proc %d, step %ld, flag %d\n",comm->me,update->ntimestep,flag); +atom_vec_body.cpp: printf("Proc %d, step %ld, flag %d\n",comm->me,update->ntimestep,flag); +atom_vec.cpp: f = memory->grow(atom->f,nmax*comm->nthreads,3,"atom:f"); +atom_vec.cpp: const int nthreads = threads[i] ? comm->nthreads : 1; +atom_vec.cpp: bytes += memory->usage(f,nmax*comm->nthreads,3); +atom_vec.cpp: const int nthreads = threads[i] ? comm->nthreads : 1; +atom_vec_hybrid.cpp: if (mass_pertype && mass_peratom && comm->me == 0) +atom_vec_hybrid.cpp: if ((ptr && strstr(ptr+1,dup)) && (comm->me == 0)) +atom_vec_line.cpp: //if (comm->me == 1 && update->ntimestep == 873) +atom_vec_line.cpp: printf("BAD vecline ptrs: %s: %d %d: %d\n",str,comm->me, +atom_vec_line.cpp: str,comm->me,update->ntimestep,count,nlocal_bonus); +balance.cpp: int *procgrid = comm->procgrid; +balance.cpp: if (style == BISECTION && comm->style == 0) +balance.cpp: // init entire system since comm->setup is done +balance.cpp: comm->setup(); +balance.cpp: comm->exchange(); +balance.cpp: if (comm->layout == Comm::LAYOUT_TILED && style != BISECTION) { +balance.cpp: if (comm->layout == Comm::LAYOUT_UNIFORM) { +balance.cpp: comm->layout = Comm::LAYOUT_NONUNIFORM; +balance.cpp: } else if (comm->layout == Comm::LAYOUT_NONUNIFORM) { +balance.cpp: comm->layout = Comm::LAYOUT_UNIFORM; +balance.cpp: } else if (comm->layout == Comm::LAYOUT_TILED) { +balance.cpp: comm->layout = Comm::LAYOUT_UNIFORM; +balance.cpp: else comm->layout = Comm::LAYOUT_NONUNIFORM; +balance.cpp: comm->xsplit[i] = i * 1.0/procgrid[0]; +balance.cpp: comm->xsplit[procgrid[0]] = 1.0; +balance.cpp: for (int i = 0; i <= procgrid[0]; i++) comm->xsplit[i] = user_xsplit[i]; +balance.cpp: comm->ysplit[i] = i * 1.0/procgrid[1]; +balance.cpp: comm->ysplit[procgrid[1]] = 1.0; +balance.cpp: for (int i = 0; i <= procgrid[1]; i++) comm->ysplit[i] = user_ysplit[i]; +balance.cpp: comm->zsplit[i] = i * 1.0/procgrid[2]; +balance.cpp: comm->zsplit[procgrid[2]] = 1.0; +balance.cpp: for (int i = 0; i <= procgrid[2]; i++) comm->zsplit[i] = user_zsplit[i]; +balance.cpp: comm->layout = Comm::LAYOUT_NONUNIFORM; +balance.cpp: comm->layout = Comm::LAYOUT_TILED; +balance.cpp: for (int i = 0; i <= comm->procgrid[0]; i++) +balance.cpp: mesg += fmt::format(" {:.8}",comm->xsplit[i]); +balance.cpp: for (int i = 0; i <= comm->procgrid[1]; i++) +balance.cpp: mesg += fmt::format(" {:.8}",comm->ysplit[i]); +balance.cpp: for (int i = 0; i <= comm->procgrid[2]; i++) +balance.cpp: mesg += fmt::format(" {:.8}",comm->zsplit[i]); +balance.cpp: if (outflag && comm->me == 0) { +balance.cpp: comm->rcbnew = 1; +balance.cpp: if (idim >= 0) comm->rcbcutfrac = (rcb->cut - boxlo[idim]) / prd[idim]; +balance.cpp: else comm->rcbcutfrac = 0.0; +balance.cpp: comm->rcbcutdim = idim; +balance.cpp: double (*mysplit)[2] = comm->mysplit; +balance.cpp: int max = MAX(comm->procgrid[0],comm->procgrid[1]); +balance.cpp: max = MAX(max,comm->procgrid[2]); +balance.cpp: if (comm->layout == Comm::LAYOUT_TILED) { +balance.cpp: int *procgrid = comm->procgrid; +balance.cpp: double *xsplit = comm->xsplit; +balance.cpp: double *ysplit = comm->ysplit; +balance.cpp: double *zsplit = comm->zsplit; +balance.cpp: int *procgrid = comm->procgrid; +balance.cpp: split = comm->xsplit; +balance.cpp: split = comm->ysplit; +balance.cpp: split = comm->zsplit; +balance.cpp: double *xsplit = comm->xsplit; +balance.cpp: double *ysplit = comm->ysplit; +balance.cpp: double *zsplit = comm->zsplit; +balance.cpp: int nx = comm->procgrid[0]; +balance.cpp: int ny = comm->procgrid[1]; +balance.cpp: int nz = comm->procgrid[2]; +bond.cpp: memory->create(eatom,comm->nthreads*maxeatom,"bond:eatom"); +bond.cpp: memory->create(vatom,comm->nthreads*maxvatom,6,"bond:vatom"); +bond.cpp: if (comm->me == 0) { +bond.cpp: if (comm->me == 0) { +bond.cpp: double bytes = comm->nthreads*maxeatom * sizeof(double); +bond.cpp: bytes += comm->nthreads*maxvatom*6 * sizeof(double); +bond_deprecated.cpp: if (lmp->comm->me == 0) +bond_fene.cpp: if (comm->me == 0) +bond_fene.cpp: if (comm->me == 0) { +bond_fene_expand.cpp: if (comm->me == 0) +bond_fene_expand.cpp: if (comm->me == 0) { +bond_gromos.cpp: if (comm->me == 0) { +bond_harmonic.cpp: if (comm->me == 0) { +bond_hybrid.cpp: const int nthreads = comm->nthreads; +bond_hybrid.cpp: int me = comm->me; +bond_morse.cpp: if (comm->me == 0) { +bond_nonlinear.cpp: if (comm->me == 0) { +bond_quartic.cpp: if (comm->me == 0) { +bond_table.cpp: if (comm->me == 0) { +bond_zero.cpp: if (comm->me == 0) { +change_box.cpp: if (comm->me == 0) utils::logmesg(lmp,"Changing box ...\n"); +change_box.cpp: if (natoms != atom->natoms && comm->me == 0) +comm_brick.cpp: if (oldcomm->layout == Comm::LAYOUT_TILED) +comm_brick.cpp: layout = oldcomm->layout; +comm.cpp: if (oldcomm->grid2proc) { +comm.cpp: memcpy(&grid2proc[0][0][0],&oldcomm->grid2proc[0][0][0], +comm.cpp: memcpy(xsplit,oldcomm->xsplit,(procgrid[0]+1)*sizeof(double)); +comm.cpp: memcpy(ysplit,oldcomm->ysplit,(procgrid[1]+1)*sizeof(double)); +comm.cpp: memcpy(zsplit,oldcomm->zsplit,(procgrid[2]+1)*sizeof(double)); +comm.cpp: if (oldcomm->cutusermulti) { +comm.cpp: memcpy(cutusermulti,oldcomm->cutusermulti,ncollections); +comm.cpp: if (oldcomm->cutusermultiold) { +comm.cpp: memcpy(cutusermultiold,oldcomm->cutusermultiold,atom->ntypes+1); +comm.cpp: int n = strlen(oldcomm->customfile) + 1; +comm.cpp: strcpy(customfile,oldcomm->customfile); +comm.cpp: int n = strlen(oldcomm->outfile) + 1; +comm.cpp: strcpy(outfile,oldcomm->outfile); +comm_tiled.cpp: layout = oldcomm->layout; +compute_adf.cpp: if (mycutneigh > comm->cutghostuser) +compute_aggregate_atom.cpp: if (count > 1 && comm->me == 0) +compute_aggregate_atom.cpp: comm->forward_comm_compute(this); +compute_aggregate_atom.cpp: comm->forward_comm_compute(this); +compute_aggregate_atom.cpp: comm->reverse_comm_compute(this); +compute_bond_local.cpp: if (velflag && !comm->ghost_velocity) ghostvelflag = 1; +compute_bond_local.cpp: if (ghostvelflag && !initflag) comm->forward_comm_compute(this); +compute_centro_atom.cpp: if (count > 1 && comm->me == 0) +compute_centroid_stress_atom.cpp: comm->reverse_comm_compute(this); +compute_chunk_atom.cpp: int nprocs = comm->nprocs; +compute_chunk_atom.cpp: comm->ring(n,sizeof(int),list,1,idring,nullptr,(void *)this,0); +compute_chunk_atom.cpp: callback from comm->ring() +compute_chunk_atom.cpp: if (flagall && comm->me == 0) +compute_cluster_atom.cpp: if (count > 1 && comm->me == 0) +compute_cluster_atom.cpp: comm->forward_comm_compute(this); +compute_cluster_atom.cpp: comm->forward_comm_compute(this); +compute_cluster_atom.cpp: comm->forward_comm_compute(this); +compute_cna_atom.cpp: comm->me == 0) +compute_cna_atom.cpp: if (count > 1 && comm->me == 0) +compute_cna_atom.cpp: if (nerrorall && comm->me == 0) +compute_cna_atom.cpp: if (nerrorall && comm->me == 0) +compute_contact_atom.cpp: if (count > 1 && comm->me == 0) +compute_contact_atom.cpp: if (force->newton_pair) comm->reverse_comm_compute(this); +compute_coord_atom.cpp: comm->forward_comm_compute(this); +compute_deprecated.cpp: if (lmp->comm->me == 0) +compute_erotate_sphere_atom.cpp: if (count > 1 && comm->me == 0) +compute_fragment_atom.cpp: if (count > 1 && comm->me == 0) +compute_fragment_atom.cpp: comm->forward_comm_compute(this); +compute_fragment_atom.cpp: comm->forward_comm_compute(this); +compute_group_group.cpp: if ((fabs(e_correction) > SMALL) && (comm->me == 0)) +compute_hexorder_atom.cpp: if (count > 1 && comm->me == 0) +compute_ke_atom.cpp: if (count > 1 && comm->me == 0) +compute_orientorder_atom.cpp: if (count > 1 && comm->me == 0) +compute_pe_atom.cpp: comm->reverse_comm_compute(this); +compute_property_atom.cpp: int me = comm->me; +compute_rdf.cpp: cutghost = MAX(force->pair->cutforce+skin,comm->cutghostuser); +compute_rdf.cpp: cutghost = comm->cutghostuser; +compute_rdf.cpp: if (comm->me == 0) +compute_stress_atom.cpp: comm->reverse_comm_compute(this); +compute_temp_deform.cpp: comm->me == 0) +compute_temp_deform.cpp: if (i == modify->nfix && comm->me == 0) +create_atoms.cpp: if (comm->layout != Comm::LAYOUT_TILED) { +create_atoms.cpp: if (comm->myloc[0] == 0) sublo[0] -= epsilon[0]; +create_atoms.cpp: if (comm->myloc[0] == comm->procgrid[0]-1) subhi[0] -= 2.0*epsilon[0]; +create_atoms.cpp: if (comm->myloc[1] == 0) sublo[1] -= epsilon[1]; +create_atoms.cpp: if (comm->myloc[1] == comm->procgrid[1]-1) subhi[1] -= 2.0*epsilon[1]; +create_atoms.cpp: if (comm->myloc[2] == 0) sublo[2] -= epsilon[2]; +create_atoms.cpp: if (comm->myloc[2] == comm->procgrid[2]-1) subhi[2] -= 2.0*epsilon[2]; +create_atoms.cpp: if (comm->mysplit[0][0] == 0.0) sublo[0] -= epsilon[0]; +create_atoms.cpp: if (comm->mysplit[0][1] == 1.0) subhi[0] -= 2.0*epsilon[0]; +create_atoms.cpp: if (comm->mysplit[1][0] == 0.0) sublo[1] -= epsilon[1]; +create_atoms.cpp: if (comm->mysplit[1][1] == 1.0) subhi[1] -= 2.0*epsilon[1]; +create_atoms.cpp: if (comm->mysplit[2][0] == 0.0) sublo[2] -= epsilon[2]; +create_atoms.cpp: if (comm->mysplit[2][1] == 1.0) subhi[2] -= 2.0*epsilon[2]; +create_bonds.cpp: // init entire system since comm->borders and neighbor->build is done +create_bonds.cpp: if (rmax > neighbor->cutneighmin && comm->me == 0) +create_bonds.cpp: comm->setup(); +create_bonds.cpp: comm->exchange(); +create_bonds.cpp: comm->borders(); +create_bonds.cpp: if (comm->me == 0) +create_box.cpp: comm->set_proc_grid(); +delete_atoms.cpp: } else if (comm->me == 0) +delete_atoms.cpp: if (comm->me == 0) { +delete_atoms.cpp: if (comm->me == 0) utils::logmesg(lmp,"System init for delete_atoms ...\n"); +delete_atoms.cpp: // init entire system since comm->borders and neighbor->build is done +delete_atoms.cpp: if (cut > neighbor->cutneighmin && comm->me == 0) +delete_atoms.cpp: comm->setup(); +delete_atoms.cpp: comm->exchange(); +delete_atoms.cpp: comm->borders(); +delete_atoms.cpp: RanMars *random = new RanMars(lmp,seed + comm->me); +delete_atoms.cpp: // pass list to all other procs via comm->ring() +delete_atoms.cpp: comm->ring(n,sizeof(tagint),list,1,bondring,nullptr,(void *)this); +delete_atoms.cpp: // pass list to all other procs via comm->ring() +delete_atoms.cpp: comm->ring(n,sizeof(tagint),list,1,molring,nullptr,(void *)this); +delete_atoms.cpp: callback from comm->ring() in delete_bond() +delete_atoms.cpp: callback from comm->ring() in delete_molecule() +delete_bonds.cpp: // init entire system since comm->borders is done +delete_bonds.cpp: if (comm->me == 0) utils::logmesg(lmp,"System init for delete_bonds ...\n"); +delete_bonds.cpp: if (comm->me == 0) utils::logmesg(lmp,"Deleting bonds ...\n"); +delete_bonds.cpp: comm->setup(); +delete_bonds.cpp: comm->exchange(); +delete_bonds.cpp: comm->borders(); +delete_bonds.cpp: if (comm->me == 0) { +deprecated.cpp: if (lmp->comm->me == 0) +deprecated.cpp: if (lmp->comm->me == 0) +dihedral_charmm.cpp: if (comm->me == 0) { +dihedral_charmmfsw.cpp: if (comm->me == 0) { +dihedral.cpp: memory->create(eatom,comm->nthreads*maxeatom,"dihedral:eatom"); +dihedral.cpp: memory->create(vatom,comm->nthreads*maxvatom,6,"dihedral:vatom"); +dihedral.cpp: memory->create(cvatom,comm->nthreads*maxcvatom,9,"dihedral:cvatom"); +dihedral.cpp: double bytes = comm->nthreads*maxeatom * sizeof(double); +dihedral.cpp: bytes += comm->nthreads*maxvatom*6 * sizeof(double); +dihedral.cpp: bytes += comm->nthreads*maxcvatom*9 * sizeof(double); +dihedral_deprecated.cpp: if (lmp->comm->me == 0) +dihedral_harmonic.cpp: if (comm->me == 0) { +dihedral_helix.cpp: if (comm->me == 0) { +dihedral_hybrid.cpp: int me = comm->me; +dihedral_multi_harmonic.cpp: if (comm->me == 0) { +dihedral_opls.cpp: if (comm->me == 0) { +displace_atoms.cpp: if (comm->me == 0) utils::logmesg(lmp,"Displacing atoms ...\n"); +displace_atoms.cpp: if (natoms != atom->natoms && comm->me == 0) +domain.cpp: else if (comm->me == 0) +domain.cpp: uses comm->xyz_split or comm->mysplit +domain.cpp: if (comm->layout != Comm::LAYOUT_TILED) { +domain.cpp: int *myloc = comm->myloc; +domain.cpp: double *xsplit = comm->xsplit; +domain.cpp: double *ysplit = comm->ysplit; +domain.cpp: double *zsplit = comm->zsplit; +domain.cpp: double (*mysplit)[2] = comm->mysplit; +domain.cpp: uses comm->xyz_split or comm->mysplit +domain.cpp: if (comm->layout != Comm::LAYOUT_TILED) { +domain.cpp: int *myloc = comm->myloc; +domain.cpp: int *procgrid = comm->procgrid; +domain.cpp: double *xsplit = comm->xsplit; +domain.cpp: double *ysplit = comm->ysplit; +domain.cpp: double *zsplit = comm->zsplit; +domain.cpp: double (*mysplit)[2] = comm->mysplit; +domain.cpp: comm->forward_comm_array(3,unwrap); +domain.cpp: if (flagall && comm->me == 0) +domain.cpp: if (all && comm->me == 0) +domain.cpp: if (all && comm->me == 0) +domain.cpp: if (flagall && comm->me == 0) +domain.cpp: since may lead to lost atoms in comm->exchange() +domain.cpp: if (flagall && comm->me == 0) +domain.cpp: if ((flag_all > 0) && (comm->me == 0)) +domain.cpp: if (comm->me == 0) { +dump_deprecated.cpp: if (lmp->comm->me == 0) +dump_image.cpp: comm->forward_comm_dump(this); +dump_movie.cpp: if ((comm->me == 0) && (fp == nullptr)) { +finish.cpp: const int nthreads = comm->nthreads; +fix_balance.cpp: if (lbstyle == BISECTION && comm->style == 0) +fix_balance.cpp: compute final imbalance factor based on nlocal after comm->exchange() +fix_balance.cpp: // invoke balancer and reset comm->uniform flag +fix_balance.cpp: comm->layout = Comm::LAYOUT_NONUNIFORM; +fix_balance.cpp: comm->layout = Comm::LAYOUT_TILED; +fix_balance.cpp: // since may lead to lost atoms in comm->exchange() +fix_balance.cpp: // else allow caller's comm->exchange() to do it +fix_balance.cpp: // b/c atoms may migrate again in comm->exchange() +fix_balance.cpp: // can only be done after atoms migrate in comm->exchange() +fix_box_relax.cpp: if (temperature->igroup != 0 && comm->me == 0) +fix_cmap.cpp: if (comm->me == 0) { +fix_cmap.cpp: if (comm->me == 0) +fix_cmap.cpp: if (comm->me == 0) fclose(fp); +fix_cmap.cpp: if (comm->me == 0) { +fix_deform.cpp: if (comm->me == 0) { +fix_deprecated.cpp: if (lmp->comm->me == 0) +fix_deprecated.cpp: if (lmp->comm->me == 0) +fix_dt_reset.cpp: strcmp(output->dump[i]->style,"xtc") == 0) && comm->me == 0) +fix_group.cpp: if (warn && comm->me == 0) +fix_halt.cpp: if (comm->me == 0 && msgflag == YESMSG) error->message(FLERR,message); +fix_langevin.cpp: random = new RanMars(lmp,seed + comm->me); +fix_langevin.cpp: if (temperature->igroup != igroup && comm->me == 0) +fix_move.cpp: if (comm->me == 0) { +fix_neigh_history.cpp: int nmypage = comm->nthreads; +fix_neigh_history.cpp: comm->reverse_comm_fix(this,0); +fix_neigh_history.cpp: comm->reverse_comm_fix_variable(this); +fix_neigh_history.cpp: int nmypage = comm->nthreads; +fix_neigh_history.cpp: if (comm->me == 0) { +fix_nh.cpp: if (comm->me == 0) { +fix_nh.cpp: if (temperature->igroup != 0 && comm->me == 0) +fix_nve_limit.cpp: if (comm->me == 0) +fix_pour.cpp: if (comm->layout != Comm::LAYOUT_TILED) { +fix_pour.cpp: if (comm->myloc[2] == comm->procgrid[2]-1 && +fix_pour.cpp: if (comm->mysplit[2][1] == 1.0 && +fix_pour.cpp: if (comm->layout != Comm::LAYOUT_TILED) { +fix_pour.cpp: if (comm->myloc[1] == comm->procgrid[1]-1 && +fix_pour.cpp: if (comm->mysplit[1][1] == 1.0 && +fix_press_berendsen.cpp: if (temperature->igroup != 0 && comm->me == 0) +fix_property_atom.cpp: if (flag && comm->me == 0) +fix_recenter.cpp: if (flag && comm->me == 0) +fix_restrain.cpp: comm->me,update->ntimestep)); +fix_restrain.cpp: comm->me,update->ntimestep)); +fix_restrain.cpp: comm->me,update->ntimestep)); +fix_restrain.cpp: comm->me,update->ntimestep)); +fix_restrain.cpp: ids[m][2],comm->me,update->ntimestep)); +fix_restrain.cpp: ids[m][2],comm->me,update->ntimestep)); +fix_restrain.cpp: ids[m][2],ids[m][3],comm->me, +fix_restrain.cpp: ids[m][2],ids[m][3],comm->me, +fix_restrain.cpp: comm->me,x[i1][0],x[i1][1],x[i1][2], +fix_restrain.cpp: comm->me,x[i2][0],x[i2][1],x[i2][2], +fix_restrain.cpp: comm->me,x[i3][0],x[i3][1],x[i3][2], +fix_restrain.cpp: comm->me,x[i4][0],x[i4][1],x[i4][2]); +fix_spring_chunk.cpp: if (comm->me == 0) { +fix_spring_chunk.cpp: if (comm->me == 0) +fix_spring_rg.cpp: if (comm->me == 0) { +fix_store.cpp: // PERATOM may be comm->exchanged before filled by caller +fix_store.cpp: if (comm->me == 0) { +fix_temp_berendsen.cpp: if (temperature->igroup != igroup && comm->me == 0) +fix_temp_berendsen.cpp: if (comm->me == 0) { +fix_temp_csld.cpp: random = new RanMars(lmp,seed + comm->me); +fix_temp_csld.cpp: if (temperature->igroup != igroup && comm->me == 0) +fix_temp_csld.cpp: int nsize = PRNGSIZE*comm->nprocs+2; // pRNG state per proc + nprocs + energy +fix_temp_csld.cpp: if (comm->me == 0) { +fix_temp_csld.cpp: list[1] = comm->nprocs; +fix_temp_csld.cpp: if (comm->me == 0) { +fix_temp_csld.cpp: if (nprocs != comm->nprocs) { +fix_temp_csld.cpp: if (comm->me == 0) +fix_temp_csld.cpp: } else random->set_state(list+2+comm->me*103); +fix_temp_csvr.cpp: random = new RanMars(lmp,seed + comm->me); +fix_temp_csvr.cpp: if (comm->me == 0) { +fix_temp_csvr.cpp: if (temperature->igroup != igroup && comm->me == 0) +fix_temp_csvr.cpp: int nsize = PRNGSIZE*comm->nprocs+2; // pRNG state per proc + nprocs + energy +fix_temp_csvr.cpp: if (comm->me == 0) { +fix_temp_csvr.cpp: list[1] = comm->nprocs; +fix_temp_csvr.cpp: if (comm->me == 0) { +fix_temp_csvr.cpp: if (nprocs != comm->nprocs) { +fix_temp_csvr.cpp: if (comm->me == 0) +fix_temp_csvr.cpp: } else random->set_state(list+2+comm->me*103); +fix_temp_rescale.cpp: if (temperature->igroup != igroup && comm->me == 0) +fix_temp_rescale.cpp: if (comm->me == 0) { +fix_wall_gran_region.cpp: if (comm->me) return; +fix_wall_reflect.cpp: if (nrigid && comm->me == 0) +force.cpp: if (comm->me == 0) { +group.cpp: // pass list to all other procs via comm->ring() +group.cpp: comm->ring(n,sizeof(tagint),list,1,molring,nullptr,(void *)this); +group.cpp: callback from comm->ring() +imbalance_neigh.cpp: if (comm->me == 0 && !did_warn) +improper.cpp: memory->create(eatom,comm->nthreads*maxeatom,"improper:eatom"); +improper.cpp: memory->create(vatom,comm->nthreads*maxvatom,6,"improper:vatom"); +improper.cpp: memory->create(cvatom,comm->nthreads*maxcvatom,9,"improper:cvatom"); +improper.cpp: double bytes = comm->nthreads*maxeatom * sizeof(double); +improper.cpp: bytes += comm->nthreads*maxvatom*6 * sizeof(double); +improper_cvff.cpp: if (comm->me == 0) { +improper_deprecated.cpp: if (lmp->comm->me == 0) +improper_harmonic.cpp: if (comm->me == 0) { +improper_hybrid.cpp: int me = comm->me; +improper_umbrella.cpp: if (comm->me == 0) { +info.cpp: if (comm->me != 0) return; +info.cpp: commstyles[comm->style], commlayout[comm->layout], +info.cpp: comm->ghost_velocity ? "yes" : "no"); +info.cpp: if (comm->mode == 0) +info.cpp: comm->get_comm_cutoff()); +info.cpp: if (comm->mode == 1) { +info.cpp: if (comm->cutusermulti) cut = MAX(cut,comm->cutusermulti[i]); +info.cpp: fmt::print(out,"Nprocs = {}, Nthreads = {}\n",comm->nprocs,comm->nthreads); +info.cpp: fmt::print(out,"Processor grid = {} x {} x {}\n",comm->procgrid[0], +info.cpp: comm->procgrid[1], comm->procgrid[2]); +info.cpp: style = commstyles[comm->style]; +info.cpp: bytes += comm->memory_usage(); +input.cpp: comm->modify_params(narg,arg); +input.cpp: if (comm->style == 0) return; +input.cpp: if (comm->style == 1) return; +input.cpp: comm->set_processors(narg,arg); +irregular.cpp: can be used in place of comm->exchange() +irregular.cpp: if (!preassign) comm->coord2proc_setup(); +irregular.cpp: mproclist[nsendatom] = comm->coord2proc(x[i],igx,igy,igz); +irregular.cpp: if not, caller can decide to use comm->exchange() instead +irregular.cpp: if (comm->layout == Comm::LAYOUT_TILED) return 1; +irregular.cpp: // cannot check via comm->procneigh since it ignores PBC +irregular.cpp: int *myloc = comm->myloc; +irregular.cpp: int *procgrid = comm->procgrid; +irregular.cpp: comm->coord2proc(x[i],igx,igy,igz); +kspace.cpp: if ((qsqsum == 0.0) && (comm->me == 0) && warn_nocharge && warning_flag) { +kspace.cpp: if (warn_nonneutral == 1 && comm->me == 0) error->warning(FLERR,message); +kspace.cpp: if (comm->me == 0) { +kspace.cpp: if ((table_accuracy > spr) && (comm->me == 0)) +kspace.cpp: if (slab_volfactor < 2.0 && comm->me == 0) +kspace_deprecated.cpp: if (lmp->comm->me == 0) +lammps.cpp: const int me = comm->me; +lammps.cpp: comm->init(); // comm must come after force, modify, neighbor, atom +lattice.cpp: if (comm->me == 0) +library.cpp: && (lmp->comm->me == 0)) { +library.cpp: && (lmp->comm->me == 0)) { +library.cpp: lmp->comm->set_proc_grid(); +library.cpp: if (strcmp(keyword,"world_rank") == 0) return lmp->comm->me; +library.cpp: if (strcmp(keyword,"world_size") == 0) return lmp->comm->nprocs; +library.cpp: if (strcmp(keyword,"nthreads") == 0) return lmp->comm->nthreads; +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: int nprocs = lmp->comm->nprocs; +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: int nprocs = lmp->comm->nprocs; +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) +library.cpp: if (lmp->comm->me == 0) lmp->error->warning(FLERR,msg); +min.cpp: if (comm->me == 0) +min.cpp: if (comm->me == 0 && screen) { +min.cpp: if (comm->me == 0 && screen) +min.cpp: comm->setup(); +min.cpp: comm->exchange(); +min.cpp: comm->borders(); +min.cpp: // atoms may have migrated in comm->exchange() +min.cpp: if (force->newton) comm->reverse_comm(); +min.cpp: comm->setup(); +min.cpp: comm->exchange(); +min.cpp: comm->borders(); +min.cpp: // atoms may have migrated in comm->exchange() +min.cpp: if (force->newton) comm->reverse_comm(); +min.cpp: comm->forward_comm(); +min.cpp: comm->setup(); +min.cpp: comm->exchange(); +min.cpp: comm->borders(); +min.cpp: comm->reverse_comm(); +min_fire.cpp: if (comm->me == 0 && logfile) { +modify.cpp: if (comm->me == 0 && checkall) +modify.cpp: if (fix[ifix]->igroup != igroup && comm->me == 0) +modify.cpp: if (comm->me == 0) +modify.cpp: if (comm->me == 0) +modify.cpp: int me = comm->me; +modify.cpp: int me = comm->me; +modify.cpp: if (flag && comm->me == 0) { +modify.cpp: if (flag && comm->me == 0) { +molecule.cpp: me = comm->me; +nbin_multi.cpp: // bsubbox lo/hi = bounding box of my subdomain extended by comm->cutghost +nbin_multi.cpp: // include dimension-dependent extension via comm->cutghost +nbin_multi.cpp: double *cutghost = comm->cutghost; +nbin_standard.cpp: // bsubbox lo/hi = bounding box of my subdomain extended by comm->cutghost +nbin_standard.cpp: // include dimension-dependent extension via comm->cutghost +nbin_standard.cpp: double *cutghost = comm->cutghost; +neighbor.cpp: if (!same && comm->me == 0) print_pairwise_info(); +neighbor.cpp: if (comm->me == 0) printf("SAME flag %d\n",same); +neighbor.cpp: const double cutghost = MAX(cutneighmax,comm->cutghostuser); +neigh_list.cpp: int nmypage = comm->nthreads; +neigh_list.cpp: if (comm->me != 0) return; +neigh_list.cpp: int nmypage = comm->nthreads; +ntopo.cpp: me = comm->me; +ntopo.cpp: nprocs = comm->nprocs; +output.cpp: if (thermo->modified && comm->me == 0) +output.cpp: if (strchr(arg[1],'%')) multiproc = comm->nprocs; +output.cpp: mbavg /= comm->nprocs; +output.cpp: if (comm->me == 0) +pair_beck.cpp: int me = comm->me; +pair_beck.cpp: int me = comm->me; +pair_born_coul_dsf.cpp: int me = comm->me; +pair_born_coul_dsf.cpp: if (comm->me == 0) { +pair_born_coul_wolf.cpp: int me = comm->me; +pair_born_coul_wolf.cpp: if (comm->me == 0) { +pair_born.cpp: int me = comm->me; +pair_born.cpp: if (comm->me == 0) { +pair_brownian.cpp: random = new RanMars(lmp,seed + comm->me); +pair_brownian.cpp: if (force->newton_pair == 0 && comm->me == 0) +pair_brownian.cpp: int me = comm->me; +pair_brownian.cpp: int me = comm->me; +pair_brownian.cpp: random = new RanMars(lmp,seed + comm->me); +pair_buck_coul_cut.cpp: int me = comm->me; +pair_buck_coul_cut.cpp: if (comm->me == 0) { +pair_buck.cpp: int me = comm->me; +pair_buck.cpp: if (comm->me == 0) { +pair_colloid.cpp: if (comm->me == 0) utils::sfread(FLERR,&setflag[i][j],sizeof(int),1,fp,nullptr,error); +pair_colloid.cpp: if (comm->me == 0) { +pair_colloid.cpp: int me = comm->me; +pair_coul_cut.cpp: int me = comm->me; +pair_coul_cut.cpp: if (comm->me == 0) { +pair_coul_debye.cpp: if (comm->me == 0) { +pair_coul_dsf.cpp: int me = comm->me; +pair_coul_dsf.cpp: if (comm->me == 0) { +pair_coul_streitz.cpp: if (comm->me == 0) { +pair_coul_streitz.cpp: if (comm->me != 0) { +pair_coul_wolf.cpp: int me = comm->me; +pair_coul_wolf.cpp: if (comm->me == 0) { +pair.cpp: if (tail_flag && domain->nonperiodic && comm->me == 0) +pair.cpp: if (!compute_flag && tail_flag && comm->me == 0) +pair.cpp: if (!compute_flag && offset_flag && comm->me == 0) +pair.cpp: if (flag && comm->me == 0) +pair.cpp: if (comm->me == 0) +pair.cpp: if (comm->me == 0) +pair.cpp: memory->create(eatom,comm->nthreads*maxeatom,"pair:eatom"); +pair.cpp: memory->create(vatom,comm->nthreads*maxvatom,6,"pair:vatom"); +pair.cpp: memory->create(cvatom,comm->nthreads*maxcvatom,9,"pair:cvatom"); +pair.cpp: if (comm->me == 0) { +pair.cpp: if ((comm->me == 0) && (epair)) +pair.cpp: if (comm->me == 0) +pair.cpp: if (comm->me == 0) fprintf(fp,"%d %.15g %.15g %.15g\n",i+1,r,e,f); +pair.cpp: if (comm->me == 0) fclose(fp); +pair.cpp: double bytes = comm->nthreads*maxeatom * sizeof(double); +pair.cpp: bytes += comm->nthreads*maxvatom*6 * sizeof(double); +pair.cpp: bytes += comm->nthreads*maxcvatom*9 * sizeof(double); +pair_deprecated.cpp: if (lmp->comm->me == 0) +pair_deprecated.cpp: if (lmp->comm->me == 0) +pair_dpd.cpp: random = new RanMars(lmp,seed + comm->me); +pair_dpd.cpp: if (comm->ghost_velocity == 0) +pair_dpd.cpp: if (force->newton_pair == 0 && comm->me == 0) error->warning(FLERR, +pair_dpd.cpp: int me = comm->me; +pair_dpd.cpp: if (comm->me == 0) { +pair_dpd.cpp: random = new RanMars(lmp,seed + comm->me); +pair_dpd_tstat.cpp: random = new RanMars(lmp,seed + comm->me); +pair_dpd_tstat.cpp: int me = comm->me; +pair_dpd_tstat.cpp: if (comm->me == 0) { +pair_dpd_tstat.cpp: random = new RanMars(lmp,seed + comm->me); +pair_gauss.cpp: int me = comm->me; +pair_gauss.cpp: if (comm->me == 0) { +pair_gran_hertz_history.cpp: comm->forward_comm_pair(this); +pair_gran_hooke.cpp: comm->forward_comm_pair(this); +pair_gran_hooke_history.cpp: comm->forward_comm_pair(this); +pair_gran_hooke_history.cpp: if (comm->ghost_velocity == 0) +pair_gran_hooke_history.cpp: int me = comm->me; +pair_gran_hooke_history.cpp: if (comm->me == 0) { +pair_granular.cpp: comm->forward_comm_pair(this); +pair_granular.cpp: if (comm->ghost_velocity == 0) +pair_granular.cpp: int me = comm->me; +pair_hybrid.cpp: int me = comm->me; +pair_lj96_cut.cpp: int me = comm->me; +pair_lj96_cut.cpp: int me = comm->me; +pair_lj_charmm_coul_charmm.cpp: int me = comm->me; +pair_lj_charmm_coul_charmm.cpp: if (comm->me == 0) { +pair_lj_charmmfsw_coul_charmmfsh.cpp: if ((comm->me == 0) && (force->qqr2e != force->qqr2e_charmm_real)) +pair_lj_charmmfsw_coul_charmmfsh.cpp: if ((comm->me == 0) && (force->qqr2e == force->qqr2e_charmm_real)) +pair_lj_charmmfsw_coul_charmmfsh.cpp: int me = comm->me; +pair_lj_charmmfsw_coul_charmmfsh.cpp: if (comm->me == 0) { +pair_lj_cubic.cpp: int me = comm->me; +pair_lj_cubic.cpp: int me = comm->me; +pair_lj_cut_coul_cut.cpp: int me = comm->me; +pair_lj_cut_coul_cut.cpp: if (comm->me == 0) { +pair_lj_cut_coul_debye.cpp: if (comm->me == 0) { +pair_lj_cut_coul_dsf.cpp: int me = comm->me; +pair_lj_cut_coul_dsf.cpp: if (comm->me == 0) { +pair_lj_cut_coul_wolf.cpp: int me = comm->me; +pair_lj_cut_coul_wolf.cpp: int me = comm->me; +pair_lj_cut.cpp: int me = comm->me; +pair_lj_cut.cpp: int me = comm->me; +pair_lj_cut_tip4p_cut.cpp: int me = comm->me; +pair_lj_cut_tip4p_cut.cpp: if (comm->me == 0) { +pair_lj_expand.cpp: int me = comm->me; +pair_lj_expand.cpp: if (comm->me == 0) { +pair_lj_gromacs_coul_gromacs.cpp: int me = comm->me; +pair_lj_gromacs_coul_gromacs.cpp: if (comm->me == 0) { +pair_lj_gromacs.cpp: int me = comm->me; +pair_lj_gromacs.cpp: int me = comm->me; +pair_lj_smooth.cpp: int me = comm->me; +pair_lj_smooth.cpp: int me = comm->me; +pair_lj_smooth_linear.cpp: int me = comm->me; +pair_lj_smooth_linear.cpp: int me = comm->me; +pair_lubricate.cpp: // no need to do this if not shearing since comm->ghost_velocity is set +pair_lubricate.cpp: comm->forward_comm_pair(this); +pair_lubricate.cpp: if (comm->ghost_velocity == 0) +pair_lubricate.cpp: int me = comm->me; +pair_lubricate.cpp: int me = comm->me; +pair_lubricate_poly.cpp: // no need to do this if not shearing since comm->ghost_velocity is set +pair_lubricate_poly.cpp: comm->forward_comm_pair(this); +pair_lubricate_poly.cpp: if (comm->ghost_velocity == 0) +pair_lubricateU.cpp: if (newton_pair) comm->reverse_comm(); +pair_lubricateU.cpp: comm->forward_comm_pair(this); +pair_lubricateU.cpp: if (newton_pair) comm->reverse_comm(); +pair_lubricateU.cpp: comm->forward_comm_pair(this); +pair_lubricateU.cpp: if (newton_pair) comm->reverse_comm(); +pair_lubricateU.cpp: comm->forward_comm_pair(this); +pair_lubricateU.cpp: if (newton_pair) comm->reverse_comm(); +pair_lubricateU.cpp: comm->forward_comm_pair(this); +pair_lubricateU.cpp: if (newton_pair) comm->reverse_comm(); +pair_lubricateU.cpp: comm->forward_comm_pair(this); +pair_lubricateU.cpp: if (newton_pair) comm->reverse_comm(); +pair_lubricateU.cpp: comm->forward_comm_pair(this); +pair_lubricateU.cpp: if (newton_pair) comm->reverse_comm(); // not really needed +pair_lubricateU.cpp: if (comm->ghost_velocity == 0) +pair_lubricateU.cpp: int me = comm->me; +pair_lubricateU.cpp: int me = comm->me; +pair_lubricateU_poly.cpp: if (newton_pair) comm->reverse_comm(); +pair_lubricateU_poly.cpp: comm->forward_comm_pair(this); +pair_lubricateU_poly.cpp: if (newton_pair) comm->reverse_comm(); +pair_lubricateU_poly.cpp: comm->forward_comm_pair(this); +pair_lubricateU_poly.cpp: if (newton_pair) comm->reverse_comm(); +pair_lubricateU_poly.cpp: comm->forward_comm_pair(this); +pair_lubricateU_poly.cpp: if (newton_pair) comm->reverse_comm(); // not really needed +pair_lubricateU_poly.cpp: if (comm->ghost_velocity == 0) +pair_lubricateU_poly.cpp: if (!comm->me) { +pair_mie_cut.cpp: int me = comm->me; +pair_mie_cut.cpp: int me = comm->me; +pair_morse.cpp: int me = comm->me; +pair_morse.cpp: if (comm->me == 0) { +pair_soft.cpp: int me = comm->me; +pair_soft.cpp: if (comm->me == 0) { +pair_table.cpp: if (comm->me == 0) { +pair_tip4p_cut.cpp: int me = comm->me; +pair_tip4p_cut.cpp: if (comm->me == 0) { +pair_ufm.cpp: int me = comm->me; +pair_ufm.cpp: int me = comm->me; +pair_yukawa.cpp: int me = comm->me; +pair_yukawa.cpp: if (comm->me == 0) { +pair_zbl.cpp: int me = comm->me; +pair_zbl.cpp: int me = comm->me; +pair_zero.cpp: int me = comm->me; +pair_zero.cpp: int me = comm->me; +potential_file_reader.cpp: if (comm->me != 0) { +random_mars.cpp: //if (comm->me == 0) printf("%d %ld %ld %g %ld\n", +read_data.cpp: if (comm->nprocs == 1) n = static_cast (atom->natoms); +read_data.cpp: else n = static_cast (LB_FACTOR * atom->natoms / comm->nprocs); +read_data.cpp: comm->set_proc_grid(); +read_data.cpp: comm->set_proc_grid(); +read_data.cpp: // do comm->init() but not comm->setup() b/c pair/neigh cutoffs not yet set +read_data.cpp: // need call to map_set() b/c comm->exchange clears atom map +read_data.cpp: if (comm->me == 0) +read_data.cpp: eof = comm->read_lines_from_file(fp,nchunk,MAXLINE,buffer); +read_data.cpp: eof = comm->read_lines_from_file(fp,nchunk,MAXLINE,buffer); +read_data.cpp: eof = comm->read_lines_from_file(fp,nchunk,MAXLINE,buffer); +read_data.cpp: eof = comm->read_lines_from_file(fp,nchunk,MAXLINE,buffer); +read_data.cpp: eof = comm->read_lines_from_file(fp,nchunk,MAXLINE,buffer); +read_data.cpp: eof = comm->read_lines_from_file(fp,nchunk,MAXLINE,buffer); +read_data.cpp: eof = comm->read_lines_from_file(fp,nchunk,MAXLINE,buffer); +read_data.cpp: int eof = comm->read_lines_from_file(fp,ntypes,MAXLINE,buf); +read_data.cpp: int eof = comm->read_lines_from_file(fp,ntypes,MAXLINE,buf); +read_data.cpp: int eof = comm->read_lines_from_file(fp,nsq,MAXLINE,buf); +read_data.cpp: int eof = comm->read_lines_from_file(fp,nbondtypes,MAXLINE,buf); +read_data.cpp: int eof = comm->read_lines_from_file(fp,nangletypes,MAXLINE,buf); +read_data.cpp: int eof = comm->read_lines_from_file(fp,ndihedraltypes,MAXLINE,buf); +read_data.cpp: int eof = comm->read_lines_from_file(fp,nimpropertypes,MAXLINE,buf); +read_data.cpp: eof = comm->read_lines_from_file(fp,nchunk,MAXLINE,buffer); +read_dump.cpp: comm->set_proc_grid(0); +read_restart.cpp: comm->set_proc_grid(); +read_restart.cpp: if (comm->me == 0) +read_restart.cpp: if (nprocs_file != comm->nprocs && me == 0) +read_restart.cpp: comm->nprocs)); +read_restart.cpp: if (comm->user_procgrid[0] != 0 && +read_restart.cpp: procgrid[0] != comm->user_procgrid[0]) flag = 1; +read_restart.cpp: if (comm->user_procgrid[1] != 0 && +read_restart.cpp: procgrid[1] != comm->user_procgrid[1]) flag = 1; +read_restart.cpp: if (comm->user_procgrid[2] != 0 && +read_restart.cpp: procgrid[2] != comm->user_procgrid[2]) flag = 1; +read_restart.cpp: if (comm->me ==0) +read_restart.cpp: comm->mode = read_int(); +read_restart.cpp: comm->cutghostuser = read_double(); +read_restart.cpp: comm->ghost_velocity = read_int(); +read_restart.cpp: if (comm->me ==0) +read_restart.cpp: if (comm->me ==0) +read_restart.cpp: if (comm->me ==0) +read_restart.cpp: if (comm->me ==0) +read_restart.cpp: if (comm->me ==0) +read_restart.cpp: if (comm->me ==0) +region_deprecated.cpp: if (lmp->comm->me == 0) +replicate.cpp: int me = comm->me; +replicate.cpp: int nprocs = comm->nprocs; +replicate.cpp: if (comm->me == 0) +replicate.cpp: comm->set_proc_grid(); +replicate.cpp: if (comm->layout != Comm::LAYOUT_TILED) { +replicate.cpp: if (comm->myloc[0] == 0) sublo[0] -= epsilon[0]; +replicate.cpp: if (comm->myloc[0] == comm->procgrid[0]-1) subhi[0] += epsilon[0]; +replicate.cpp: if (comm->myloc[1] == 0) sublo[1] -= epsilon[1]; +replicate.cpp: if (comm->myloc[1] == comm->procgrid[1]-1) subhi[1] += epsilon[1]; +replicate.cpp: if (comm->myloc[2] == 0) sublo[2] -= epsilon[2]; +replicate.cpp: if (comm->myloc[2] == comm->procgrid[2]-1) subhi[2] += epsilon[2]; +replicate.cpp: if (comm->mysplit[0][0] == 0.0) sublo[0] -= epsilon[0]; +replicate.cpp: if (comm->mysplit[0][1] == 1.0) subhi[0] += epsilon[0]; +replicate.cpp: if (comm->mysplit[1][0] == 0.0) sublo[1] -= epsilon[1]; +replicate.cpp: if (comm->mysplit[1][1] == 1.0) subhi[1] += epsilon[1]; +replicate.cpp: if (comm->mysplit[2][0] == 0.0) sublo[2] -= epsilon[2]; +replicate.cpp: if (comm->mysplit[2][1] == 1.0) subhi[2] += epsilon[2]; +reset_atom_ids.cpp: if (comm->me == 0) utils::logmesg(lmp,"Resetting atom IDs ...\n"); +reset_atom_ids.cpp: // initialize system since comm->borders() will be invoked +reset_atom_ids.cpp: comm->setup(); +reset_atom_ids.cpp: comm->exchange(); +reset_atom_ids.cpp: comm->borders(); +reset_atom_ids.cpp: comm->forward_comm_array(1,newIDs); +reset_atom_ids.cpp: int me = comm->me; +reset_atom_ids.cpp: int nprocs = comm->nprocs; +reset_atom_ids.cpp: int nreturn = comm->rendezvous(1,nlocal,(char *) atombuf,sizeof(AtomRvous), +reset_atom_ids.cpp: // rptr->comm->me,i,n,in[i].ibin,binlo,binhi); +reset_atom_ids.cpp: // pass outbuf of IDRvous datums back to comm->rendezvous +reset_mol_ids.cpp: if (comm->me == 0) utils::logmesg(lmp,"Resetting molecule IDs ...\n"); +reset_mol_ids.cpp: // initialize system since comm->borders() will be invoked +reset_mol_ids.cpp: comm->setup(); +reset_mol_ids.cpp: comm->exchange(); +reset_mol_ids.cpp: comm->borders(); +reset_mol_ids.cpp: if (comm->me == 0) { +respa.cpp: if (comm->me == 0) { +respa.cpp: if (flag && comm->me == 0) +respa.cpp: if (modify->nfix == 0 && comm->me == 0) +respa.cpp: if (comm->me == 0 && screen) { +respa.cpp: comm->setup(); +respa.cpp: comm->exchange(); +respa.cpp: comm->borders(); +respa.cpp: if (newton[ilevel]) comm->reverse_comm(); +respa.cpp: comm->setup(); +respa.cpp: comm->exchange(); +respa.cpp: comm->borders(); +respa.cpp: if (newton[ilevel]) comm->reverse_comm(); +respa.cpp: comm->setup(); +respa.cpp: comm->exchange(); +respa.cpp: comm->borders(); +respa.cpp: comm->forward_comm(); +respa.cpp: comm->forward_comm(); +respa.cpp: comm->reverse_comm(); +set.cpp: if (comm->me == 0) utils::logmesg(lmp,"Setting atom values ...\n"); +set.cpp: if (comm->me == 0) { +set.cpp: RanMars *ranmars = new RanMars(lmp,seed + comm->me); +set.cpp: // init entire system since comm->exchange is done +set.cpp: if (comm->me == 0) utils::logmesg(lmp," system init for set ...\n"); +set.cpp: comm->setup(); +set.cpp: comm->exchange(); +set.cpp: comm->borders(); +special.cpp: comm->rendezvous(RVOUS,nlocal,(char *) idbuf,sizeof(IDRvous),0,proclist, +special.cpp: int nreturn = comm->rendezvous(RVOUS,nsend,(char *) inbuf,sizeof(PairRvous), +special.cpp: int nreturn = comm->rendezvous(RVOUS,nsend,(char *) inbuf,sizeof(PairRvous), +special.cpp: int nreturn = comm->rendezvous(RVOUS,nsend,(char *) inbuf,sizeof(PairRvous), +special.cpp: int nreturn = comm->rendezvous(RVOUS,nsend,(char *) inbuf,sizeof(PairRvous), +special.cpp: int nreturn = comm->rendezvous(RVOUS,nsend,(char *) inbuf,sizeof(PairRvous), +special.cpp: if (comm->me == 0) diff --git a/src/nbin.cpp b/src/nbin.cpp index c28c22882c..1323a9042a 100644 --- a/src/nbin.cpp +++ b/src/nbin.cpp @@ -41,11 +41,8 @@ NBin::NBin(LAMMPS *lmp) : Pointers(lmp) bininvx_multi = nullptr; bininvy_multi = nullptr; bininvz_multi = nullptr; binhead_multi = nullptr; maxbins_multi = nullptr; - - map_type_multi = nullptr; - cutmultisq = nullptr; - maxgroups = 0; + maxcollections = 0; neighbor->last_setup_bins = -1; @@ -85,7 +82,7 @@ NBin::~NBin() memory->destroy(bininvy_multi); memory->destroy(bininvz_multi); - for (int n = 0; n < maxgroups; n++) { + for (int n = 0; n < maxcollections; n++) { memory->destroy(binhead_multi[n]); } delete [] binhead_multi; @@ -115,9 +112,8 @@ void NBin::copy_neighbor_info() bboxlo = neighbor->bboxlo; bboxhi = neighbor->bboxhi; - n_multi_groups = neighbor->n_multi_groups; - map_type_multi = neighbor->map_type_multi; - cutmultisq = neighbor->cutmultisq; + ncollections = neighbor->ncollections; + cutcollectionsq = neighbor->cutcollectionsq; // overwrite Neighbor cutoff with custom value set by requestor // only works for style = BIN (checked by Neighbor class) @@ -175,10 +171,10 @@ int NBin::coord2bin(double *x) /* ---------------------------------------------------------------------- - convert atom coords into local bin # for a particular grouping + convert atom coords into local bin # for a particular collection ------------------------------------------------------------------------- */ -int NBin::coord2bin_multi(double *x, int ig) +int NBin::coord2bin_multi(double *x, int ic) { int ix,iy,iz; int ibin; @@ -187,33 +183,33 @@ int NBin::coord2bin_multi(double *x, int ig) error->one(FLERR,"Non-numeric positions - simulation unstable"); if (x[0] >= bboxhi[0]) - ix = static_cast ((x[0]-bboxhi[0])*bininvx_multi[ig]) + nbinx_multi[ig]; + ix = static_cast ((x[0]-bboxhi[0])*bininvx_multi[ic]) + nbinx_multi[ic]; else if (x[0] >= bboxlo[0]) { - ix = static_cast ((x[0]-bboxlo[0])*bininvx_multi[ig]); - ix = MIN(ix,nbinx_multi[ig]-1); + ix = static_cast ((x[0]-bboxlo[0])*bininvx_multi[ic]); + ix = MIN(ix,nbinx_multi[ic]-1); } else - ix = static_cast ((x[0]-bboxlo[0])*bininvx_multi[ig]) - 1; + ix = static_cast ((x[0]-bboxlo[0])*bininvx_multi[ic]) - 1; if (x[1] >= bboxhi[1]) - iy = static_cast ((x[1]-bboxhi[1])*bininvy_multi[ig]) + nbiny_multi[ig]; + iy = static_cast ((x[1]-bboxhi[1])*bininvy_multi[ic]) + nbiny_multi[ic]; else if (x[1] >= bboxlo[1]) { - iy = static_cast ((x[1]-bboxlo[1])*bininvy_multi[ig]); - iy = MIN(iy,nbiny_multi[ig]-1); + iy = static_cast ((x[1]-bboxlo[1])*bininvy_multi[ic]); + iy = MIN(iy,nbiny_multi[ic]-1); } else - iy = static_cast ((x[1]-bboxlo[1])*bininvy_multi[ig]) - 1; + iy = static_cast ((x[1]-bboxlo[1])*bininvy_multi[ic]) - 1; if (x[2] >= bboxhi[2]) - iz = static_cast ((x[2]-bboxhi[2])*bininvz_multi[ig]) + nbinz_multi[ig]; + iz = static_cast ((x[2]-bboxhi[2])*bininvz_multi[ic]) + nbinz_multi[ic]; else if (x[2] >= bboxlo[2]) { - iz = static_cast ((x[2]-bboxlo[2])*bininvz_multi[ig]); - iz = MIN(iz,nbinz_multi[ig]-1); + iz = static_cast ((x[2]-bboxlo[2])*bininvz_multi[ic]); + iz = MIN(iz,nbinz_multi[ic]-1); } else - iz = static_cast ((x[2]-bboxlo[2])*bininvz_multi[ig]) - 1; + iz = static_cast ((x[2]-bboxlo[2])*bininvz_multi[ic]) - 1; - ibin = (iz-mbinzlo_multi[ig])*mbiny_multi[ig]*mbinx_multi[ig] - + (iy-mbinylo_multi[ig])*mbinx_multi[ig] - + (ix-mbinxlo_multi[ig]); + ibin = (iz-mbinzlo_multi[ic])*mbiny_multi[ic]*mbinx_multi[ic] + + (iy-mbinylo_multi[ic])*mbinx_multi[ic] + + (ix-mbinxlo_multi[ic]); return ibin; } diff --git a/src/nbin.h b/src/nbin.h index e31e5bab99..97299398fa 100644 --- a/src/nbin.h +++ b/src/nbin.h @@ -73,9 +73,8 @@ class NBin : protected Pointers { int binsizeflag; double binsize_user; double *bboxlo,*bboxhi; - int n_multi_groups; - int *map_type_multi; - double **cutmultisq; + int ncollections; + double **cutcollectionsq; // data common to all NBin variants @@ -89,7 +88,7 @@ class NBin : protected Pointers { // data for multi NBin - int maxgroups; // size of multi arrays + int maxcollections; // size of multi arrays int * maxbins_multi; // size of 2nd dimension of binhead_multi array // methods diff --git a/src/nbin_multi.cpp b/src/nbin_multi.cpp index 04756ff803..db1e8d5682 100644 --- a/src/nbin_multi.cpp +++ b/src/nbin_multi.cpp @@ -38,7 +38,7 @@ void NBinMulti::bin_atoms_setup(int nall) { // binhead_multi[n] = per-bin vector mbins in length mbins_multi[n] - for (int n = 0; n < maxgroups; n++) { + for (int n = 0; n < maxcollections; n++) { if (mbins_multi[n] > maxbins_multi[n]) { maxbins_multi[n] = mbins_multi[n]; memory->destroy(binhead_multi[n]); @@ -75,81 +75,80 @@ void NBinMulti::bin_atoms_setup(int nall) for triclinic boxes: tilted simulation box cannot contain integer # of bins stencil & neigh list built differently to account for this - mbinlo_multi = lowest global bin any of my ghost atoms could fall into for each grouping - mbinhi_multi = highest global bin any of my ghost atoms could fall into for each grouping - mbin_multi = number of bins I need in a dimension for each grouping + mbinlo_multi = lowest global bin any of my ghost atoms could fall into for each collection + mbinhi_multi = highest global bin any of my ghost atoms could fall into for each collection + mbin_multi = number of bins I need in a dimension for each collection ------------------------------------------------------------------------- */ -void NBinMulti::setup_bins(int style) +void NBinMulti::setup_bins(int /*style*/) { int n; // Initialize arrays - if (n_multi_groups > maxgroups) { + if (ncollections > maxcollections) { // Clear any/all memory for existing groupings - for (n = 0; n < maxgroups; n++) + for (n = 0; n < maxcollections; n++) memory->destroy(binhead_multi[n]); delete [] binhead_multi; - // Realloacte at updated maxgroups - maxgroups = n_multi_groups; + maxcollections = ncollections; - binhead_multi = new int*[maxgroups](); + binhead_multi = new int*[maxcollections](); memory->destroy(nbinx_multi); memory->destroy(nbiny_multi); memory->destroy(nbinz_multi); - memory->create(nbinx_multi, maxgroups, "neigh:nbinx_multi"); - memory->create(nbiny_multi, maxgroups, "neigh:nbiny_multi"); - memory->create(nbinz_multi, maxgroups, "neigh:nbinz_multi"); + memory->create(nbinx_multi, maxcollections, "neigh:nbinx_multi"); + memory->create(nbiny_multi, maxcollections, "neigh:nbiny_multi"); + memory->create(nbinz_multi, maxcollections, "neigh:nbinz_multi"); memory->destroy(mbins_multi); memory->destroy(mbinx_multi); memory->destroy(mbiny_multi); memory->destroy(mbinz_multi); - memory->create(mbins_multi, maxgroups, "neigh:mbins_multi"); - memory->create(mbinx_multi, maxgroups, "neigh:mbinx_multi"); - memory->create(mbiny_multi, maxgroups, "neigh:mbiny_multi"); - memory->create(mbinz_multi, maxgroups, "neigh:mbinz_multi"); + memory->create(mbins_multi, maxcollections, "neigh:mbins_multi"); + memory->create(mbinx_multi, maxcollections, "neigh:mbinx_multi"); + memory->create(mbiny_multi, maxcollections, "neigh:mbiny_multi"); + memory->create(mbinz_multi, maxcollections, "neigh:mbinz_multi"); memory->destroy(mbinxlo_multi); memory->destroy(mbinylo_multi); memory->destroy(mbinzlo_multi); - memory->create(mbinxlo_multi, maxgroups, "neigh:mbinxlo_multi"); - memory->create(mbinylo_multi, maxgroups, "neigh:mbinylo_multi"); - memory->create(mbinzlo_multi, maxgroups, "neigh:mbinzlo_multi"); + memory->create(mbinxlo_multi, maxcollections, "neigh:mbinxlo_multi"); + memory->create(mbinylo_multi, maxcollections, "neigh:mbinylo_multi"); + memory->create(mbinzlo_multi, maxcollections, "neigh:mbinzlo_multi"); memory->destroy(binsizex_multi); memory->destroy(binsizey_multi); memory->destroy(binsizez_multi); - memory->create(binsizex_multi, maxgroups, "neigh:binsizex_multi"); - memory->create(binsizey_multi, maxgroups, "neigh:binsizey_multi"); - memory->create(binsizez_multi, maxgroups, "neigh:binsizez_multi"); + memory->create(binsizex_multi, maxcollections, "neigh:binsizex_multi"); + memory->create(binsizey_multi, maxcollections, "neigh:binsizey_multi"); + memory->create(binsizez_multi, maxcollections, "neigh:binsizez_multi"); memory->destroy(bininvx_multi); memory->destroy(bininvy_multi); memory->destroy(bininvz_multi); - memory->create(bininvx_multi, maxgroups, "neigh:bininvx_multi"); - memory->create(bininvy_multi, maxgroups, "neigh:bininvy_multi"); - memory->create(bininvz_multi, maxgroups, "neigh:bininvz_multi"); + memory->create(bininvx_multi, maxcollections, "neigh:bininvx_multi"); + memory->create(bininvy_multi, maxcollections, "neigh:bininvy_multi"); + memory->create(bininvz_multi, maxcollections, "neigh:bininvz_multi"); memory->destroy(maxbins_multi); - memory->create(maxbins_multi, maxgroups, "neigh:maxbins_multi"); + memory->create(maxbins_multi, maxcollections, "neigh:maxbins_multi"); // ensure reallocation occurs in bin_atoms_setup() - for (n = 0; n < maxgroups; n++) { + for (n = 0; n < maxcollections; n++) { maxbins_multi[n] = 0; } maxatom = 0; } - // Identify smallest group - int igroupmin = 0; - for (n = 0; n < maxgroups; n++) - if (cutmultisq[n][n] < cutmultisq[igroupmin][igroupmin]) - igroupmin = n; + // Identify smallest collection + int icollectionmin = 0; + for (n = 0; n < maxcollections; n++) + if (cutcollectionsq[n][n] < cutcollectionsq[icollectionmin][icollectionmin]) + icollectionmin = n; // bbox = size of bbox of entire domain // bsubbox lo/hi = bounding box of my subdomain extended by comm->cutghost @@ -189,13 +188,13 @@ void NBinMulti::setup_bins(int style) double binsize_optimal, binsizeinv, coord; int mbinxhi,mbinyhi,mbinzhi; - for (n = 0; n < maxgroups; n++) { - // binsize_user only relates to smallest group - // optimal bin size is roughly 1/2 the group-group cutoff + for (n = 0; n < maxcollections; n++) { + // binsize_user only relates to smallest collection + // optimal bin size is roughly 1/2 the collection-collection cutoff // special case of all cutoffs = 0.0, binsize = box size - if (n == igroupmin && binsizeflag) binsize_optimal = binsize_user; - else binsize_optimal = 0.5*sqrt(cutmultisq[n][n]); + if (n == icollectionmin && binsizeflag) binsize_optimal = binsize_user; + else binsize_optimal = 0.5*sqrt(cutcollectionsq[n][n]); if (binsize_optimal == 0.0) binsize_optimal = bbox[0]; binsizeinv = 1.0/binsize_optimal; @@ -296,16 +295,15 @@ void NBinMulti::bin_atoms() int i,ibin,n; last_bin = update->ntimestep; - for (n = 0; n < maxgroups; n++) { + for (n = 0; n < maxcollections; n++) { for (i = 0; i < mbins_multi[n]; i++) binhead_multi[n][i] = -1; } - // bin in reverse order so linked list will be in forward order // also puts ghost atoms at end of list, which is necessary + int *collection = neighbor->collection; double **x = atom->x; int *mask = atom->mask; - int *type = atom->type; int nlocal = atom->nlocal; int nall = nlocal + atom->nghost; @@ -313,7 +311,7 @@ void NBinMulti::bin_atoms() int bitmask = group->bitmask[includegroup]; for (i = nall-1; i >= nlocal; i--) { if (mask[i] & bitmask) { - n = map_type_multi[type[i]]; + n = collection[i]; ibin = coord2bin_multi(x[i], n); atom2bin[i] = ibin; bins[i] = binhead_multi[n][ibin]; @@ -321,15 +319,15 @@ void NBinMulti::bin_atoms() } } for (i = atom->nfirst-1; i >= 0; i--) { - n = map_type_multi[type[i]]; + n = collection[i]; ibin = coord2bin_multi(x[i], n); atom2bin[i] = ibin; bins[i] = binhead_multi[n][ibin]; binhead_multi[n][ibin] = i; } } else { - for (i = nall-1; i >= 0; i--) { - n = map_type_multi[type[i]]; + for (i = nall-1; i >= 0; i--) { + n = collection[i]; ibin = coord2bin_multi(x[i], n); atom2bin[i] = ibin; bins[i] = binhead_multi[n][ibin]; @@ -343,7 +341,7 @@ void NBinMulti::bin_atoms() double NBinMulti::memory_usage() { double bytes = 0; - for (int m = 0; m < maxgroups; m++) + for (int m = 0; m < maxcollections; m++) bytes += maxbins_multi[m]*sizeof(int); bytes += 2*maxatom*sizeof(int); return bytes; diff --git a/src/neighbor.cpp b/src/neighbor.cpp index 5a1255edbe..16aca434f1 100644 --- a/src/neighbor.cpp +++ b/src/neighbor.cpp @@ -43,6 +43,7 @@ #include "style_npair.h" #include "style_nstencil.h" #include "style_ntopo.h" +#include "tokenizer.h" #include "update.h" #include @@ -53,6 +54,7 @@ using namespace NeighConst; #define RQDELTA 1 #define EXDELTA 1 +#define DELTA_PERATOM 64 #define BIG 1.0e20 @@ -195,9 +197,13 @@ pairclass(nullptr), pairnames(nullptr), pairmasks(nullptr) // Multi data - map_type_multi = nullptr; - cutmultisq = nullptr; - multi_groups = 0; + type2collection = nullptr; + collection2cut = nullptr; + collection = nullptr; + cutcollectionsq = nullptr; + custom_collection_flag = 0; + interval_collection_flag = 0; + nmax_collection = 0; // Kokkos setting @@ -265,8 +271,10 @@ Neighbor::~Neighbor() delete [] ex_mol_bit; memory->destroy(ex_mol_intra); - memory->destroy(map_type_multi); - memory->destroy(cutmultisq); + memory->destroy(type2collection); + memory->destroy(collection2cut); + memory->destroy(collection); + memory->destroy(cutcollectionsq); } /* ---------------------------------------------------------------------- */ @@ -350,33 +358,96 @@ void Neighbor::init() } cutneighmaxsq = cutneighmax * cutneighmax; - // multi cutoffs + // Define cutoffs for multi if(style == Neighbor::MULTI){ - int igroup, jgroup; + int icollection, jcollection; - // If not defined from custom grouping, create default map - if(not map_type_multi) { - n_multi_groups = n; - memory->create(map_type_multi,n+1,"neigh:map_type_multi"); + // If collections not yet defined, create default map using types + if(not custom_collection_flag) { + ncollections = n; + interval_collection_flag = 0; + memory->create(type2collection,n+1,"neigh:type2collection"); for(i = 1; i <= n; i++) - map_type_multi[i] = i-1; + type2collection[i] = i-1; } - // Define maximum interaction distance for each pair of groups - memory->grow(cutmultisq, n_multi_groups, n_multi_groups, "neigh:cutmultisq"); - for(i = 0; i < n_multi_groups; i++) - for(j = 0; j < n_multi_groups; j++) - cutmultisq[i][j] = 0.0; + memory->grow(cutcollectionsq, ncollections, ncollections, "neigh:cutcollectionsq"); - for(i = 1; i <= n; i++){ - igroup = map_type_multi[i]; - for(j = 1; j <= n; j++){ - jgroup = map_type_multi[j]; - if(cutneighsq[i][j] > cutmultisq[igroup][jgroup]) { - cutmultisq[igroup][jgroup] = cutneighsq[i][j]; - cutmultisq[jgroup][igroup] = cutneighsq[i][j]; + // 3 possible ways of defining collections + // 1) Types are used to define collections + // Each collection loops through its owned types, and uses cutneighsq to calculate its cutoff + // 2) Collections are defined by intervals, point particles + // Types are first sorted into collections based on cutneighsq[i][i] + // Each collection loops through its owned types, and uses cutneighsq to calculate its cutoff + // 3) Collections are defined by intervals, finite particles + // + + // Define collection cutoffs + for(i = 0; i < ncollections; i++) + for(j = 0; j < ncollections; j++) + cutcollectionsq[i][j] = 0.0; + + if(not interval_collection_flag){ + finite_cut_flag = 0; + for(i = 1; i <= n; i++){ + icollection = type2collection[i]; + for(j = 1; j <= n; j++){ + jcollection = type2collection[j]; + if(cutneighsq[i][j] > cutcollectionsq[icollection][jcollection]) { + cutcollectionsq[icollection][jcollection] = cutneighsq[i][j]; + cutcollectionsq[jcollection][icollection] = cutneighsq[i][j]; + } } } + } else { + if(force->pair->finitecutflag){ + finite_cut_flag = 1; + // If cutoffs depend on finite atom sizes, use radii of intervals to find cutoffs + double ri, rj, tmp; + for(i = 0; i < ncollections; i++){ + ri = collection2cut[i]*0.5; + for(j = 0; j < ncollections; j++){ + rj = collection2cut[j]*0.5; + tmp = force->pair->radii2cut(ri, rj); + cutcollectionsq[i][j] = tmp*tmp; + } + } + } else { + finite_cut_flag = 0; + + // Map types to collections + if(not type2collection) + memory->create(type2collection,n+1,"neigh:type2collection"); + + for(i = 1; i <= n; i++) + type2collection[i] = -1; + + double cuttmp; + for(i = 1; i <= n; i++){ + cuttmp = sqrt(cutneighsq[i][i]); + for(icollection = 0; icollection < ncollections; icollection ++){ + if(collection2cut[icollection] > cuttmp) { + type2collection[i] = icollection; + break; + } + } + + if(type2collection[i] == -1) + error->all(FLERR, "Pair cutoff exceeds interval cutoffs for multi"); + } + + // Define cutoffs + for(i = 1; i <= n; i++){ + icollection = type2collection[i]; + for(j = 1; j <= n; j++){ + jcollection = type2collection[j]; + if(cutneighsq[i][j] > cutcollectionsq[icollection][jcollection]) { + cutcollectionsq[icollection][jcollection] = cutneighsq[i][j]; + cutcollectionsq[jcollection][icollection] = cutneighsq[i][j]; + } + } + } + } } } @@ -2100,6 +2171,8 @@ void Neighbor::build(int topoflag) int nlocal = atom->nlocal; int nall = nlocal + atom->nghost; + // rebuild collection array from scratch + if(style == Neighbor::MULTI) build_collection(0); // check that using special bond flags will not overflow neigh lists @@ -2140,7 +2213,7 @@ void Neighbor::build(int topoflag) } } } - + // bin atoms for all NBin instances // not just NBin associated with perpetual lists, also occasional lists // b/c cannot wait to bin occasional lists in build_one() call @@ -2423,53 +2496,91 @@ void Neighbor::modify_params(int narg, char **arg) iarg += 2; } else error->all(FLERR,"Illegal neigh_modify command"); - } else if (strcmp(arg[iarg],"multi/custom") == 0) { + } else if (strcmp(arg[iarg],"collection/interval") == 0) { if(style != Neighbor::MULTI) - error->all(FLERR,"Cannot use multi/custom command without multi setting"); + error->all(FLERR,"Cannot use collection/interval command without multi setting"); if(iarg+2 > narg) - error->all(FLERR,"Invalid multi/custom command"); - int nextra = utils::inumeric(FLERR,arg[iarg+1],false,lmp); - if(iarg+1+nextra > narg) - error->all(FLERR,"Invalid multi/custom command"); + error->all(FLERR,"Invalid collection/interval command"); + ncollections = utils::inumeric(FLERR,arg[iarg+1],false,lmp); + if(iarg+1+ncollections > narg) + error->all(FLERR,"Invalid collection/interval command"); int ntypes = atom->ntypes; - int n, nlo, nhi, i, j, igroup, jgroup; + int n, nlo, nhi, i, j; + + interval_collection_flag = 1; + custom_collection_flag = 1; + if(not collection2cut) + memory->create(collection2cut,ncollections,"neigh:collection2cut"); - if(not map_type_multi) - memory->create(map_type_multi,ntypes+1,"neigh:map_type_multi"); + // Set upper cutoff for each collection + char *id; + double cut_interval; + for(i = 0; i < ncollections; i++){ + cut_interval = utils::numeric(FLERR,arg[iarg+2+i],false,lmp); + collection2cut[i] = cut_interval; + + if(i != 0) + if(collection2cut[i-1] >= collection2cut[i]) + error->all(FLERR,"Nonsequential interval cutoffs in collection/interval setting"); + } + + iarg += 2 + ncollections; + } else if (strcmp(arg[iarg],"collection/type") == 0) { + if(style != Neighbor::MULTI) + error->all(FLERR,"Cannot use collection/type command without multi setting"); + + if(iarg+2 > narg) + error->all(FLERR,"Invalid collection/type command"); + ncollections = utils::inumeric(FLERR,arg[iarg+1],false,lmp); + if(iarg+1+ncollections > narg) + error->all(FLERR,"Invalid collection/type command"); + + int ntypes = atom->ntypes; + int n, nlo, nhi, i, j, k; + + interval_collection_flag = 0; + custom_collection_flag = 1; + if(not type2collection) + memory->create(type2collection,ntypes+1,"neigh:type2collection"); // Erase previous mapping for(i = 1; i <= ntypes; i++) - map_type_multi[i] = -1; + type2collection[i] = -1; // For each custom range, define mapping for types in interval - n_multi_groups = 0; - char *id; - for(i = 0; i < nextra; i++){ + int nfield; + char *str; + for(i = 0; i < ncollections; i++){ n = strlen(arg[iarg+2+i]) + 1; - id = new char[n]; - strcpy(id,arg[iarg+2+i]); - utils::bounds(FLERR,id,1,ntypes,nlo,nhi,error); - delete [] id; + str = new char[n]; + strcpy(str,arg[iarg+2+i]); + std::vector words = Tokenizer(str, ",").as_vector(); + nfield = words.size(); + + for (j = 0; j < nfield; j++) { + const char * field = words[j].c_str(); + utils::bounds(FLERR,field,1,ntypes,nlo,nhi,error); - for (j = nlo; j <= nhi; j++) { - if(map_type_multi[j] != -1) - error->all(FLERR,"Type specified more than once in multi/custom commnd"); - map_type_multi[j] = n_multi_groups; + for (k = nlo; k <= nhi; k++) { + if(type2collection[k] != -1) + error->all(FLERR,"Type specified more than once in collection/type commnd"); + type2collection[k] = i; + } } - n_multi_groups += 1; + + delete [] str; } - // Create separate group for each undefined atom type + // Check for undefined atom type for(i = 1; i <= ntypes; i++){ - if(map_type_multi[i] == -1){ - map_type_multi[i] = n_multi_groups; - n_multi_groups += 1; + if(type2collection[i] == -1){ + error->all(FLERR,"Type missing in collection/type commnd"); } } - iarg += 2 + nextra; + iarg += 2 + ncollections; } else error->all(FLERR,"Illegal neigh_modify command"); } } @@ -2521,6 +2632,46 @@ int Neighbor::any_full() return any_full; } +/* ---------------------------------------------------------------------- + populate collection array for multi starting at the index istart +------------------------------------------------------------------------- */ + +void Neighbor::build_collection(int istart) +{ + if (style != Neighbor::MULTI) + error->all(FLERR, "Cannot define atom collections without neighbor style multi"); + + int nmax = atom->nlocal+atom->nghost; + if(nmax > nmax_collection){ + nmax_collection = nmax+DELTA_PERATOM; + memory->grow(collection, nmax_collection, "neigh:collection"); + } + + if(finite_cut_flag){ + double cut; + int icollection; + for(int i = istart; i < nmax; i++){ + cut = force->pair->atom2cut(i); + collection[i] = -1; + + for(icollection = 0; icollection < ncollections; icollection++){ + if(collection2cut[icollection] > cut) { + collection[i] = icollection; + break; + } + } + + if(collection[i] == -1) + error->one(FLERR, "Atom cutoff exceeds interval cutoffs for multi"); + } + } else { + int *type = atom->type; + for(int i = istart; i < nmax; i++){ + collection[i] = type2collection[type[i]]; + } + } +} + /* ---------------------------------------------------------------------- return # of bytes of allocated memory ------------------------------------------------------------------------- */ diff --git a/src/neighbor.h b/src/neighbor.h index b8e1189b35..7b6246dec7 100644 --- a/src/neighbor.h +++ b/src/neighbor.h @@ -104,11 +104,16 @@ class Neighbor : protected Pointers { // optional type grouping for multi - int multi_groups; // 1 if custom groupings are defined for multi - int n_multi_groups; // # of custom groupings - int *map_type_multi; // ntype array mapping types to custom groupings - double **cutmultisq; // cutoffs for each combination of custom groupings - + int custom_collection_flag; // 1 if custom collections are defined for multi + int interval_collection_flag; // 1 if custom collections use intervals + int finite_cut_flag; // 1 if multi considers finite atom size + int ncollections; // # of custom collections + int nmax_collection; // maximum atoms stored in collection array + int *type2collection; // ntype array mapping types to custom collections + double *collection2cut; // ncollection array with upper bounds on cutoff intervals + double **cutcollectionsq; // cutoffs for each combination of collections + int *collection; // local per-atom array to store collection id + // public methods Neighbor(class LAMMPS *); @@ -130,6 +135,7 @@ class Neighbor : protected Pointers { int exclude_setting(); // return exclude value to accelerator pkg class NeighRequest *find_request(void *); // find a neighbor request int any_full(); // Check if any old requests had full neighbor lists + void build_collection(int); // build peratom collection array starting at the given index double memory_usage(); diff --git a/src/npair.cpp b/src/npair.cpp index 95a61fbefa..eee8610958 100644 --- a/src/npair.cpp +++ b/src/npair.cpp @@ -97,9 +97,8 @@ void NPair::copy_neighbor_info() // multi info - n_multi_groups = neighbor->n_multi_groups; - map_type_multi = neighbor->map_type_multi; - cutmultisq = neighbor->cutmultisq; + ncollections = neighbor->ncollections; + cutcollectionsq = neighbor->cutcollectionsq; // overwrite per-type Neighbor cutoffs with custom value set by requestor // only works for style = BIN (checked by Neighbor class) @@ -184,9 +183,8 @@ void NPair::build_setup() { if (nb) copy_bin_info(); if (ns) copy_stencil_info(); - + // set here, since build_setup() always called before build() - last_build = update->ntimestep; } @@ -269,10 +267,10 @@ int NPair::coord2bin(double *x, int &ix, int &iy, int &iz) /* ---------------------------------------------------------------------- - multi version of coord2bin for a given grouping + multi version of coord2bin for a given collection ------------------------------------------------------------------------- */ -int NPair::coord2bin(double *x, int ig) +int NPair::coord2bin(double *x, int ic) { int ix,iy,iz; int ibin; @@ -281,32 +279,32 @@ int NPair::coord2bin(double *x, int ig) error->one(FLERR,"Non-numeric positions - simulation unstable"); if (x[0] >= bboxhi[0]) - ix = static_cast ((x[0]-bboxhi[0])*bininvx_multi[ig]) + nbinx_multi[ig]; + ix = static_cast ((x[0]-bboxhi[0])*bininvx_multi[ic]) + nbinx_multi[ic]; else if (x[0] >= bboxlo[0]) { - ix = static_cast ((x[0]-bboxlo[0])*bininvx_multi[ig]); - ix = MIN(ix,nbinx_multi[ig]-1); + ix = static_cast ((x[0]-bboxlo[0])*bininvx_multi[ic]); + ix = MIN(ix,nbinx_multi[ic]-1); } else - ix = static_cast ((x[0]-bboxlo[0])*bininvx_multi[ig]) - 1; + ix = static_cast ((x[0]-bboxlo[0])*bininvx_multi[ic]) - 1; if (x[1] >= bboxhi[1]) - iy = static_cast ((x[1]-bboxhi[1])*bininvy_multi[ig]) + nbiny_multi[ig]; + iy = static_cast ((x[1]-bboxhi[1])*bininvy_multi[ic]) + nbiny_multi[ic]; else if (x[1] >= bboxlo[1]) { - iy = static_cast ((x[1]-bboxlo[1])*bininvy_multi[ig]); - iy = MIN(iy,nbiny_multi[ig]-1); + iy = static_cast ((x[1]-bboxlo[1])*bininvy_multi[ic]); + iy = MIN(iy,nbiny_multi[ic]-1); } else - iy = static_cast ((x[1]-bboxlo[1])*bininvy_multi[ig]) - 1; + iy = static_cast ((x[1]-bboxlo[1])*bininvy_multi[ic]) - 1; if (x[2] >= bboxhi[2]) - iz = static_cast ((x[2]-bboxhi[2])*bininvz_multi[ig]) + nbinz_multi[ig]; + iz = static_cast ((x[2]-bboxhi[2])*bininvz_multi[ic]) + nbinz_multi[ic]; else if (x[2] >= bboxlo[2]) { - iz = static_cast ((x[2]-bboxlo[2])*bininvz_multi[ig]); - iz = MIN(iz,nbinz_multi[ig]-1); + iz = static_cast ((x[2]-bboxlo[2])*bininvz_multi[ic]); + iz = MIN(iz,nbinz_multi[ic]-1); } else - iz = static_cast ((x[2]-bboxlo[2])*bininvz_multi[ig]) - 1; + iz = static_cast ((x[2]-bboxlo[2])*bininvz_multi[ic]) - 1; - ix -= mbinxlo_multi[ig]; - iy -= mbinylo_multi[ig]; - iz -= mbinzlo_multi[ig]; - ibin = iz*mbiny_multi[ig]*mbinx_multi[ig] + iy*mbinx_multi[ig] + ix; + ix -= mbinxlo_multi[ic]; + iy -= mbinylo_multi[ic]; + iz -= mbinzlo_multi[ic]; + ibin = iz*mbiny_multi[ic]*mbinx_multi[ic] + iy*mbinx_multi[ic] + ix; return ibin; } diff --git a/src/npair.h b/src/npair.h index c857ec0924..b8416208ca 100644 --- a/src/npair.h +++ b/src/npair.h @@ -48,9 +48,8 @@ class NPair : protected Pointers { double cut_middle_sq; double cut_middle_inside_sq; double *bboxlo,*bboxhi; - int n_multi_groups; - int *map_type_multi; - double **cutmultisq; + int ncollections; + double **cutcollectionsq; // exclusion data from Neighbor class diff --git a/src/npair_full_multi.cpp b/src/npair_full_multi.cpp index c78f7b2d1d..0ce717f722 100644 --- a/src/npair_full_multi.cpp +++ b/src/npair_full_multi.cpp @@ -12,6 +12,7 @@ ------------------------------------------------------------------------- */ #include "npair_full_multi.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -28,18 +29,19 @@ NPairFullMulti::NPairFullMulti(LAMMPS *lmp) : NPair(lmp) {} /* ---------------------------------------------------------------------- binned neighbor list construction for all neighbors - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent every neighbor pair appears in list of both atoms i and j ------------------------------------------------------------------------- */ void NPairFullMulti::build(NeighList *list) { - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,which,ns,imol,iatom,moltemplate; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,which,ns,imol,iatom,moltemplate; tagint tagprev; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; int *neighptr,*s; int js; + int *collection = neighbor->collection; double **x = atom->x; int *type = atom->type; int *mask = atom->mask; @@ -68,7 +70,7 @@ void NPairFullMulti::build(NeighList *list) n = 0; neighptr = ipage->vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -80,22 +82,22 @@ void NPairFullMulti::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + 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 group combinations + // use full stencil for all collection combinations - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >= 0; j = bins[j]) { if (i == j) continue; diff --git a/src/npair_half_multi_newtoff.cpp b/src/npair_half_multi_newtoff.cpp index 11bab91e6a..5f050429ac 100755 --- a/src/npair_half_multi_newtoff.cpp +++ b/src/npair_half_multi_newtoff.cpp @@ -12,6 +12,7 @@ ------------------------------------------------------------------------- */ #include "npair_half_multi_newtoff.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -28,7 +29,7 @@ NPairHalfMultiNewtoff::NPairHalfMultiNewtoff(LAMMPS *lmp) : NPair(lmp) {} /* ---------------------------------------------------------------------- binned neighbor list construction with partial Newton's 3rd law - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks own bin and other bins in stencil pair stored once if i,j are both owned and i < j pair stored by me if j is ghost (also stored by proc owning j) @@ -36,12 +37,13 @@ NPairHalfMultiNewtoff::NPairHalfMultiNewtoff(LAMMPS *lmp) : NPair(lmp) {} void NPairHalfMultiNewtoff::build(NeighList *list) { - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,which,ns,imol,iatom,moltemplate; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,which,ns,imol,iatom,moltemplate; tagint tagprev; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; int *neighptr,*s; int js; + int *collection = neighbor->collection; double **x = atom->x; int *type = atom->type; int *mask = atom->mask; @@ -70,7 +72,7 @@ void NPairHalfMultiNewtoff::build(NeighList *list) n = 0; neighptr = ipage->vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -82,24 +84,24 @@ void NPairHalfMultiNewtoff::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // loop over all atoms in other bins in stencil including self // only store pair if i < j // stores own/own pairs only once // stores own/ghost pairs on both procs - // use full stencil for all group combinations + // use full stencil for all collection combinations - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >=0; j = bins[j]) { if (j <= i) continue; diff --git a/src/npair_half_multi_newton.cpp b/src/npair_half_multi_newton.cpp index aa8f606242..18e66d4adb 100755 --- a/src/npair_half_multi_newton.cpp +++ b/src/npair_half_multi_newton.cpp @@ -12,6 +12,7 @@ ------------------------------------------------------------------------- */ #include "npair_half_multi_newton.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -28,19 +29,20 @@ NPairHalfMultiNewton::NPairHalfMultiNewton(LAMMPS *lmp) : NPair(lmp) {} /* ---------------------------------------------------------------------- binned neighbor list construction with full Newton's 3rd law - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks its own bin and other bins in Newton stencil every pair stored exactly once by some processor ------------------------------------------------------------------------- */ void NPairHalfMultiNewton::build(NeighList *list) { - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,which,ns,imol,iatom,moltemplate; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,which,ns,imol,iatom,moltemplate; tagint tagprev; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; int *neighptr,*s; int js; + int *collection = neighbor->collection; double **x = atom->x; int *type = atom->type; int *mask = atom->mask; @@ -69,7 +71,7 @@ void NPairHalfMultiNewton::build(NeighList *list) n = 0; neighptr = ipage->vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -81,29 +83,29 @@ void NPairHalfMultiNewton::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // if same size: uses half stencil so check central bin - if(cutmultisq[igroup][igroup] == cutmultisq[jgroup][jgroup]){ + if(cutcollectionsq[icollection][icollection] == cutcollectionsq[jcollection][jcollection]){ - if(igroup == jgroup) js = bins[i]; - else js = binhead_multi[jgroup][jbin]; + if(icollection == jcollection) js = bins[i]; + else js = binhead_multi[jcollection][jbin]; - // if same group, + // 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 different groups, + // 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 for (j = js; j >= 0; j = bins[j]) { - if(igroup != jgroup and j < i) continue; + if(icollection != jcollection and j < i) continue; if (j >= nlocal) { if (x[j][2] < ztmp) continue; @@ -139,16 +141,16 @@ void NPairHalfMultiNewton::build(NeighList *list) } } - // for all groups, loop over all atoms in other bins in stencil, store every pair + // for all collections, loop over all atoms in other bins in stencil, store every pair // stencil is empty if i larger than j // stencil is half if i same size as j // stencil is full if i smaller than j - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >= 0; j = bins[j]) { jtype = type[j]; diff --git a/src/npair_half_multi_newton_tri.cpp b/src/npair_half_multi_newton_tri.cpp index 6969c9cfab..90c6cf714d 100755 --- a/src/npair_half_multi_newton_tri.cpp +++ b/src/npair_half_multi_newton_tri.cpp @@ -12,6 +12,7 @@ ------------------------------------------------------------------------- */ #include "npair_half_multi_newton_tri.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -28,19 +29,20 @@ NPairHalfMultiNewtonTri::NPairHalfMultiNewtonTri(LAMMPS *lmp) : NPair(lmp) {} /* ---------------------------------------------------------------------- binned neighbor list construction with Newton's 3rd law for triclinic - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks its own bin and other bins in triclinic stencil every pair stored exactly once by some processor ------------------------------------------------------------------------- */ void NPairHalfMultiNewtonTri::build(NeighList *list) { - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,which,ns,imol,iatom,moltemplate; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,which,ns,imol,iatom,moltemplate; tagint tagprev; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; int *neighptr,*s; int js; - + + int *collection = neighbor->collection; double **x = atom->x; int *type = atom->type; int *mask = atom->mask; @@ -69,7 +71,7 @@ void NPairHalfMultiNewtonTri::build(NeighList *list) n = 0; neighptr = ipage->vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -81,12 +83,12 @@ void NPairHalfMultiNewtonTri::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // loop over all atoms in bins in stencil // stencil is empty if i larger than j @@ -97,15 +99,15 @@ void NPairHalfMultiNewtonTri::build(NeighList *list) // (equal zyx and j <= i) // latter excludes self-self interaction but allows superposed atoms - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >= 0; j = bins[j]) { - // if same size (same group), use half stencil - if(cutmultisq[igroup][igroup] == cutmultisq[jgroup][jgroup]){ + // if same size (same collection), use half stencil + if(cutcollectionsq[icollection][icollection] == cutcollectionsq[jcollection][jcollection]){ if (x[j][2] < ztmp) continue; if (x[j][2] == ztmp) { if (x[j][1] < ytmp) continue; diff --git a/src/npair_half_size_multi_newtoff.cpp b/src/npair_half_size_multi_newtoff.cpp index e361bbce76..77aa4f746b 100644 --- a/src/npair_half_size_multi_newtoff.cpp +++ b/src/npair_half_size_multi_newtoff.cpp @@ -13,6 +13,7 @@ es certain rights in this software. This software is distributed under #include #include "npair_half_size_multi_newtoff.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -28,7 +29,7 @@ NPairHalfSizeMultiNewtoff::NPairHalfSizeMultiNewtoff(LAMMPS *lmp) : NPair(lmp) { /* ---------------------------------------------------------------------- size particles binned neighbor list construction with partial Newton's 3rd law - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks own bin and other bins in stencil pair stored once if i,j are both owned and i < j pair stored by me if j is ghost (also stored by proc owning j) @@ -36,12 +37,13 @@ NPairHalfSizeMultiNewtoff::NPairHalfSizeMultiNewtoff(LAMMPS *lmp) : NPair(lmp) { void NPairHalfSizeMultiNewtoff::build(NeighList *list) { - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,ns; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,ns; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; double radi,radsum,cutdistsq; int *neighptr,*s; int js; + int *collection = neighbor->collection; double **x = atom->x; double *radius = atom->radius; int *type = atom->type; @@ -65,7 +67,7 @@ void NPairHalfSizeMultiNewtoff::build(NeighList *list) n = 0; neighptr = ipage->vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -73,24 +75,24 @@ void NPairHalfSizeMultiNewtoff::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // loop over all atoms in other bins in stencil including self // only store pair if i < j // stores own/own pairs only once // stores own/ghost pairs on both procs - // use full stencil for all group combinations + // use full stencil for all collection combinations - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >=0; j = bins[j]) { if (j <= i) continue; diff --git a/src/npair_half_size_multi_newton.cpp b/src/npair_half_size_multi_newton.cpp index 84571a3b02..67715ebe68 100755 --- a/src/npair_half_size_multi_newton.cpp +++ b/src/npair_half_size_multi_newton.cpp @@ -12,6 +12,7 @@ ------------------------------------------------------------------------- */ #include "npair_half_size_multi_newton.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -27,18 +28,19 @@ NPairHalfSizeMultiNewton::NPairHalfSizeMultiNewton(LAMMPS *lmp) : NPair(lmp) {} /* ---------------------------------------------------------------------- size particles binned neighbor list construction with full Newton's 3rd law - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks its own bin and other bins in Newton stencil every pair stored exactly once by some processor ------------------------------------------------------------------------- */ void NPairHalfSizeMultiNewton::build(NeighList *list) { - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,ns,js; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,ns,js; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; double radi,radsum,cutdistsq; int *neighptr,*s; + int *collection = neighbor->collection; double **x = atom->x; double *radius = atom->radius; int *type = atom->type; @@ -62,7 +64,7 @@ void NPairHalfSizeMultiNewton::build(NeighList *list) n = 0; neighptr = ipage->vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -70,29 +72,29 @@ void NPairHalfSizeMultiNewton::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // if same size: uses half stencil so check central bin - if(cutmultisq[igroup][igroup] == cutmultisq[jgroup][jgroup]){ + if(cutcollectionsq[icollection][icollection] == cutcollectionsq[jcollection][jcollection]){ - if(igroup == jgroup) js = bins[i]; - else js = binhead_multi[jgroup][jbin]; + if(icollection == jcollection) js = bins[i]; + else js = binhead_multi[jcollection][jbin]; - // if same group, + // 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 different groups, + // 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 for (j = js; j >= 0; j = bins[j]) { - if(igroup != jgroup and j < i) continue; + if(icollection != jcollection and j < i) continue; if (j >= nlocal) { if (x[j][2] < ztmp) continue; @@ -121,16 +123,16 @@ void NPairHalfSizeMultiNewton::build(NeighList *list) } } - // for all groups, loop over all atoms in other bins in stencil, store every pair + // for all collections, loop over all atoms in other bins in stencil, store every pair // stencil is empty if i larger than j // stencil is half if i same size as j // stencil is full if i smaller than j - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >= 0; j = bins[j]) { jtype = type[j]; diff --git a/src/npair_half_size_multi_newton_tri.cpp b/src/npair_half_size_multi_newton_tri.cpp index 2000927c63..2bf033781e 100755 --- a/src/npair_half_size_multi_newton_tri.cpp +++ b/src/npair_half_size_multi_newton_tri.cpp @@ -12,6 +12,7 @@ ------------------------------------------------------------------------- */ #include "npair_half_size_multi_newton_tri.h" +#include "neighbor.h" #include "neigh_list.h" #include "atom.h" #include "atom_vec.h" @@ -27,18 +28,19 @@ NPairHalfSizeMultiNewtonTri::NPairHalfSizeMultiNewtonTri(LAMMPS *lmp) : NPair(lm /* ---------------------------------------------------------------------- size particles binned neighbor list construction with Newton's 3rd law for triclinic - multi stencil is igroup-jgroup dependent + multi stencil is icollection-jcollection dependent each owned atom i checks its own bin and other bins in triclinic stencil every pair stored exactly once by some processor ------------------------------------------------------------------------- */ void NPairHalfSizeMultiNewtonTri::build(NeighList *list) { - int i,j,k,n,itype,jtype,igroup,jgroup,ibin,jbin,ns,js; + int i,j,k,n,itype,jtype,icollection,jcollection,ibin,jbin,ns,js; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; double radi,radsum,cutdistsq; int *neighptr,*s; + int *collection = neighbor->collection; double **x = atom->x; double *radius = atom->radius; int *type = atom->type; @@ -62,7 +64,7 @@ void NPairHalfSizeMultiNewtonTri::build(NeighList *list) n = 0; neighptr = ipage->vget(); itype = type[i]; - igroup = map_type_multi[itype]; + icollection = collection[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -70,12 +72,12 @@ void NPairHalfSizeMultiNewtonTri::build(NeighList *list) ibin = atom2bin[i]; - // loop through stencils for all groups - for (jgroup = 0; jgroup < n_multi_groups; jgroup++) { + // loop through stencils for all collections + for (jcollection = 0; jcollection < ncollections; jcollection++) { - // if same group use own bin - if(igroup == jgroup) jbin = ibin; - else jbin = coord2bin(x[i], jgroup); + // if same collection use own bin + if(icollection == jcollection) jbin = ibin; + else jbin = coord2bin(x[i], jcollection); // loop over all atoms in bins in stencil // stencil is empty if i larger than j @@ -86,15 +88,15 @@ void NPairHalfSizeMultiNewtonTri::build(NeighList *list) // (equal zyx and j <= i) // latter excludes self-self interaction but allows superposed atoms - s = stencil_multi[igroup][jgroup]; - ns = nstencil_multi[igroup][jgroup]; + s = stencil_multi[icollection][jcollection]; + ns = nstencil_multi[icollection][jcollection]; for (k = 0; k < ns; k++) { - js = binhead_multi[jgroup][jbin + s[k]]; + js = binhead_multi[jcollection][jbin + s[k]]; for (j = js; j >= 0; j = bins[j]) { - // if same size (same group), use half stencil - if(cutmultisq[igroup][igroup] == cutmultisq[jgroup][jgroup]){ + // if same size (same collection), use half stencil + if(cutcollectionsq[icollection][icollection] == cutcollectionsq[jcollection][jcollection]){ if (x[j][2] < ztmp) continue; if (x[j][2] == ztmp) { if (x[j][1] < ytmp) continue; diff --git a/src/npair_half_size_multi_old_newtoff.cpp b/src/npair_half_size_multi_old_newtoff.cpp index 35095b4856..fd22412f00 100644 --- a/src/npair_half_size_multi_old_newtoff.cpp +++ b/src/npair_half_size_multi_old_newtoff.cpp @@ -35,7 +35,7 @@ NPairHalfSizeMultiOldNewtoff::NPairHalfSizeMultiOldNewtoff(LAMMPS *lmp) : NPair( void NPairHalfSizeMultiOldNewtoff::build(NeighList *list) { - int i,j,k,m,n,itype,jtype,ibin,ns; + int i,j,k,n,itype,jtype,ibin,ns; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; double radi,radsum,cutdistsq; int *neighptr,*s; diff --git a/src/npair_half_size_multi_old_newton.cpp b/src/npair_half_size_multi_old_newton.cpp index 0112c0593f..1f6918625b 100644 --- a/src/npair_half_size_multi_old_newton.cpp +++ b/src/npair_half_size_multi_old_newton.cpp @@ -34,7 +34,7 @@ NPairHalfSizeMultiOldNewton::NPairHalfSizeMultiOldNewton(LAMMPS *lmp) : NPair(lm void NPairHalfSizeMultiOldNewton::build(NeighList *list) { - int i,j,k,m,n,itype,jtype,ibin,ns; + int i,j,k,n,itype,jtype,ibin,ns; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; double radi,radsum,cutdistsq; int *neighptr,*s; diff --git a/src/npair_half_size_multi_old_newton_tri.cpp b/src/npair_half_size_multi_old_newton_tri.cpp index f29c323944..e7477fada4 100644 --- a/src/npair_half_size_multi_old_newton_tri.cpp +++ b/src/npair_half_size_multi_old_newton_tri.cpp @@ -33,7 +33,7 @@ NPairHalfSizeMultiOldNewtonTri::NPairHalfSizeMultiOldNewtonTri(LAMMPS *lmp) : NP void NPairHalfSizeMultiOldNewtonTri::build(NeighList *list) { - int i,j,k,m,n,itype,jtype,ibin,ns; + int i,j,k,n,itype,jtype,ibin,ns; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; double radi,radsum,cutdistsq; int *neighptr,*s; diff --git a/src/nstencil.cpp b/src/nstencil.cpp index 3d86b57b9e..c6ed76b6eb 100644 --- a/src/nstencil.cpp +++ b/src/nstencil.cpp @@ -52,14 +52,14 @@ using namespace LAMMPS_NS; cutoff is not cutneighmaxsq, but max cutoff for that atom type no versions that allow ghost on (any need for it?) for multi: - create one stencil for each igroup-jgroup pairing - the full/half stencil label refers to the same-group stencil - a half list with newton on has a half same-group stencil - a full list or half list with newton off has a full same-group stencil - cross group stencils are always full to allow small-to-large lookups + create one stencil for each icollection-jcollection pairing + the full/half stencil label refers to the same-collection stencil + a half list with newton on has a half same-collection stencil + a full list or half list with newton off has a full same-collection stencil + cross collection stencils are always full to allow small-to-large lookups for orthogonal boxes, a half stencil includes bins to the "upper right" of central bin for triclinic, a half stencil includes bins in the z (3D) or y (2D) plane of self and above - cutoff is not cutneighmaxsq, but max cutoff for that atom group + cutoff is not cutneighmaxsq, but max cutoff for that atom collection no versions that allow ghost on (any need for it?) ------------------------------------------------------------------------- */ @@ -81,7 +81,7 @@ NStencil::NStencil(LAMMPS *lmp) : Pointers(lmp) flag_half_multi = nullptr; flag_skip_multi = nullptr; - bin_group_multi = nullptr; + bin_collection_multi = nullptr; dimension = domain->dimension; } @@ -107,7 +107,7 @@ NStencil::~NStencil() if (stencil_multi) { - int n = n_multi_groups; + int n = ncollections; memory->destroy(nstencil_multi); for (int i = 0; i < n; i++) { for (int j = 0; j < n; j++) { @@ -120,7 +120,7 @@ NStencil::~NStencil() memory->destroy(maxstencil_multi); memory->destroy(flag_half_multi); memory->destroy(flag_skip_multi); - memory->destroy(bin_group_multi); + memory->destroy(bin_collection_multi); memory->destroy(stencil_sx_multi); memory->destroy(stencil_sy_multi); @@ -156,9 +156,9 @@ void NStencil::copy_neighbor_info() cuttypesq = neighbor->cuttypesq; cutneighsq = neighbor->cutneighsq; - n_multi_groups = neighbor->n_multi_groups; - map_type_multi = neighbor->map_type_multi; - cutmultisq = neighbor->cutmultisq; + ncollections = neighbor->ncollections; + collection = neighbor->collection; + cutcollectionsq = neighbor->cutcollectionsq; // overwrite Neighbor cutoff with custom value set by requestor // only works for style = BIN (checked by Neighbor class) @@ -269,9 +269,9 @@ void NStencil::create_setup() } } } else { - int i, j, bin_group, smax; + int i, j, bin_collection, smax; double stencil_range; - int n = n_multi_groups; + int n = ncollections; if(nb) copy_bin_info_multi(); @@ -280,8 +280,8 @@ void NStencil::create_setup() "neighstencil:flag_half_multi"); memory->create(flag_skip_multi, n, n, "neighstencil:flag_skip_multi"); - memory->create(bin_group_multi, n, n, - "neighstencil:bin_group_multi"); + memory->create(bin_collection_multi, n, n, + "neighstencil:bin_collection_multi"); memory->create(stencil_sx_multi, n, n, "neighstencil:stencil_sx_multi"); @@ -335,25 +335,25 @@ void NStencil::create_setup() // Skip creation of unused stencils if (flag_skip_multi[i][j]) continue; - // Copy bin info for this pair of atom types - bin_group = bin_group_multi[i][j]; + // Copy bin info for this pair of atom collections + bin_collection = bin_collection_multi[i][j]; - stencil_binsizex_multi[i][j] = binsizex_multi[bin_group]; - stencil_binsizey_multi[i][j] = binsizey_multi[bin_group]; - stencil_binsizez_multi[i][j] = binsizez_multi[bin_group]; + stencil_binsizex_multi[i][j] = binsizex_multi[bin_collection]; + stencil_binsizey_multi[i][j] = binsizey_multi[bin_collection]; + stencil_binsizez_multi[i][j] = binsizez_multi[bin_collection]; - stencil_mbinx_multi[i][j] = mbinx_multi[bin_group]; - stencil_mbiny_multi[i][j] = mbiny_multi[bin_group]; - stencil_mbinz_multi[i][j] = mbinz_multi[bin_group]; + stencil_mbinx_multi[i][j] = mbinx_multi[bin_collection]; + stencil_mbiny_multi[i][j] = mbiny_multi[bin_collection]; + stencil_mbinz_multi[i][j] = mbinz_multi[bin_collection]; - stencil_range = sqrt(cutmultisq[i][j]); + stencil_range = sqrt(cutcollectionsq[i][j]); - sx = static_cast (stencil_range*bininvx_multi[bin_group]); - if (sx*binsizex_multi[bin_group] < stencil_range) sx++; - sy = static_cast (stencil_range*bininvy_multi[bin_group]); - if (sy*binsizey_multi[bin_group] < stencil_range) sy++; - sz = static_cast (stencil_range*bininvz_multi[bin_group]); - if (sz*binsizez_multi[bin_group] < stencil_range) sz++; + sx = static_cast (stencil_range*bininvx_multi[bin_collection]); + if (sx*binsizex_multi[bin_collection] < stencil_range) sx++; + sy = static_cast (stencil_range*bininvy_multi[bin_collection]); + if (sy*binsizey_multi[bin_collection] < stencil_range) sy++; + sz = static_cast (stencil_range*bininvz_multi[bin_collection]); + if (sz*binsizez_multi[bin_collection] < stencil_range) sz++; stencil_sx_multi[i][j] = sx; stencil_sy_multi[i][j] = sy; @@ -397,24 +397,24 @@ double NStencil::bin_distance(int i, int j, int k) /* ---------------------------------------------------------------------- - compute closest distance for a given atom grouping + compute closest distance for a given atom collection ------------------------------------------------------------------------- */ -double NStencil::bin_distance_multi(int i, int j, int k, int ig) +double NStencil::bin_distance_multi(int i, int j, int k, int ic) { double delx,dely,delz; - if (i > 0) delx = (i-1)*binsizex_multi[ig]; + if (i > 0) delx = (i-1)*binsizex_multi[ic]; else if (i == 0) delx = 0.0; - else delx = (i+1)*binsizex_multi[ig]; + else delx = (i+1)*binsizex_multi[ic]; - if (j > 0) dely = (j-1)*binsizey_multi[ig]; + if (j > 0) dely = (j-1)*binsizey_multi[ic]; else if (j == 0) dely = 0.0; - else dely = (j+1)*binsizey_multi[ig]; + else dely = (j+1)*binsizey_multi[ic]; - if (k > 0) delz = (k-1)*binsizez_multi[ig]; + if (k > 0) delz = (k-1)*binsizez_multi[ic]; else if (k == 0) delz = 0.0; - else delz = (k+1)*binsizez_multi[ig]; + else delz = (k+1)*binsizez_multi[ic]; return (delx*delx + dely*dely + delz*delz); } @@ -431,7 +431,7 @@ double NStencil::memory_usage() bytes += atom->ntypes*maxstencil_multi_old * sizeof(int); bytes += atom->ntypes*maxstencil_multi_old * sizeof(double); } else if (neighstyle == Neighbor::MULTI) { - int n = n_multi_groups; + int n = ncollections; for (int i = 0; i < n; i++) { for (int j = 0; j < n; j++) { bytes += maxstencil_multi[i][j] * sizeof(int); diff --git a/src/nstencil.h b/src/nstencil.h index b38276b45f..37cedec045 100644 --- a/src/nstencil.h +++ b/src/nstencil.h @@ -43,8 +43,8 @@ class NStencil : protected Pointers { // Arrays to store options for multi itype-jtype stencils bool **flag_half_multi; // flag creation of a half stencil for igroup-jgroup - bool **flag_skip_multi; // skip creation of igroup-jgroup stencils (for newton on) - int **bin_group_multi; // what group to use for bin information + bool **flag_skip_multi; // skip creation of icollection-jcollection stencils (for newton on) + int **bin_collection_multi; // what collection to use for bin information NStencil(class LAMMPS *); virtual ~NStencil(); @@ -66,9 +66,9 @@ class NStencil : protected Pointers { double cutneighmaxsq; double *cuttypesq; double **cutneighsq; - double **cutmultisq; - int n_multi_groups; - int *map_type_multi; + double **cutcollectionsq; + int ncollections; + int *collection; // data from NBin class @@ -112,7 +112,7 @@ class NStencil : protected Pointers { // methods for multi NStencil - double bin_distance_multi(int, int, int, int); // distance between bin corners for different types + double bin_distance_multi(int, int, int, int); // distance between bin corners for different collections void copy_bin_info_multi(); // copy multi info from NBin class virtual void set_stencil_properties(){} // determine which stencils to build and how }; diff --git a/src/nstencil_full_multi_2d.cpp b/src/nstencil_full_multi_2d.cpp index 4c96f41a56..b27e3c2a9c 100644 --- a/src/nstencil_full_multi_2d.cpp +++ b/src/nstencil_full_multi_2d.cpp @@ -29,7 +29,7 @@ NStencilFullMulti2d::NStencilFullMulti2d(LAMMPS *lmp) : NStencil(lmp) {} void NStencilFullMulti2d::set_stencil_properties() { - int n = n_multi_groups; + int n = ncollections; int i, j; // Always look up neighbor using full stencil and neighbor's bin @@ -38,7 +38,7 @@ void NStencilFullMulti2d::set_stencil_properties() for (j = 0; j < n; j++) { flag_half_multi[i][j] = 0; flag_skip_multi[i][j] = 0; - bin_group_multi[i][j] = j; + bin_collection_multi[i][j] = j; } } } @@ -49,33 +49,33 @@ void NStencilFullMulti2d::set_stencil_properties() void NStencilFullMulti2d::create() { - int igroup, jgroup, bin_group, i, j, k, ns; - int n = n_multi_groups; + int icollection, jcollection, bin_collection, i, j, k, ns; + int n = ncollections; double cutsq; - for (igroup = 0; igroup < n; igroup++) { - for (jgroup = 0; jgroup < n; jgroup++) { - if (flag_skip_multi[igroup][jgroup]) continue; + for (icollection = 0; icollection < n; icollection++) { + for (jcollection = 0; jcollection < n; jcollection++) { + if (flag_skip_multi[icollection][jcollection]) continue; ns = 0; - sx = stencil_sx_multi[igroup][jgroup]; - sy = stencil_sy_multi[igroup][jgroup]; + sx = stencil_sx_multi[icollection][jcollection]; + sy = stencil_sy_multi[icollection][jcollection]; - mbinx = stencil_mbinx_multi[igroup][jgroup]; - mbiny = stencil_mbiny_multi[igroup][jgroup]; + mbinx = stencil_mbinx_multi[icollection][jcollection]; + mbiny = stencil_mbiny_multi[icollection][jcollection]; - bin_group = bin_group_multi[igroup][jgroup]; + bin_collection = bin_collection_multi[icollection][jcollection]; - cutsq = cutmultisq[igroup][jgroup]; + cutsq = cutcollectionsq[icollection][jcollection]; for (j = -sy; j <= sy; j++) for (i = -sx; i <= sx; i++) - if (bin_distance_multi(i,j,0,bin_group) < cutsq) - stencil_multi[igroup][jgroup][ns++] = j*mbinx + i; + if (bin_distance_multi(i,j,0,bin_collection) < cutsq) + stencil_multi[icollection][jcollection][ns++] = j*mbinx + i; - nstencil_multi[igroup][jgroup] = ns; + nstencil_multi[icollection][jcollection] = ns; } } } diff --git a/src/nstencil_full_multi_3d.cpp b/src/nstencil_full_multi_3d.cpp index bc76433723..3d97bc10ad 100644 --- a/src/nstencil_full_multi_3d.cpp +++ b/src/nstencil_full_multi_3d.cpp @@ -29,7 +29,7 @@ NStencilFullMulti3d::NStencilFullMulti3d(LAMMPS *lmp) : NStencil(lmp) {} void NStencilFullMulti3d::set_stencil_properties() { - int n = n_multi_groups; + int n = ncollections; int i, j; // Always look up neighbor using full stencil and neighbor's bin @@ -39,7 +39,7 @@ void NStencilFullMulti3d::set_stencil_properties() for (j = 0; j < n; j++) { flag_half_multi[i][j] = 0; flag_skip_multi[i][j] = 0; - bin_group_multi[i][j] = j; + bin_collection_multi[i][j] = j; } } } @@ -50,37 +50,37 @@ void NStencilFullMulti3d::set_stencil_properties() void NStencilFullMulti3d::create() { - int igroup, jgroup, bin_group, i, j, k, ns; - int n = n_multi_groups; + int icollection, jcollection, bin_collection, i, j, k, ns; + int n = ncollections; double cutsq; - for (igroup = 0; igroup < n; igroup++) { - for (jgroup = 0; jgroup < n; jgroup++) { - if (flag_skip_multi[igroup][jgroup]) continue; + for (icollection = 0; icollection < n; icollection++) { + for (jcollection = 0; jcollection < n; jcollection++) { + if (flag_skip_multi[icollection][jcollection]) continue; ns = 0; - sx = stencil_sx_multi[igroup][jgroup]; - sy = stencil_sy_multi[igroup][jgroup]; - sz = stencil_sz_multi[igroup][jgroup]; + sx = stencil_sx_multi[icollection][jcollection]; + sy = stencil_sy_multi[icollection][jcollection]; + sz = stencil_sz_multi[icollection][jcollection]; - mbinx = stencil_mbinx_multi[igroup][jgroup]; - mbiny = stencil_mbiny_multi[igroup][jgroup]; - mbinz = stencil_mbinz_multi[igroup][jgroup]; + mbinx = stencil_mbinx_multi[icollection][jcollection]; + mbiny = stencil_mbiny_multi[icollection][jcollection]; + mbinz = stencil_mbinz_multi[icollection][jcollection]; - bin_group = bin_group_multi[igroup][jgroup]; + bin_collection = bin_collection_multi[icollection][jcollection]; - cutsq = cutmultisq[igroup][jgroup]; + cutsq = cutcollectionsq[icollection][jcollection]; for (k = -sz; k <= sz; k++) for (j = -sy; j <= sy; j++) for (i = -sx; i <= sx; i++) - if (bin_distance_multi(i,j,k,bin_group) < cutsq) - stencil_multi[igroup][jgroup][ns++] = + if (bin_distance_multi(i,j,k,bin_collection) < cutsq) + stencil_multi[icollection][jcollection][ns++] = k*mbiny*mbinx + j*mbinx + i; - nstencil_multi[igroup][jgroup] = ns; + nstencil_multi[icollection][jcollection] = ns; } } } diff --git a/src/nstencil_half_multi_2d.cpp b/src/nstencil_half_multi_2d.cpp index 8be8a1c151..57352c7c6f 100644 --- a/src/nstencil_half_multi_2d.cpp +++ b/src/nstencil_half_multi_2d.cpp @@ -30,26 +30,26 @@ NStencilHalfMulti2d::NStencilHalfMulti2d(LAMMPS *lmp) : void NStencilHalfMulti2d::set_stencil_properties() { - int n = n_multi_groups; + int n = ncollections; int i, j; - // Cross groups: use full stencil, looking one way through hierarchy + // Cross collections: use full stencil, looking one way through hierarchy // smaller -> larger => use full stencil in larger bin // larger -> smaller => no nstencil required // If cut offs are same, use half stencil for (i = 0; i < n; i++) { for (j = 0; j < n; j++) { - if(cutmultisq[i][i] > cutmultisq[j][j]) continue; + if(cutcollectionsq[i][i] > cutcollectionsq[j][j]) continue; flag_skip_multi[i][j] = 0; - if(cutmultisq[i][i] == cutmultisq[j][j]){ + if(cutcollectionsq[i][i] == cutcollectionsq[j][j]){ flag_half_multi[i][j] = 1; - bin_group_multi[i][j] = i; + bin_collection_multi[i][j] = i; } else { flag_half_multi[i][j] = 0; - bin_group_multi[i][j] = j; + bin_collection_multi[i][j] = j; } } } @@ -61,42 +61,42 @@ void NStencilHalfMulti2d::set_stencil_properties() void NStencilHalfMulti2d::create() { - int igroup, jgroup, bin_group, i, j, ns; - int n = n_multi_groups; + int icollection, jcollection, bin_collection, i, j, ns; + int n = ncollections; double cutsq; - for (igroup = 0; igroup < n; igroup++) { - for (jgroup = 0; jgroup < n; jgroup++) { - if (flag_skip_multi[igroup][jgroup]) continue; + for (icollection = 0; icollection < n; icollection++) { + for (jcollection = 0; jcollection < n; jcollection++) { + if (flag_skip_multi[icollection][jcollection]) continue; ns = 0; - sx = stencil_sx_multi[igroup][jgroup]; - sy = stencil_sy_multi[igroup][jgroup]; + sx = stencil_sx_multi[icollection][jcollection]; + sy = stencil_sy_multi[icollection][jcollection]; - mbinx = stencil_mbinx_multi[igroup][jgroup]; - mbiny = stencil_mbiny_multi[igroup][jgroup]; + mbinx = stencil_mbinx_multi[icollection][jcollection]; + mbiny = stencil_mbiny_multi[icollection][jcollection]; - bin_group = bin_group_multi[igroup][jgroup]; + bin_collection = bin_collection_multi[icollection][jcollection]; - cutsq = cutmultisq[igroup][jgroup]; + cutsq = cutcollectionsq[icollection][jcollection]; - if (flag_half_multi[igroup][jgroup]) { + if (flag_half_multi[icollection][jcollection]) { for (j = 0; j <= sy; j++) for (i = -sx; i <= sx; i++) if (j > 0 || (j == 0 && i > 0)) { - if (bin_distance_multi(i,j,0,bin_group) < cutsq) - stencil_multi[igroup][jgroup][ns++] = j*mbinx + i; + if (bin_distance_multi(i,j,0,bin_collection) < cutsq) + stencil_multi[icollection][jcollection][ns++] = j*mbinx + i; } } else { for (j = -sy; j <= sy; j++) for (i = -sx; i <= sx; i++) - if (bin_distance_multi(i,j,0,bin_group) < cutsq) - stencil_multi[igroup][jgroup][ns++] = j*mbinx + i; + if (bin_distance_multi(i,j,0,bin_collection) < cutsq) + stencil_multi[icollection][jcollection][ns++] = j*mbinx + i; } - nstencil_multi[igroup][jgroup] = ns; + nstencil_multi[icollection][jcollection] = ns; } } } diff --git a/src/nstencil_half_multi_2d_tri.cpp b/src/nstencil_half_multi_2d_tri.cpp index 4a4e74c4a0..6588248109 100755 --- a/src/nstencil_half_multi_2d_tri.cpp +++ b/src/nstencil_half_multi_2d_tri.cpp @@ -30,26 +30,26 @@ NStencilHalfMulti2dTri::NStencilHalfMulti2dTri(LAMMPS *lmp) : void NStencilHalfMulti2dTri::set_stencil_properties() { - int n = n_multi_groups; + int n = ncollections; int i, j; - // Cross groups: use full stencil, looking one way through hierarchy + // Cross collections: use full stencil, looking one way through hierarchy // smaller -> larger => use full stencil in larger bin // larger -> smaller => no nstencil required // If cut offs are same, use half stencil for (i = 0; i < n; i++) { for (j = 0; j < n; j++) { - if(cutmultisq[i][i] > cutmultisq[j][j]) continue; + if(cutcollectionsq[i][i] > cutcollectionsq[j][j]) continue; flag_skip_multi[i][j] = 0; - if(cutmultisq[i][i] == cutmultisq[j][j]){ + if(cutcollectionsq[i][i] == cutcollectionsq[j][j]){ flag_half_multi[i][j] = 1; - bin_group_multi[i][j] = i; + bin_collection_multi[i][j] = i; } else { flag_half_multi[i][j] = 0; - bin_group_multi[i][j] = j; + bin_collection_multi[i][j] = j; } } } @@ -61,40 +61,40 @@ void NStencilHalfMulti2dTri::set_stencil_properties() void NStencilHalfMulti2dTri::create() { - int igroup, jgroup, bin_group, i, j, ns; - int n = n_multi_groups; + int icollection, jcollection, bin_collection, i, j, ns; + int n = ncollections; double cutsq; - for (igroup = 0; igroup < n; igroup++) { - for (jgroup = 0; jgroup < n; jgroup++) { - if (flag_skip_multi[igroup][jgroup]) continue; + for (icollection = 0; icollection < n; icollection++) { + for (jcollection = 0; jcollection < n; jcollection++) { + if (flag_skip_multi[icollection][jcollection]) continue; ns = 0; - sx = stencil_sx_multi[igroup][jgroup]; - sy = stencil_sy_multi[igroup][jgroup]; + sx = stencil_sx_multi[icollection][jcollection]; + sy = stencil_sy_multi[icollection][jcollection]; - mbinx = stencil_mbinx_multi[igroup][jgroup]; - mbiny = stencil_mbiny_multi[igroup][jgroup]; + mbinx = stencil_mbinx_multi[icollection][jcollection]; + mbiny = stencil_mbiny_multi[icollection][jcollection]; - bin_group = bin_group_multi[igroup][jgroup]; + bin_collection = bin_collection_multi[icollection][jcollection]; - cutsq = cutmultisq[igroup][jgroup]; + cutsq = cutcollectionsq[icollection][jcollection]; - if (flag_half_multi[igroup][jgroup]) { + if (flag_half_multi[icollection][jcollection]) { for (j = 0; j <= sy; j++) for (i = -sx; i <= sx; i++) - if (bin_distance_multi(i,j,0,bin_group) < cutsq) - stencil_multi[igroup][jgroup][ns++] = j*mbinx + i; + if (bin_distance_multi(i,j,0,bin_collection) < cutsq) + stencil_multi[icollection][jcollection][ns++] = j*mbinx + i; } else { for (j = -sy; j <= sy; j++) for (i = -sx; i <= sx; i++) - if (bin_distance_multi(i,j,0,bin_group) < cutsq) - stencil_multi[igroup][jgroup][ns++] = j*mbinx + i; + if (bin_distance_multi(i,j,0,bin_collection) < cutsq) + stencil_multi[icollection][jcollection][ns++] = j*mbinx + i; } - nstencil_multi[igroup][jgroup] = ns; + nstencil_multi[icollection][jcollection] = ns; } } } diff --git a/src/nstencil_half_multi_3d.cpp b/src/nstencil_half_multi_3d.cpp index ad4970603c..19cd051392 100644 --- a/src/nstencil_half_multi_3d.cpp +++ b/src/nstencil_half_multi_3d.cpp @@ -30,26 +30,26 @@ NStencilHalfMulti3d::NStencilHalfMulti3d(LAMMPS *lmp) : void NStencilHalfMulti3d::set_stencil_properties() { - int n = n_multi_groups; + int n = ncollections; int i, j; - // Cross groups: use full stencil, looking one way through hierarchy + // Cross collections: use full stencil, looking one way through hierarchy // smaller -> larger => use full stencil in larger bin // larger -> smaller => no nstencil required // If cut offs are same, use half stencil for (i = 0; i < n; i++) { for (j = 0; j < n; j++) { - if(cutmultisq[i][i] > cutmultisq[j][j]) continue; + if(cutcollectionsq[i][i] > cutcollectionsq[j][j]) continue; flag_skip_multi[i][j] = 0; - if(cutmultisq[i][i] == cutmultisq[j][j]){ + if(cutcollectionsq[i][i] == cutcollectionsq[j][j]){ flag_half_multi[i][j] = 1; - bin_group_multi[i][j] = i; + bin_collection_multi[i][j] = i; } else { flag_half_multi[i][j] = 0; - bin_group_multi[i][j] = j; + bin_collection_multi[i][j] = j; } } } @@ -61,48 +61,48 @@ void NStencilHalfMulti3d::set_stencil_properties() void NStencilHalfMulti3d::create() { - int igroup, jgroup, bin_group, i, j, k, ns; - int n = n_multi_groups; + int icollection, jcollection, bin_collection, i, j, k, ns; + int n = ncollections; double cutsq; - for (igroup = 0; igroup < n; igroup++) { - for (jgroup = 0; jgroup < n; jgroup++) { - if (flag_skip_multi[igroup][jgroup]) continue; + for (icollection = 0; icollection < n; icollection++) { + for (jcollection = 0; jcollection < n; jcollection++) { + if (flag_skip_multi[icollection][jcollection]) continue; ns = 0; - sx = stencil_sx_multi[igroup][jgroup]; - sy = stencil_sy_multi[igroup][jgroup]; - sz = stencil_sz_multi[igroup][jgroup]; + sx = stencil_sx_multi[icollection][jcollection]; + sy = stencil_sy_multi[icollection][jcollection]; + sz = stencil_sz_multi[icollection][jcollection]; - mbinx = stencil_mbinx_multi[igroup][jgroup]; - mbiny = stencil_mbiny_multi[igroup][jgroup]; - mbinz = stencil_mbinz_multi[igroup][jgroup]; + mbinx = stencil_mbinx_multi[icollection][jcollection]; + mbiny = stencil_mbiny_multi[icollection][jcollection]; + mbinz = stencil_mbinz_multi[icollection][jcollection]; - bin_group = bin_group_multi[igroup][jgroup]; + bin_collection = bin_collection_multi[icollection][jcollection]; - cutsq = cutmultisq[igroup][jgroup]; + cutsq = cutcollectionsq[icollection][jcollection]; - if (flag_half_multi[igroup][jgroup]) { + if (flag_half_multi[icollection][jcollection]) { for (k = 0; k <= sz; k++) for (j = -sy; j <= sy; j++) for (i = -sx; i <= sx; i++) if (k > 0 || j > 0 || (j == 0 && i > 0)) { - if (bin_distance_multi(i,j,k,bin_group) < cutsq) - stencil_multi[igroup][jgroup][ns++] = + if (bin_distance_multi(i,j,k,bin_collection) < cutsq) + stencil_multi[icollection][jcollection][ns++] = k*mbiny*mbinx + j*mbinx + i; } } else { for (k = -sz; k <= sz; k++) for (j = -sy; j <= sy; j++) for (i = -sx; i <= sx; i++) - if (bin_distance_multi(i,j,k,bin_group) < cutsq) - stencil_multi[igroup][jgroup][ns++] = + if (bin_distance_multi(i,j,k,bin_collection) < cutsq) + stencil_multi[icollection][jcollection][ns++] = k*mbiny*mbinx + j*mbinx + i; } - nstencil_multi[igroup][jgroup] = ns; + nstencil_multi[icollection][jcollection] = ns; } } } diff --git a/src/nstencil_half_multi_3d_tri.cpp b/src/nstencil_half_multi_3d_tri.cpp index 779f0865a5..3c5efac463 100755 --- a/src/nstencil_half_multi_3d_tri.cpp +++ b/src/nstencil_half_multi_3d_tri.cpp @@ -30,26 +30,26 @@ NStencilHalfMulti3dTri::NStencilHalfMulti3dTri(LAMMPS *lmp) : void NStencilHalfMulti3dTri::set_stencil_properties() { - int n = n_multi_groups; + int n = ncollections; int i, j; - // Cross groups: use full stencil, looking one way through hierarchy + // Cross collections: use full stencil, looking one way through hierarchy // smaller -> larger => use full stencil in larger bin // larger -> smaller => no nstencil required // If cut offs are same, use half stencil for (i = 0; i < n; i++) { for (j = 0; j < n; j++) { - if(cutmultisq[i][i] > cutmultisq[j][j]) continue; + if(cutcollectionsq[i][i] > cutcollectionsq[j][j]) continue; flag_skip_multi[i][j] = 0; - if(cutmultisq[i][i] == cutmultisq[j][j]){ + if(cutcollectionsq[i][i] == cutcollectionsq[j][j]){ flag_half_multi[i][j] = 1; - bin_group_multi[i][j] = i; + bin_collection_multi[i][j] = i; } else { flag_half_multi[i][j] = 0; - bin_group_multi[i][j] = j; + bin_collection_multi[i][j] = j; } } } @@ -61,46 +61,46 @@ void NStencilHalfMulti3dTri::set_stencil_properties() void NStencilHalfMulti3dTri::create() { - int igroup, jgroup, bin_group, i, j, k, ns; - int n = n_multi_groups; + int icollection, jcollection, bin_collection, i, j, k, ns; + int n = ncollections; double cutsq; - for (igroup = 0; igroup < n; igroup++) { - for (jgroup = 0; jgroup < n; jgroup++) { - if (flag_skip_multi[igroup][jgroup]) continue; + for (icollection = 0; icollection < n; icollection++) { + for (jcollection = 0; jcollection < n; jcollection++) { + if (flag_skip_multi[icollection][jcollection]) continue; ns = 0; - sx = stencil_sx_multi[igroup][jgroup]; - sy = stencil_sy_multi[igroup][jgroup]; - sz = stencil_sz_multi[igroup][jgroup]; + sx = stencil_sx_multi[icollection][jcollection]; + sy = stencil_sy_multi[icollection][jcollection]; + sz = stencil_sz_multi[icollection][jcollection]; - mbinx = stencil_mbinx_multi[igroup][jgroup]; - mbiny = stencil_mbiny_multi[igroup][jgroup]; - mbinz = stencil_mbinz_multi[igroup][jgroup]; + mbinx = stencil_mbinx_multi[icollection][jcollection]; + mbiny = stencil_mbiny_multi[icollection][jcollection]; + mbinz = stencil_mbinz_multi[icollection][jcollection]; - bin_group = bin_group_multi[igroup][jgroup]; + bin_collection = bin_collection_multi[icollection][jcollection]; - cutsq = cutmultisq[igroup][jgroup]; + cutsq = cutcollectionsq[icollection][jcollection]; - if (flag_half_multi[igroup][jgroup]) { + if (flag_half_multi[icollection][jcollection]) { for (k = 0; k <= sz; k++) for (j = -sy; j <= sy; j++) for (i = -sx; i <= sx; i++) - if (bin_distance_multi(i,j,k,bin_group) < cutsq) - stencil_multi[igroup][jgroup][ns++] = + if (bin_distance_multi(i,j,k,bin_collection) < cutsq) + stencil_multi[icollection][jcollection][ns++] = k*mbiny*mbinx + j*mbinx + i; } else { for (k = -sz; k <= sz; k++) for (j = -sy; j <= sy; j++) for (i = -sx; i <= sx; i++) - if (bin_distance_multi(i,j,k,bin_group) < cutsq) - stencil_multi[igroup][jgroup][ns++] = + if (bin_distance_multi(i,j,k,bin_collection) < cutsq) + stencil_multi[icollection][jcollection][ns++] = k*mbiny*mbinx + j*mbinx + i; } - nstencil_multi[igroup][jgroup] = ns; + nstencil_multi[icollection][jcollection] = ns; } } } diff --git a/src/pair.cpp b/src/pair.cpp index afc5eb2785..04aafbe863 100644 --- a/src/pair.cpp +++ b/src/pair.cpp @@ -63,6 +63,7 @@ Pair::Pair(LAMMPS *lmp) : Pointers(lmp) one_coeff = 0; no_virial_fdotr_compute = 0; writedata = 0; + finitecutflag = 0; ghostneigh = 0; unit_convert_flag = utils::NOCONVERT; diff --git a/src/pair.h b/src/pair.h index 5801941458..3c6710229e 100644 --- a/src/pair.h +++ b/src/pair.h @@ -56,6 +56,7 @@ class Pair : protected Pointers { int unit_convert_flag; // value != 0 indicates support for unit conversion. int no_virial_fdotr_compute; // 1 if does not invoke virial_fdotr_compute() int writedata; // 1 if writes coeffs to data file + int finitecutflag; // 1 if cut depends on finite atom size int ghostneigh; // 1 if pair style needs neighbors of ghosts double **cutghost; // cutoff for each ghost pair @@ -205,6 +206,8 @@ class Pair : protected Pointers { virtual void min_xf_get(int) {} virtual void min_x_set(int) {} virtual void transfer_history(double *, double*) {} + virtual double atom2cut(int) {return 0.0;} + virtual double radii2cut(double,double) {return 0.0;} // management of callbacks to be run from ev_tally() diff --git a/src/pair_hybrid.cpp b/src/pair_hybrid.cpp index 3026c9c741..e559e81f84 100644 --- a/src/pair_hybrid.cpp +++ b/src/pair_hybrid.cpp @@ -401,6 +401,7 @@ void PairHybrid::flags() if (styles[m]->dispersionflag) dispersionflag = 1; if (styles[m]->tip4pflag) tip4pflag = 1; if (styles[m]->compute_flag) compute_flag = 1; + if (styles[m]->finitecutflag) finitecutflag = 1; } single_enable = (single_enable == nstyles) ? 1 : 0; respa_enable = (respa_enable == nstyles) ? 1 : 0; @@ -1072,6 +1073,42 @@ int PairHybrid::check_ijtype(int itype, int jtype, char *substyle) return 0; } +/* ---------------------------------------------------------------------- + check if substyles calculate self-interaction range of particle +------------------------------------------------------------------------- */ + +double PairHybrid::atom2cut(int i) +{ + double temp, cut; + + cut = 0.0; + for (int m = 0; m < nstyles; m++) { + if (styles[m]->finitecutflag){ + temp = styles[m]->atom2cut(i); + if(temp > cut) cut = temp; + } + } + return cut; +} + +/* ---------------------------------------------------------------------- + check if substyles calculate maximum interaction range for two finite particles +------------------------------------------------------------------------- */ + +double PairHybrid::radii2cut(double r1, double r2) +{ + double temp, cut; + + cut = 0.0; + for (int m = 0; m < nstyles; m++) { + if (styles[m]->finitecutflag){ + temp = styles[m]->radii2cut(r1,r2); + if(temp > cut) cut = temp; + } + } + return cut; +} + /* ---------------------------------------------------------------------- memory usage of each sub-style ------------------------------------------------------------------------- */ diff --git a/src/pair_hybrid.h b/src/pair_hybrid.h index 7717d1fd51..9d10e91684 100644 --- a/src/pair_hybrid.h +++ b/src/pair_hybrid.h @@ -58,6 +58,8 @@ class PairHybrid : public Pair { virtual void add_tally_callback(class Compute *); virtual void del_tally_callback(class Compute *); + double atom2cut(int); + double radii2cut(double,double); protected: int nstyles; // # of sub-styles