Commit e42678ed authored by Tim Mattox's avatar Tim Mattox
Browse files

USER-DPD: track & use the extent of the local atoms in the bins

parent ff2786c8
Loading
Loading
Loading
Loading
+47 −26
Original line number Diff line number Diff line
@@ -30,22 +30,24 @@ using namespace LAMMPS_NS;
NBinSSA::NBinSSA(LAMMPS *lmp) : NBinStandard(lmp)
{
  maxbin_ssa = 0;
  bins_ssa = NULL;
  maxhead_ssa = 0;
  binhead_ssa = NULL;
  for (int i = 0; i < 9; i++) gairhead_ssa[i] = -1;
  binlist_ssa = NULL;
  binct_ssa = NULL;
  for (int i = 0; i < 9; i++) {
    gairhead_ssa[i] = -1;
    gairct_ssa[i] = 0;
  }
}

NBinSSA::~NBinSSA()
{
  memory->destroy(bins_ssa);
  memory->destroy(binhead_ssa);
  memory->destroy(binlist_ssa);
  memory->destroy(binct_ssa);
}

/* ----------------------------------------------------------------------
   bin owned and ghost atoms for the Shardlow Splitting Algorithm (SSA)
   local atoms are in distinct bins (binhead_ssa) from the ghosts
   ghost atoms are in distinct bins (gbinhead_ssa) from the locals
   local atoms are in distinct bins (binhead[]) from the ghosts
   ghost atoms are "binned" in gairhead_ssa[] instead
     ghosts which are not in an Active Interaction Region (AIR) are skipped
------------------------------------------------------------------------- */

@@ -58,6 +60,7 @@ void NBinSSA::bin_atoms()
  double **x = atom->x;
  int *mask = atom->mask;
  int *ssaAIR = atom->ssaAIR;
  int xbin,ybin,zbin;

  last_bin = update->ntimestep;

@@ -66,10 +69,13 @@ void NBinSSA::bin_atoms()

  for (i = 0; i < 9; i++) {
    gairhead_ssa[i] = -1;
    gairct_ssa[i] = 0;
  }

  for (i = 0; i < mbins; i++) {
    binhead_ssa[i] = -1;
    binhead[i] = -1;
    binlist_ssa[i] = -1;
    binct_ssa[i] = 0;
  }

  // bin in reverse order so linked list will be in forward order
@@ -81,23 +87,34 @@ void NBinSSA::bin_atoms()
      ibin = ssaAIR[i];
      if (ibin < 2) continue; // skip ghost atoms not in AIR
      if (mask[i] & bitmask) {
        bins_ssa[i] = gairhead_ssa[ibin];
        bins[i] = gairhead_ssa[ibin];
        gairhead_ssa[ibin] = i;
        ++(gairct_ssa[ibin]);
      }
    }
  } else {
    for (i = nall-1; i >= nlocal; i--) {
      ibin = ssaAIR[i];
      if (ibin < 2) continue; // skip ghost atoms not in AIR
      bins_ssa[i] = gairhead_ssa[ibin];
      bins[i] = gairhead_ssa[ibin];
      gairhead_ssa[ibin] = i;
      ++(gairct_ssa[ibin]);
    }
  }
  for (i = nlocal-1; i >= 0; i--) {
    ibin = coord2bin(x[i][0], x[i][1], x[i][2]);
    bins_ssa[i] = binhead_ssa[ibin];
    binhead_ssa[ibin] = i;
    ibin = coord2bin(x[i][0], x[i][1], x[i][2], xbin, ybin, zbin);
    // Find the bounding box of the local atoms in the bins
    if (xbin < lbinxlo) lbinxlo = xbin;
    if (xbin >= lbinxhi) lbinxhi = xbin + 1;
    if (ybin < lbinylo) lbinylo = ybin;
    if (ybin >= lbinyhi) lbinyhi = ybin + 1;
    if (zbin < lbinzlo) lbinzlo = zbin;
    if (zbin >= lbinzhi) lbinzhi = zbin + 1;
    bins[i] = binhead[ibin];
    binhead[ibin] = i;
    ++(binct_ssa[ibin]);
  }

}

/* ---------------------------------------------------------------------- */
@@ -106,17 +123,21 @@ void NBinSSA::bin_atoms_setup(int nall)
{
  NBinStandard::bin_atoms_setup(nall); // Setup the parent class's data too

  if (mbins > maxhead_ssa) {
    maxhead_ssa = mbins;
    memory->destroy(binhead_ssa);
    memory->create(binhead_ssa,maxhead_ssa,"binhead_ssa");
  if (mbins > maxbin_ssa) {
    maxbin_ssa = mbins;
    memory->destroy(binlist_ssa);
    memory->destroy(binct_ssa);
    memory->create(binlist_ssa,maxbin_ssa,"binlist_ssa");
    memory->create(binct_ssa,maxbin_ssa,"binct_ssa");
  }

  if (nall > maxbin_ssa) {
    maxbin_ssa = nall;
    memory->destroy(bins_ssa);
    memory->create(bins_ssa,maxbin_ssa,"bins_ssa");
  }
  // Clear the local bin extent bounding box.
  lbinxlo = mbinx - 1; // Safe to = stencil->sx + 1
  lbinylo = mbiny - 1; // Safe to = stencil->sy + 1
  lbinzlo = mbinz - 1; // Safe to = stencil->sz + 1
  lbinxhi = 0; // Safe to = mbinx - stencil->sx - 1
  lbinyhi = 0; // Safe to = mbiny - stencil->sy - 1
  lbinzhi = 0; // Safe to = mbinz - stencil->sz - 1
}

/* ---------------------------------------------------------------------- */
@@ -125,9 +146,9 @@ bigint NBinSSA::memory_usage()
{
  bigint bytes = NBinStandard::memory_usage(); // Count the parent's usage too

  if (maxbin_ssa) bytes += memory->usage(bins_ssa,maxbin_ssa);
  if (maxhead_ssa) {
    bytes += memory->usage(binhead_ssa,maxhead_ssa);
  if (maxbin_ssa) {
    bytes += memory->usage(binlist_ssa,maxbin_ssa);
    bytes += memory->usage(binct_ssa,maxbin_ssa);
  }
  return bytes;
}
+12 −4
Original line number Diff line number Diff line
@@ -29,11 +29,19 @@ namespace LAMMPS_NS {
class NBinSSA : public NBinStandard {
 public:

  int *bins_ssa;             // index of next atom in each bin
  int maxbin_ssa;            // size of bins_ssa array
  int *binhead_ssa;          // index of 1st local atom in each bin
  int *binlist_ssa;          // index in neighlist of 1st local atom in each bin
  int *binct_ssa;            // count of local atoms in each bin
  int gairhead_ssa[9];       // index of 1st ghost atom in each AIR
  int maxhead_ssa;           // size of binhead_ssa and gbinhead_ssa arrays
  int gairct_ssa[9];         // count of ghost atoms in each AIR
  int maxbin_ssa;            // size of binlist_ssa and binct_ssa arrays

  // Bounds of the local atoms in the binhead array
  int lbinxlo;               // lowest local bin x-dim coordinate
  int lbinylo;               // lowest local bin y-dim coordinate
  int lbinzlo;               // lowest local bin z-dim coordinate
  int lbinxhi;               // highest local bin x-dim coordinate
  int lbinyhi;               // highest local bin y-dim coordinate
  int lbinzhi;               // highest local bin z-dim coordinate

  NBinSSA(class LAMMPS *);
  ~NBinSSA();
+33 −16
Original line number Diff line number Diff line
@@ -79,20 +79,32 @@ void NPairHalfBinNewtonSSA::build(NeighList *list)

  NBinSSA *nb_ssa = dynamic_cast<NBinSSA*>(nb);
  if (!nb_ssa) error->one(FLERR, "NBin wasn't a NBinSSA object");
  int *bins_ssa = nb_ssa->bins_ssa;
  int *binhead_ssa = nb_ssa->binhead_ssa;
  int *bins = nb_ssa->bins;
  int *binhead = nb_ssa->binhead;
  int *binlist_ssa = nb_ssa->binlist_ssa;
  int *binct_ssa = nb_ssa->binct_ssa;
  int *gairhead_ssa = &(nb_ssa->gairhead_ssa[0]);

  int inum = 0;
  int gnum = 0;
  int xbin,ybin,zbin,xbin2,ybin2,zbin2;
  int **stencilxyz = ns_ssa->stencilxyz;
  int lbinxlo = nb_ssa->lbinxlo;
  int lbinxhi = nb_ssa->lbinxhi;
  int lbinylo = nb_ssa->lbinylo;
  int lbinyhi = nb_ssa->lbinyhi;
  int lbinzlo = nb_ssa->lbinzlo;
  int lbinzhi = nb_ssa->lbinzhi;

  ipage->reset();

  // loop over owned atoms, storing half of the neighbors

  for (i = 0; i < nlocal; i++) {
  // loop over bins with local atoms, storing half of the neighbors
  for (zbin = lbinzlo; zbin < lbinzhi; zbin++) {
  for (ybin = lbinylo; ybin < lbinyhi; ybin++) {
  for (xbin = lbinxlo; xbin < lbinxhi; xbin++) {
  ibin = zbin*mbiny*mbinx + ybin*mbinx + xbin;
  binlist_ssa[ibin] = inum; // record where ibin starts in ilist
  for (i = binhead[ibin]; i >= 0; i = bins[i]) {
    n = 0;
    neighptr = ipage->vget();

@@ -109,7 +121,7 @@ void NPairHalfBinNewtonSSA::build(NeighList *list)
    // loop over rest of local atoms in i's bin
    // just store them, since j is beyond i in linked list

    for (j = bins_ssa[i]; j >= 0; j = bins_ssa[j]) {
    for (j = bins[i]; j >= 0; j = bins[j]) {

      jtype = type[j];
      if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
@@ -136,13 +148,11 @@ void NPairHalfBinNewtonSSA::build(NeighList *list)
      }
    }

    ibin = coord2bin(x[i]);

    // loop over all local atoms in other bins in "half" stencil

    for (k = 0; k < nstencil_half; k++) {
      for (j = binhead_ssa[ibin+stencil[k]]; j >= 0;
           j = bins_ssa[j]) {
      for (j = binhead[ibin+stencil[k]]; j >= 0;
           j = bins[j]) {

        jtype = type[j];
        if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
@@ -177,13 +187,20 @@ void NPairHalfBinNewtonSSA::build(NeighList *list)
    if (ipage->status())
      error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
  }
  // verify count of atoms in ibin
  if (binct_ssa[ibin] != (inum - binlist_ssa[ibin]))
    error->one(FLERR,"binct_ssa didn't agree with lenght in ilist");
  }
  }
  }

  list->inum = inum;

  // loop over AIR ghost atoms, storing their local neighbors
  // since these are ghosts, must check if stencil bin is out of bounds
  for (int airnum = 2; airnum <= 8; airnum++) {
    for (i = gairhead_ssa[airnum]; i >= 0; i = bins_ssa[i]) {
    list->AIRct_ssa[airnum - 1] = nb_ssa->gairct_ssa[airnum];
    for (i = gairhead_ssa[airnum]; i >= 0; i = bins[i]) {
      n = 0;
      neighptr = ipage->vget();

@@ -205,11 +222,11 @@ void NPairHalfBinNewtonSSA::build(NeighList *list)
        xbin2 = xbin + stencilxyz[k][0];
        ybin2 = ybin + stencilxyz[k][1];
        zbin2 = zbin + stencilxyz[k][2];
        // since we only care about ghost to local neighbors, these "bounds" could be inset
        if (xbin2 < 0 || xbin2 >= mbinx ||
            ybin2 < 0 || ybin2 >= mbiny ||
            zbin2 < 0 || zbin2 >= mbinz) continue;
        for (j = binhead_ssa[ibin+stencil[k]]; j >= 0; j = bins_ssa[j]) {
        // Skip it if this bin is outside the extent of local bins
        if (xbin2 < lbinxlo || xbin2 >= lbinxhi ||
            ybin2 < lbinylo || ybin2 >= lbinyhi ||
            zbin2 < lbinzlo || zbin2 >= lbinzhi) continue;
        for (j = binhead[ibin+stencil[k]]; j >= 0; j = bins[j]) {

          jtype = type[j];
          if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
+1 −0
Original line number Diff line number Diff line
@@ -80,6 +80,7 @@ class NeighList : protected Pointers {

  // USER-DPD package and Shardlow Splitting Algorithm (SSA) support

  int AIRct_ssa[8]; // count of how many atoms in each AIR
  uint16_t (*ndxAIR_ssa)[8]; // for each atom, last neighbor index of each AIR

  // methods