diff --git a/doc/src/neighbor.rst b/doc/src/neighbor.rst index 0cfedfc090..663170ef47 100644 --- a/doc/src/neighbor.rst +++ b/doc/src/neighbor.rst @@ -49,29 +49,27 @@ sometimes be faster. Either style should give the same answers. The *multi* style is a modified binning algorithm that is useful for systems with a wide range of cutoff distances, e.g. due to different -size particles. For granular pair styles, cutoffs are set to the -sum of the maximum atomic radii for each atom type. -For the *bin* style, the bin size is set to 1/2 of -the largest cutoff distance between any pair of atom types and a -single set of bins is defined to search over for all atom types. This -can be inefficient if one pair of types has a very long cutoff, but -other type pairs have a much shorter cutoff. The *multi* style uses -different sized bins for collections of different sized particles, where -"size" may mean the physical size of the particle or its cutoff -distance for interacting with other particles. Different +size particles. For granular pair styles, cutoffs are set to the sum of +the maximum atomic radii for each atom type. For the *bin* style, the +bin size is set to 1/2 of the largest cutoff distance between any pair +of atom types and a single set of bins is defined to search over for all +atom types. This can be inefficient if one pair of types has a very +long cutoff, but other type pairs have a much shorter cutoff. The +*multi* style uses different sized bins for collections of different +sized particles, where "size" may mean the physical size of the particle +or its cutoff distance for interacting with other particles. Different sets of bins are then used to construct the neighbor lists as as further described by Shire, Hanley, and Stratford :ref:`(Shire) `. -This imposes some extra setup overhead, but the searches themselves -may be much faster. By default, each atom type defines a separate -collection of particles. For systems where two or more atom types -have the same size (either physical size or cutoff distance), the -definition of collections can be customized, which can result in less -overhead and faster performance. See the :doc:`neigh_modify ` -command for how to define custom collections. Whether the collection -definition is customized or not, also see the -:doc:`comm_modify mode multi ` command for communication -options that further improve performance in a manner consistent with -neighbor style multi. +This imposes some extra setup overhead, but the searches themselves may +be much faster. By default, each atom type defines a separate collection +of particles. For systems where two or more atom types have the same +size (either physical size or cutoff distance), the definition of +collections can be customized, which can result in less overhead and +faster performance. See the :doc:`neigh_modify ` command +for how to define custom collections. Whether the collection definition +is customized or not, also see the :doc:`comm_modify mode multi +` command for communication options that further improve +performance in a manner consistent with neighbor style multi. An alternate style, *multi/old*, sets the bin size to 1/2 of the shortest cutoff distance and multiple sets of bins are defined to search over for @@ -80,6 +78,16 @@ algorithm in LAMMPS but was found to be significantly slower than the new approach. For now we are keeping the old option in case there are use cases where multi/old outperforms the new multi style. +.. note:: + + If there are multiple sub-styles in a :doc:`hybrid/overlay pair style + ` that cover the same atom types, but have significantly + different cutoffs, the *multi* style does not apply. Instead, the + :doc:`pair_modify neigh/trim ` setting applies (which is + *yes* by default). Please check the neighbor list summary printed at + the beginning of a calculation to verify that the desired set of + neighbor list builds is performed. + The :doc:`neigh_modify ` command has additional options that control how often neighbor lists are built and which pairs are diff --git a/doc/src/pair_modify.rst b/doc/src/pair_modify.rst index 4941693fbd..6d4171fbc9 100644 --- a/doc/src/pair_modify.rst +++ b/doc/src/pair_modify.rst @@ -13,7 +13,7 @@ Syntax * one or more keyword/value pairs may be listed * keyword = *pair* or *shift* or *mix* or *table* or *table/disp* or *tabinner* or *tabinner/disp* or *tail* or *compute* or *nofdotr* or *special* or - *compute/tally* + *compute/tally* or *neigh/trim* .. parsed-literal:: @@ -37,6 +37,7 @@ Syntax which = *lj/coul* or *lj* or *coul* w1,w2,w3 = 1-2, 1-3, 1-4 weights from 0.0 to 1.0 inclusive *compute/tally* value = *yes* or *no* + *neigh/trim* value = *yes* or *no* Examples """""""" @@ -283,6 +284,31 @@ the *pair* keyword. Use *no* to disable, or *yes* to enable. The "pair_modify pair compute/tally" command must be issued **before** the corresponding compute style is defined. +The *neigh/trim* keyword controls whether an explicit cutoff is set for +each neighbor list request issued by individual pair sub-styles when +using :doc:`pair hybrid/overlay `. When this keyword is +set to *no*, then the cutoff of each pair sub-style neighbor list will +be set equal to the largest cutoff, even if a shorter cutoff is +specified for a particular sub-style. If possible the neighbor list +will be copied directly from another list. When this keyword is set to +*yes* then the cutoff of the neighbor list will be explicitly set to the +value requested by the pair sub-style, and if possible the list will be +created by trimming neighbors from another list with a longer cutoff, +otherwise a new neighbor list will be created with the specified cutoff. +The *yes* option can be faster when there are multiple pair styles with +different cutoffs since the number of pair-wise distance checks between +neighbors is reduced (but the time required to build the neighbor lists +is increased). The *no* option could be faster when two or more neighbor +lists have similar (but not exactly the same) cutoffs. + +.. note:: + + The "pair_modify neigh/trim" command *only* applies when there are + multiple pair sub-styles for the same atoms with different cutoffs, + i.e. when using pair style hybrid/overlay. If you have different + cutoffs for different pairs for atoms type, the :doc:`neighbor style + multi ` should be used to create optimized neighbor lists. + ---------- Restrictions @@ -298,13 +324,13 @@ Related commands :doc:`pair_style `, :doc:`pair_style hybrid `, :doc:`pair_coeff `, :doc:`thermo_style `, -:doc:`compute \*/tally ` +:doc:`compute \*/tally `, :doc:`neighbor multi ` Default """"""" The option defaults are mix = geometric, shift = no, table = 12, -tabinner = sqrt(2.0), tail = no, and compute = yes. +tabinner = sqrt(2.0), tail = no, compute = yes, and neigh/trim yes. Note that some pair styles perform mixing, but only a certain style of mixing. See the doc pages for individual pair styles for details. diff --git a/doc/utils/sphinx-config/false_positives.txt b/doc/utils/sphinx-config/false_positives.txt index d5cd578da0..ade1f8794f 100644 --- a/doc/utils/sphinx-config/false_positives.txt +++ b/doc/utils/sphinx-config/false_positives.txt @@ -3219,6 +3219,7 @@ Souza sp spacings Spearot +specieslist specular spellcheck Spellmeyer diff --git a/src/Depend.sh b/src/Depend.sh index 98f3e5de4f..0d49512c3e 100755 --- a/src/Depend.sh +++ b/src/Depend.sh @@ -121,6 +121,7 @@ fi if (test $1 = "MANYBODY") then depend ATC depend GPU + depend INTEL depend KOKKOS depend OPT depend QEQ diff --git a/src/INTEL/npair_halffull_newtoff_trim_intel.h b/src/INTEL/npair_halffull_newtoff_trim_intel.h new file mode 100644 index 0000000000..923a0c8c95 --- /dev/null +++ b/src/INTEL/npair_halffull_newtoff_trim_intel.h @@ -0,0 +1,44 @@ +// clang-format off +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +/* ---------------------------------------------------------------------- + Contributing author: Stan Moore (SNL) +------------------------------------------------------------------------- */ + +// Only used for hybrid to generate list for non-intel style. Use +// standard routines. + +#ifdef NPAIR_CLASS +// clang-format off +NPairStyle(halffull/newtoff/trim/intel, + NPairHalffullNewtoffTrim, + NP_HALF_FULL | NP_NEWTOFF | NP_NSQ | NP_BIN | NP_MULTI | NP_HALF | + NP_ORTHO | NP_TRI | NP_TRIM | NP_INTEL); + +NPairStyle(halffull/newtoff/skip/trim/intel, + NPairHalffullNewtoffTrim, + NP_HALF_FULL | NP_NEWTOFF | NP_NSQ | NP_BIN | NP_MULTI | NP_HALF | + NP_ORTHO | NP_TRI | NP_SKIP | NP_TRIM | NP_INTEL); + +NPairStyle(halffull/newtoff/ghost/trim/intel, + NPairHalffullNewtoffTrim, + NP_HALF_FULL | NP_NEWTOFF | NP_NSQ | NP_BIN | NP_MULTI | NP_HALF | + NP_ORTHO | NP_TRI | NP_GHOST | NP_TRIM | NP_INTEL); + +NPairStyle(halffull/newtoff/skip/ghost/trim/intel, + NPairHalffullNewtoffTrim, + NP_HALF_FULL | NP_NEWTOFF | NP_NSQ | NP_BIN | NP_MULTI | NP_HALF | + NP_ORTHO | NP_TRI | NP_SKIP | NP_GHOST | NP_TRIM | NP_INTEL); +// clang-format on +#endif diff --git a/src/INTEL/npair_halffull_newton_trim_intel.cpp b/src/INTEL/npair_halffull_newton_trim_intel.cpp new file mode 100644 index 0000000000..e0c01086da --- /dev/null +++ b/src/INTEL/npair_halffull_newton_trim_intel.cpp @@ -0,0 +1,258 @@ +// clang-format off +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +/* ---------------------------------------------------------------------- + Contributing author: Stan Moore (SNL) +------------------------------------------------------------------------- */ + +#include "npair_halffull_newton_trim_intel.h" + +#include "atom.h" +#include "comm.h" +#include "error.h" +#include "modify.h" +#include "my_page.h" +#include "neigh_list.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NPairHalffullNewtonTrimIntel::NPairHalffullNewtonTrimIntel(LAMMPS *lmp) : NPair(lmp) { + _fix = static_cast(modify->get_fix_by_id("package_intel")); + if (!_fix) error->all(FLERR, "The 'package intel' command is required for /intel styles"); +} + +/* ---------------------------------------------------------------------- + build half list from full list and trim to shorter cutoff + pair stored once if i,j are both owned and i < j + if j is ghost, only store if j coords are "above and to the right" of i + works if full list is a skip list +------------------------------------------------------------------------- */ + +template +void NPairHalffullNewtonTrimIntel::build_t(NeighList *list, + IntelBuffers *buffers) +{ + const int inum_full = list->listfull->inum; + const int nlocal = atom->nlocal; + const int e_nall = nlocal + atom->nghost; + const ATOM_T * _noalias const x = buffers->get_x(); + int * _noalias const ilist = list->ilist; + int * _noalias const numneigh = list->numneigh; + int ** _noalias const firstneigh = list->firstneigh; + const int * _noalias const ilist_full = list->listfull->ilist; + const int * _noalias const numneigh_full = list->listfull->numneigh; + const int ** _noalias const firstneigh_full = (const int ** const)list->listfull->firstneigh; // NOLINT + + const flt_t cutsq_custom = cutoff_custom * cutoff_custom; + + #if defined(_OPENMP) + #pragma omp parallel + #endif + { + int tid, ifrom, ito; + IP_PRE_omp_range_id(ifrom, ito, tid, inum_full, comm->nthreads); + + // each thread has its own page allocator + MyPage &ipage = list->ipage[tid]; + ipage.reset(); + + // loop over parent full list + for (int ii = ifrom; ii < ito; ii++) { + int n = 0; + int *neighptr = ipage.vget(); + + const int i = ilist_full[ii]; + const flt_t xtmp = x[i].x; + const flt_t ytmp = x[i].y; + const flt_t ztmp = x[i].z; + + // loop over full neighbor list + + const int * _noalias const jlist = firstneigh_full[i]; + const int jnum = numneigh_full[i]; + + #if defined(LMP_SIMD_COMPILER) + #pragma vector aligned + #pragma ivdep + #endif + for (int jj = 0; jj < jnum; jj++) { + const int joriginal = jlist[jj]; + const int j = joriginal & NEIGHMASK; + int addme = 1; + if (j < nlocal) { + if (i > j) addme = 0; + } else { + if (x[j].z < ztmp) addme = 0; + if (x[j].z == ztmp) { + if (x[j].y < ytmp) addme = 0; + if (x[j].y == ytmp && x[j].x < xtmp) addme = 0; + } + } + + // trim to shorter cutoff + + const flt_t delx = xtmp - x[j].x; + const flt_t dely = ytmp - x[j].y; + const flt_t delz = ztmp - x[j].z; + const flt_t rsq = delx * delx + dely * dely + delz * delz; + + if (rsq > cutsq_custom) addme = 0; + + if (addme) + neighptr[n++] = joriginal; + } + + ilist[ii] = i; + firstneigh[i] = neighptr; + numneigh[i] = n; + + int pad_end = n; + IP_PRE_neighbor_pad(pad_end, 0); + #if defined(LMP_SIMD_COMPILER) + #pragma vector aligned + #pragma loop_count min=1, max=INTEL_COMPILE_WIDTH-1, \ + avg=INTEL_COMPILE_WIDTH/2 + #endif + for ( ; n < pad_end; n++) + neighptr[n] = e_nall; + + ipage.vgot(n); + if (ipage.status()) + error->one(FLERR,"Neighbor list overflow, boost neigh_modify one"); + } + } + list->inum = inum_full; +} + +/* ---------------------------------------------------------------------- + build half list from full 3-body list and trim to shorter cutoff + half list is already stored as first part of 3-body list +------------------------------------------------------------------------- */ + +template +void NPairHalffullNewtonTrimIntel::build_t3(NeighList *list, int *numhalf, + IntelBuffers *buffers) +{ + const int inum_full = list->listfull->inum; + const int e_nall = atom->nlocal + atom->nghost; + const ATOM_T * _noalias const x = buffers->get_x(); + int * _noalias const ilist = list->ilist; + int * _noalias const numneigh = list->numneigh; + int ** _noalias const firstneigh = list->firstneigh; + const int * _noalias const ilist_full = list->listfull->ilist; + const int * _noalias const numneigh_full = numhalf; + const int ** _noalias const firstneigh_full = (const int ** const)list->listfull->firstneigh; // NOLINT + + const flt_t cutsq_custom = cutoff_custom * cutoff_custom; + + int packthreads = 1; + if (comm->nthreads > INTEL_HTHREADS) packthreads = comm->nthreads; + + #if defined(_OPENMP) + #pragma omp parallel if (packthreads > 1) + #endif + { + int tid, ifrom, ito; + IP_PRE_omp_range_id(ifrom, ito, tid, inum_full, packthreads); + + // each thread has its own page allocator + MyPage &ipage = list->ipage[tid]; + ipage.reset(); + + // loop over parent full list + for (int ii = ifrom; ii < ito; ii++) { + int n = 0; + int *neighptr = ipage.vget(); + + const int i = ilist_full[ii]; + const flt_t xtmp = x[i].x; + const flt_t ytmp = x[i].y; + const flt_t ztmp = x[i].z; + + // loop over full neighbor list + + const int * _noalias const jlist = firstneigh_full[i]; + const int jnum = numneigh_full[ii]; + + #if defined(LMP_SIMD_COMPILER) + #pragma vector aligned + #pragma ivdep + #endif + for (int jj = 0; jj < jnum; jj++) { + const int joriginal = jlist[jj]; + const int j = joriginal & NEIGHMASK; + int addme = 1; + + // trim to shorter cutoff + + const flt_t delx = xtmp - x[j].x; + const flt_t dely = ytmp - x[j].y; + const flt_t delz = ztmp - x[j].z; + const flt_t rsq = delx * delx + dely * dely + delz * delz; + + if (rsq > cutsq_custom) addme = 0; + + if (addme) + neighptr[n++] = joriginal; + } + + ilist[ii] = i; + firstneigh[i] = neighptr; + numneigh[i] = n; + + int pad_end = n; + IP_PRE_neighbor_pad(pad_end, 0); + #if defined(LMP_SIMD_COMPILER) + #pragma vector aligned + #pragma loop_count min=1, max=INTEL_COMPILE_WIDTH-1, \ + avg=INTEL_COMPILE_WIDTH/2 + #endif + for ( ; n < pad_end; n++) + neighptr[n] = e_nall; + + ipage.vgot(n); + if (ipage.status()) + error->one(FLERR,"Neighbor list overflow, boost neigh_modify one"); + } + } + list->inum = inum_full; +} + +/* ---------------------------------------------------------------------- */ + +void NPairHalffullNewtonTrimIntel::build(NeighList *list) +{ + if (_fix->three_body_neighbor() == 0) { + if (_fix->precision() == FixIntel::PREC_MODE_MIXED) + build_t(list, _fix->get_mixed_buffers()); + else if (_fix->precision() == FixIntel::PREC_MODE_DOUBLE) + build_t(list, _fix->get_double_buffers()); + else + build_t(list, _fix->get_single_buffers()); + } else { + int *nhalf, *cnum; + if (_fix->precision() == FixIntel::PREC_MODE_MIXED) { + _fix->get_mixed_buffers()->get_list_data3(list->listfull, nhalf, cnum); + build_t3(list, nhalf, _fix->get_mixed_buffers()); + } else if (_fix->precision() == FixIntel::PREC_MODE_DOUBLE) { + _fix->get_double_buffers()->get_list_data3(list->listfull, nhalf, cnum); + build_t3(list, nhalf, _fix->get_double_buffers()); + } else { + _fix->get_single_buffers()->get_list_data3(list->listfull, nhalf, cnum); + build_t3(list, nhalf, _fix->get_single_buffers()); + } + } +} diff --git a/src/INTEL/npair_halffull_newton_trim_intel.h b/src/INTEL/npair_halffull_newton_trim_intel.h new file mode 100644 index 0000000000..096cf3bd66 --- /dev/null +++ b/src/INTEL/npair_halffull_newton_trim_intel.h @@ -0,0 +1,61 @@ +// clang-format off +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +/* ---------------------------------------------------------------------- + Contributing author: Stan Moore (SNL) +------------------------------------------------------------------------- */ + +#ifdef NPAIR_CLASS +// clang-format off +NPairStyle(halffull/newton/trim/intel, + NPairHalffullNewtonTrimIntel, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI| NP_TRIM | NP_INTEL); + +NPairStyle(halffull/newton/skip/trim/intel, + NPairHalffullNewtonTrimIntel, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_SKIP | NP_TRIM | NP_INTEL); +// clang-format on +#else + +#ifndef LMP_NPAIR_HALFFULL_NEWTON_TRIM_INTEL_H +#define LMP_NPAIR_HALFFULL_NEWTON_TRIM_INTEL_H + +#include "fix_intel.h" +#include "npair.h" + +#if defined(_OPENMP) +#include +#endif + +namespace LAMMPS_NS { + +class NPairHalffullNewtonTrimIntel : public NPair { + public: + NPairHalffullNewtonTrimIntel(class LAMMPS *); + void build(class NeighList *) override; + + protected: + FixIntel *_fix; + + template void build_t(NeighList *, IntelBuffers *); + + template void build_t3(NeighList *, int *, IntelBuffers *); +}; + +} // namespace LAMMPS_NS + +#endif +#endif diff --git a/src/INTEL/npair_trim_intel.cpp b/src/INTEL/npair_trim_intel.cpp new file mode 100644 index 0000000000..61c20db547 --- /dev/null +++ b/src/INTEL/npair_trim_intel.cpp @@ -0,0 +1,138 @@ +// clang-format off +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +/* ---------------------------------------------------------------------- + Contributing author: Stan Moore (SNL) +------------------------------------------------------------------------- */ + +#include "npair_trim_intel.h" + +#include "atom.h" +#include "comm.h" +#include "error.h" +#include "modify.h" +#include "my_page.h" +#include "neigh_list.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NPairTrimIntel::NPairTrimIntel(LAMMPS *lmp) : NPair(lmp) { + _fix = static_cast(modify->get_fix_by_id("package_intel")); + if (!_fix) error->all(FLERR, "The 'package intel' command is required for /intel styles"); +} + +/* ---------------------------------------------------------------------- + trim from copy list to shorter cutoff +------------------------------------------------------------------------- */ + +template +void NPairTrimIntel::build_t(NeighList *list, + IntelBuffers *buffers) +{ + const int inum_copy = list->listcopy->inum; + const int nlocal = atom->nlocal; + const int e_nall = nlocal + atom->nghost; + const ATOM_T * _noalias const x = buffers->get_x(); + int * _noalias const ilist = list->ilist; + int * _noalias const numneigh = list->numneigh; + int ** _noalias const firstneigh = list->firstneigh; + const int * _noalias const ilist_copy = list->listcopy->ilist; + const int * _noalias const numneigh_copy = list->listcopy->numneigh; + const int ** _noalias const firstneigh_copy = (const int ** const)list->listcopy->firstneigh; // NOLINT + + const flt_t cutsq_custom = cutoff_custom * cutoff_custom; + + #if defined(_OPENMP) + #pragma omp parallel + #endif + { + int tid, ifrom, ito; + IP_PRE_omp_range_id(ifrom, ito, tid, inum_copy, comm->nthreads); + + // each thread has its own page allocator + MyPage &ipage = list->ipage[tid]; + ipage.reset(); + + // loop over parent copy list + for (int ii = ifrom; ii < ito; ii++) { + int n = 0; + int *neighptr = ipage.vget(); + + const int i = ilist_copy[ii]; + const flt_t xtmp = x[i].x; + const flt_t ytmp = x[i].y; + const flt_t ztmp = x[i].z; + + // loop over copy neighbor list + + const int * _noalias const jlist = firstneigh_copy[i]; + const int jnum = numneigh_copy[i]; + + #if defined(LMP_SIMD_COMPILER) + #pragma vector aligned + #pragma ivdep + #endif + for (int jj = 0; jj < jnum; jj++) { + const int joriginal = jlist[jj]; + const int j = joriginal & NEIGHMASK; + int addme = 1; + + // trim to shorter cutoff + + const flt_t delx = xtmp - x[j].x; + const flt_t dely = ytmp - x[j].y; + const flt_t delz = ztmp - x[j].z; + const flt_t rsq = delx * delx + dely * dely + delz * delz; + + if (rsq > cutsq_custom) addme = 0; + + if (addme) + neighptr[n++] = joriginal; + } + + ilist[ii] = i; + firstneigh[i] = neighptr; + numneigh[i] = n; + + int pad_end = n; + IP_PRE_neighbor_pad(pad_end, 0); + #if defined(LMP_SIMD_COMPILER) + #pragma vector aligned + #pragma loop_count min=1, max=INTEL_COMPILE_WIDTH-1, \ + avg=INTEL_COMPILE_WIDTH/2 + #endif + for ( ; n < pad_end; n++) + neighptr[n] = e_nall; + + ipage.vgot(n); + if (ipage.status()) + error->one(FLERR,"Neighbor list overflow, boost neigh_modify one"); + } + } + list->inum = inum_copy; +} + +/* ---------------------------------------------------------------------- */ + +void NPairTrimIntel::build(NeighList *list) +{ + if (_fix->precision() == FixIntel::PREC_MODE_MIXED) + build_t(list, _fix->get_mixed_buffers()); + else if (_fix->precision() == FixIntel::PREC_MODE_DOUBLE) + build_t(list, _fix->get_double_buffers()); + else + build_t(list, _fix->get_single_buffers()); +} diff --git a/src/INTEL/npair_trim_intel.h b/src/INTEL/npair_trim_intel.h new file mode 100644 index 0000000000..6a68c19b26 --- /dev/null +++ b/src/INTEL/npair_trim_intel.h @@ -0,0 +1,53 @@ +// clang-format off +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +/* ---------------------------------------------------------------------- + Contributing author: Stan Moore (SNL) +------------------------------------------------------------------------- */ + +#ifdef NPAIR_CLASS +// clang-format off +NPairStyle(trim/intel, + NPairTrimIntel, + NP_COPY | NP_TRIM | NP_INTEL); +// clang-format on +#else + +#ifndef LMP_NPAIR_TRIM_INTEL_H +#define LMP_NPAIR_TRIM_INTEL_H + +#include "fix_intel.h" +#include "npair.h" + +#if defined(_OPENMP) +#include +#endif + +namespace LAMMPS_NS { + +class NPairTrimIntel : public NPair { + public: + NPairTrimIntel(class LAMMPS *); + void build(class NeighList *) override; + + protected: + FixIntel *_fix; + + template void build_t(NeighList *, IntelBuffers *); +}; + +} // namespace LAMMPS_NS + +#endif +#endif diff --git a/src/KOKKOS/Install.sh b/src/KOKKOS/Install.sh index 314205ea3e..6ad96be466 100755 --- a/src/KOKKOS/Install.sh +++ b/src/KOKKOS/Install.sh @@ -204,6 +204,8 @@ action npair_halffull_kokkos.cpp action npair_halffull_kokkos.h action npair_skip_kokkos.cpp action npair_skip_kokkos.h +action npair_trim_kokkos.cpp +action npair_trim_kokkos.h action npair_kokkos.cpp action npair_kokkos.h action npair_ssa_kokkos.cpp npair_half_bin_newton_ssa.cpp diff --git a/src/KOKKOS/npair_halffull_kokkos.cpp b/src/KOKKOS/npair_halffull_kokkos.cpp index 53bd132cac..3005f0e463 100644 --- a/src/KOKKOS/npair_halffull_kokkos.cpp +++ b/src/KOKKOS/npair_halffull_kokkos.cpp @@ -26,8 +26,8 @@ using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -template -NPairHalffullKokkos::NPairHalffullKokkos(LAMMPS *lmp) : NPair(lmp) { +template +NPairHalffullKokkos::NPairHalffullKokkos(LAMMPS *lmp) : NPair(lmp) { atomKK = (AtomKokkos *) atom; execution_space = ExecutionSpaceFromDevice::space; } @@ -41,15 +41,17 @@ NPairHalffullKokkos::NPairHalffullKokkos(LAMMPS *lmp) : NPair if ghost, also store neighbors of ghost atoms & set inum,gnum correctly ------------------------------------------------------------------------- */ -template -void NPairHalffullKokkos::build(NeighList *list) +template +void NPairHalffullKokkos::build(NeighList *list) { - if (NEWTON) { + if (NEWTON || TRIM) { x = atomKK->k_x.view(); atomKK->sync(execution_space,X_MASK); } nlocal = atom->nlocal; + cutsq_custom = cutoff_custom*cutoff_custom; + NeighListKokkos* k_list_full = static_cast*>(list->listfull); d_ilist_full = k_list_full->d_ilist; d_numneigh_full = k_list_full->d_numneigh; @@ -76,14 +78,14 @@ void NPairHalffullKokkos::build(NeighList *list) k_list->k_ilist.template modify(); } -template +template KOKKOS_INLINE_FUNCTION -void NPairHalffullKokkos::operator()(TagNPairHalffullCompute, const int &ii) const { +void NPairHalffullKokkos::operator()(TagNPairHalffullCompute, const int &ii) const { int n = 0; const int i = d_ilist_full(ii); F_FLOAT xtmp,ytmp,ztmp; - if (NEWTON) { + if (NEWTON || TRIM) { xtmp = x(i,0); ytmp = x(i,1); ztmp = x(i,2); @@ -108,9 +110,29 @@ void NPairHalffullKokkos::operator()(TagNPairHalffullCompute, if (x(j,1) == ytmp && x(j,0) < xtmp) continue; } } + + if (TRIM) { + const double delx = xtmp - x(j,0); + const double dely = ytmp - x(j,1); + const double delz = ztmp - x(j,2); + const double rsq = delx*delx + dely*dely + delz*delz; + + if (rsq > cutsq_custom) continue; + } + + neighbors_i(n++) = joriginal; + } else if (j > i) { + + if (TRIM) { + const double delx = xtmp - x(j,0); + const double dely = ytmp - x(j,1); + const double delz = ztmp - x(j,2); + const double rsq = delx*delx + dely*dely + delz*delz; + + if (rsq > cutsq_custom) continue; + } + neighbors_i(n++) = joriginal; - } else { - if (j > i) neighbors_i(n++) = joriginal; } } @@ -119,10 +141,14 @@ void NPairHalffullKokkos::operator()(TagNPairHalffullCompute, } namespace LAMMPS_NS { -template class NPairHalffullKokkos; -template class NPairHalffullKokkos; +template class NPairHalffullKokkos; +template class NPairHalffullKokkos; +template class NPairHalffullKokkos; +template class NPairHalffullKokkos; #ifdef LMP_KOKKOS_GPU -template class NPairHalffullKokkos; -template class NPairHalffullKokkos; +template class NPairHalffullKokkos; +template class NPairHalffullKokkos; +template class NPairHalffullKokkos; +template class NPairHalffullKokkos; #endif } diff --git a/src/KOKKOS/npair_halffull_kokkos.h b/src/KOKKOS/npair_halffull_kokkos.h index 6d4e722098..ce2321f6fd 100644 --- a/src/KOKKOS/npair_halffull_kokkos.h +++ b/src/KOKKOS/npair_halffull_kokkos.h @@ -13,27 +13,30 @@ #ifdef NPAIR_CLASS // clang-format off + +// Trim off + // Newton -typedef NPairHalffullKokkos NPairKokkosHalffullNewtonDevice; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonDevice; NPairStyle(halffull/newton/kk/device, NPairKokkosHalffullNewtonDevice, NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_KOKKOS_DEVICE); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtonHost; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonHost; NPairStyle(halffull/newton/kk/host, NPairKokkosHalffullNewtonHost, NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_KOKKOS_HOST); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtonDevice; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonDevice; NPairStyle(halffull/newton/skip/kk/device, NPairKokkosHalffullNewtonDevice, NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_SKIP | NP_KOKKOS_DEVICE); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtonHost; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonHost; NPairStyle(halffull/newton/skip/kk/host, NPairKokkosHalffullNewtonHost, NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | @@ -41,25 +44,25 @@ NPairStyle(halffull/newton/skip/kk/host, // Newtoff -typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffDevice; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffDevice; NPairStyle(halffull/newtoff/kk/device, NPairKokkosHalffullNewtoffDevice, NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_KOKKOS_DEVICE); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffHost; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffHost; NPairStyle(halffull/newtoff/kk/host, NPairKokkosHalffullNewtoffHost, NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_KOKKOS_HOST); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffDevice; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffDevice; NPairStyle(halffull/newtoff/skip/kk/device, NPairKokkosHalffullNewtoffDevice, NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_SKIP | NP_KOKKOS_DEVICE); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffHost; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffHost; NPairStyle(halffull/newtoff/skip/kk/host, NPairKokkosHalffullNewtoffHost, NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | @@ -69,25 +72,25 @@ NPairStyle(halffull/newtoff/skip/kk/host, // Newton -typedef NPairHalffullKokkos NPairKokkosHalffullNewtonGhostDevice; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonGhostDevice; NPairStyle(halffull/newton/ghost/kk/device, NPairKokkosHalffullNewtonGhostDevice, NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_GHOST | NP_KOKKOS_DEVICE); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtonHost; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonHost; NPairStyle(halffull/newton/ghost/kk/host, NPairKokkosHalffullNewtonHost, NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_GHOST | NP_KOKKOS_HOST); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtonGhostDevice; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonGhostDevice; NPairStyle(halffull/newton/skip/ghost/kk/device, NPairKokkosHalffullNewtonGhostDevice, NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_GHOST | NP_SKIP | NP_KOKKOS_DEVICE); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtonHost; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonHost; NPairStyle(halffull/newton/skip/ghost/kk/host, NPairKokkosHalffullNewtonHost, NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | @@ -95,29 +98,138 @@ NPairStyle(halffull/newton/skip/ghost/kk/host, // Newtoff -typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffGhostDevice; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffGhostDevice; NPairStyle(halffull/newtoff/ghost/kk/device, NPairKokkosHalffullNewtoffGhostDevice, NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_GHOST | NP_KOKKOS_DEVICE); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffHost; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffHost; NPairStyle(halffull/newtoff/ghost/kk/host, NPairKokkosHalffullNewtoffHost, NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_GHOST | NP_KOKKOS_HOST); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffGhostDevice; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffGhostDevice; NPairStyle(halffull/newtoff/skip/ghost/kk/device, NPairKokkosHalffullNewtoffGhostDevice, NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_GHOST | NP_SKIP | NP_KOKKOS_DEVICE); -typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffHost; +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffHost; NPairStyle(halffull/newtoff/skip/ghost/kk/host, NPairKokkosHalffullNewtoffHost, NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_ORTHO | NP_TRI | NP_GHOST | NP_SKIP | NP_KOKKOS_HOST); + + +//************ Trim ************** + +// Newton + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonTrimDevice; +NPairStyle(halffull/newton/trim/kk/device, + NPairKokkosHalffullNewtonTrimDevice, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_TRIM | NP_KOKKOS_DEVICE); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonTrimHost; +NPairStyle(halffull/newton/trim/kk/host, + NPairKokkosHalffullNewtonTrimHost, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_TRIM | NP_KOKKOS_HOST); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonTrimDevice; +NPairStyle(halffull/newton/skip/trim/kk/device, + NPairKokkosHalffullNewtonTrimDevice, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_SKIP | NP_TRIM | NP_KOKKOS_DEVICE); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonTrimHost; +NPairStyle(halffull/newton/skip/trim/kk/host, + NPairKokkosHalffullNewtonTrimHost, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_SKIP | NP_TRIM | NP_KOKKOS_HOST); + +// Newtoff + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffTrimDevice; +NPairStyle(halffull/newtoff/trim/kk/device, + NPairKokkosHalffullNewtoffTrimDevice, + NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_TRIM | NP_KOKKOS_DEVICE); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffTrimHost; +NPairStyle(halffull/newtoff/trim/kk/host, + NPairKokkosHalffullNewtoffTrimHost, + NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_TRIM | NP_KOKKOS_HOST); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffTrimDevice; +NPairStyle(halffull/newtoff/skip/trim/kk/device, + NPairKokkosHalffullNewtoffTrimDevice, + NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_SKIP | NP_TRIM | NP_KOKKOS_DEVICE); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffTrimHost; +NPairStyle(halffull/newtoff/skip/trim/kk/host, + NPairKokkosHalffullNewtoffTrimHost, + NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_SKIP | NP_TRIM | NP_KOKKOS_HOST); + +//************ Ghost ************** + +// Newton + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonGhostTrimDevice; +NPairStyle(halffull/newton/ghost/trim/kk/device, + NPairKokkosHalffullNewtonGhostTrimDevice, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_GHOST | NP_TRIM | NP_KOKKOS_DEVICE); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonTrimHost; +NPairStyle(halffull/newton/ghost/trim/kk/host, + NPairKokkosHalffullNewtonTrimHost, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_GHOST | NP_TRIM | NP_KOKKOS_HOST); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonGhostTrimDevice; +NPairStyle(halffull/newton/skip/ghost/trim/kk/device, + NPairKokkosHalffullNewtonGhostTrimDevice, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_GHOST | NP_SKIP | NP_TRIM | NP_KOKKOS_DEVICE); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtonTrimHost; +NPairStyle(halffull/newton/skip/ghost/trim/kk/host, + NPairKokkosHalffullNewtonTrimHost, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_GHOST | NP_SKIP | NP_TRIM | NP_KOKKOS_HOST); + +// Newtoff + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffGhostTrimDevice; +NPairStyle(halffull/newtoff/ghost/trim/kk/device, + NPairKokkosHalffullNewtoffGhostTrimDevice, + NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_GHOST | NP_TRIM | NP_KOKKOS_DEVICE); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffTrimHost; +NPairStyle(halffull/newtoff/ghost/trim/kk/host, + NPairKokkosHalffullNewtoffTrimHost, + NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_GHOST | NP_TRIM | NP_KOKKOS_HOST); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffGhostTrimDevice; +NPairStyle(halffull/newtoff/skip/ghost/trim/kk/device, + NPairKokkosHalffullNewtoffGhostTrimDevice, + NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_GHOST | NP_SKIP | NP_TRIM | NP_KOKKOS_DEVICE); + +typedef NPairHalffullKokkos NPairKokkosHalffullNewtoffTrimHost; +NPairStyle(halffull/newtoff/skip/ghost/trim/kk/host, + NPairKokkosHalffullNewtoffTrimHost, + NP_HALF_FULL | NP_NEWTOFF | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | + NP_ORTHO | NP_TRI | NP_GHOST | NP_SKIP | NP_TRIM | NP_KOKKOS_HOST); // clang-format on #else @@ -132,7 +244,7 @@ namespace LAMMPS_NS { struct TagNPairHalffullCompute{}; -template +template class NPairHalffullKokkos : public NPair { public: typedef DeviceType device_type; @@ -146,6 +258,7 @@ class NPairHalffullKokkos : public NPair { private: int nlocal; + double cutsq_custom; typename AT::t_x_array_randomread x; diff --git a/src/KOKKOS/npair_trim_kokkos.cpp b/src/KOKKOS/npair_trim_kokkos.cpp new file mode 100644 index 0000000000..bd2761f40c --- /dev/null +++ b/src/KOKKOS/npair_trim_kokkos.cpp @@ -0,0 +1,196 @@ +// clang-format off +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Trimright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#include "npair_trim_kokkos.h" +#include "atom_kokkos.h" +#include "atom_masks.h" +#include "neigh_list_kokkos.h" +#include "my_page.h" +#include "error.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +template +NPairTrimKokkos::NPairTrimKokkos(LAMMPS *lmp) : NPair(lmp) {} + +/* ---------------------------------------------------------------------- + create list which is simply a copy of parent list +------------------------------------------------------------------------- */ + +template +void NPairTrimKokkos::build(NeighList *list) +{ + NeighList *listcopy = list->listcopy; + + cutsq_custom = cutoff_custom*cutoff_custom; + + if (list->kokkos) { + if (!listcopy->kokkos) + error->all(FLERR,"Cannot trim non-Kokkos neighbor list to Kokkos neighbor list"); + trim_to_kokkos(list); + } else { + if (!listcopy->kokkos) + error->all(FLERR,"Missing Kokkos neighbor list for trim"); + trim_to_cpu(list); + } +} + +/* ---------------------------------------------------------------------- */ + +template +void NPairTrimKokkos::trim_to_kokkos(NeighList *list) +{ + x = atomKK->k_x.view(); + atomKK->sync(execution_space,X_MASK); + + cutsq_custom = cutoff_custom*cutoff_custom; + + NeighListKokkos* k_list_copy = static_cast*>(list->listcopy); + d_ilist_copy = k_list_copy->d_ilist; + d_numneigh_copy = k_list_copy->d_numneigh; + d_neighbors_copy = k_list_copy->d_neighbors; + int inum_copy = list->listcopy->inum; + if (list->ghost) inum_copy += list->listcopy->gnum; + + NeighListKokkos* k_list = static_cast*>(list); + k_list->maxneighs = k_list_copy->maxneighs; // simple, but could be made more memory efficient + k_list->grow(atom->nmax); + d_ilist = k_list->d_ilist; + d_numneigh = k_list->d_numneigh; + d_neighbors = k_list->d_neighbors; + + // loop over parent list and trim + + copymode = 1; + Kokkos::parallel_for(Kokkos::RangePolicy(0,inum_copy),*this); + copymode = 0; + + list->inum = k_list_copy->inum; + list->gnum = k_list_copy->gnum; + + k_list->k_ilist.template modify(); +} + +template +KOKKOS_INLINE_FUNCTION +void NPairTrimKokkos::operator()(TagNPairTrim, const int &ii) const { + int n = 0; + + const int i = d_ilist_copy(ii); + const double xtmp = x(i,0); + const double ytmp = x(i,1); + const double ztmp = x(i,2); + + // loop over copy neighbor list + + const int jnum = d_numneigh_copy(i); + const AtomNeighbors neighbors_i = AtomNeighbors(&d_neighbors(i,0),d_numneigh(i), + &d_neighbors(i,1)-&d_neighbors(i,0)); + + for (int jj = 0; jj < jnum; jj++) { + const int joriginal = d_neighbors_copy(i,jj); + const int j = joriginal & NEIGHMASK; + + const double delx = xtmp - x(j,0); + const double dely = ytmp - x(j,1); + const double delz = ztmp - x(j,2); + const double rsq = delx*delx + dely*dely + delz*delz; + + if (rsq > cutsq_custom) continue; + + neighbors_i(n++) = joriginal; + } + + d_numneigh(i) = n; + d_ilist(ii) = i; +} + +/* ---------------------------------------------------------------------- */ + +template +void NPairTrimKokkos::trim_to_cpu(NeighList *list) +{ + NeighList *listcopy = list->listcopy; + NeighListKokkos* listcopy_kk = (NeighListKokkos*) listcopy; + + listcopy_kk->k_ilist.template sync(); + + double** x = atom->x; + + int inum = listcopy->inum; + int gnum = listcopy->gnum; + int inum_all = inum; + if (list->ghost) inum_all += gnum; + auto h_ilist = listcopy_kk->k_ilist.h_view; + auto h_numneigh = Kokkos::create_mirror_view_and_copy(LMPHostType(),listcopy_kk->d_numneigh); + auto h_neighbors = Kokkos::create_mirror_view_and_copy(LMPHostType(),listcopy_kk->d_neighbors); + + list->inum = inum; + list->gnum = gnum; + auto ilist = list->ilist; + auto numneigh = list->numneigh; + + // Kokkos neighbor data is stored differently than regular CPU, + // must loop over lists + + int *neighptr; + int **firstneigh = list->firstneigh; + MyPage *ipage = list->ipage; + ipage->reset(); + + for (int ii = 0; ii < inum_all; ii++) { + int n = 0; + neighptr = ipage->vget(); + + const int i = h_ilist[ii]; + ilist[ii] = i; + const double xtmp = x[i][0]; + const double ytmp = x[i][1]; + const double ztmp = x[i][2]; + + // loop over Kokkos neighbor list + + const int jnum = h_numneigh[i]; + + for (int jj = 0; jj < jnum; jj++) { + const int joriginal = h_neighbors(i,jj); + + const int j = joriginal & NEIGHMASK; + + const double delx = xtmp - x[j][0]; + const double dely = ytmp - x[j][1]; + const double delz = ztmp - x[j][2]; + const double rsq = delx*delx + dely*dely + delz*delz; + + if (rsq > cutsq_custom) continue; + + neighptr[n++] = joriginal; + } + + firstneigh[i] = neighptr; + numneigh[i] = n; + ipage->vgot(n); + if (ipage->status()) + error->one(FLERR,"Neighbor list overflow, boost neigh_modify one"); + } +} + +namespace LAMMPS_NS { +template class NPairTrimKokkos; +#ifdef LMP_KOKKOS_GPU +template class NPairTrimKokkos; +#endif +} diff --git a/src/KOKKOS/npair_trim_kokkos.h b/src/KOKKOS/npair_trim_kokkos.h new file mode 100644 index 0000000000..67b60f7337 --- /dev/null +++ b/src/KOKKOS/npair_trim_kokkos.h @@ -0,0 +1,70 @@ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Trimright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#ifdef NPAIR_CLASS +// clang-format off +NPairStyle(trim/kk/device, + NPairTrimKokkos, + NP_COPY | NP_TRIM | NP_KOKKOS_DEVICE); + +NPairStyle(trim/kk/host, + NPairTrimKokkos, + NP_COPY | NP_TRIM | NP_KOKKOS_HOST); +// clang-format on +#else + +// clang-format off +#ifndef LMP_NPAIR_TRIM_KOKKOS_H +#define LMP_NPAIR_TRIM_KOKKOS_H + +#include "npair.h" +#include "kokkos_type.h" + +namespace LAMMPS_NS { + +struct TagNPairTrim{}; + +template +class NPairTrimKokkos : public NPair { + public: + typedef DeviceType device_type; + typedef ArrayTypes AT; + + NPairTrimKokkos(class LAMMPS *); + void build(class NeighList *) override; + + KOKKOS_INLINE_FUNCTION + void operator()(TagNPairTrim, const int&) const; + + private: + double cutsq_custom; + + typename AT::t_x_array_randomread x; + + typename AT::t_neighbors_2d_const d_neighbors_copy; + typename AT::t_int_1d_const d_ilist_copy; + typename AT::t_int_1d_const d_numneigh_copy; + + typename AT::t_neighbors_2d d_neighbors; + typename AT::t_int_1d d_ilist; + typename AT::t_int_1d d_numneigh; + + void trim_to_kokkos(class NeighList *); + void trim_to_cpu(class NeighList *); +}; + +} + +#endif +#endif + diff --git a/src/OPENMP/npair_halffull_newtoff_trim_omp.cpp b/src/OPENMP/npair_halffull_newtoff_trim_omp.cpp new file mode 100644 index 0000000000..20d497d53e --- /dev/null +++ b/src/OPENMP/npair_halffull_newtoff_trim_omp.cpp @@ -0,0 +1,111 @@ +// clang-format off +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#include "npair_halffull_newtoff_trim_omp.h" + +#include "atom.h" +#include "error.h" +#include "my_page.h" +#include "neigh_list.h" +#include "npair_omp.h" + +#include "omp_compat.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NPairHalffullNewtoffTrimOmp::NPairHalffullNewtoffTrimOmp(LAMMPS *lmp) : NPair(lmp) {} + +/* ---------------------------------------------------------------------- + build half list from full list and trim to shorter cutoff + 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) + works if full list is a skip list +------------------------------------------------------------------------- */ + +void NPairHalffullNewtoffTrimOmp::build(NeighList *list) +{ + const int inum_full = list->listfull->inum; + + NPAIR_OMP_INIT; + +#if defined(_OPENMP) +#pragma omp parallel LMP_DEFAULT_NONE LMP_SHARED(list) +#endif + NPAIR_OMP_SETUP(inum_full); + + int i,j,ii,jj,n,jnum,joriginal; + int *neighptr,*jlist; + double xtmp,ytmp,ztmp; + double delx,dely,delz,rsq; + + double **x = atom->x; + + int *ilist = list->ilist; + int *numneigh = list->numneigh; + int **firstneigh = list->firstneigh; + int *ilist_full = list->listfull->ilist; + int *numneigh_full = list->listfull->numneigh; + int **firstneigh_full = list->listfull->firstneigh; + + // each thread has its own page allocator + MyPage &ipage = list->ipage[tid]; + ipage.reset(); + + double cutsq_custom = cutoff_custom * cutoff_custom; + + // loop over atoms in full list + + for (ii = ifrom; ii < ito; ii++) { + + n = 0; + neighptr = ipage.vget(); + + // loop over parent full list + + i = ilist_full[ii]; + xtmp = x[i][0]; + ytmp = x[i][1]; + ztmp = x[i][2]; + + jlist = firstneigh_full[i]; + jnum = numneigh_full[i]; + + for (jj = 0; jj < jnum; jj++) { + joriginal = jlist[jj]; + j = joriginal & NEIGHMASK; + + // trim to shorter cutoff + + delx = xtmp - x[j][0]; + dely = ytmp - x[j][1]; + delz = ztmp - x[j][2]; + rsq = delx * delx + dely * dely + delz * delz; + + if (rsq > cutsq_custom) continue; + + if (j > i) neighptr[n++] = joriginal; + } + + ilist[ii] = i; + firstneigh[i] = neighptr; + numneigh[i] = n; + ipage.vgot(n); + if (ipage.status()) + error->one(FLERR,"Neighbor list overflow, boost neigh_modify one"); + } + NPAIR_OMP_CLOSE; + list->inum = inum_full; +} diff --git a/src/OPENMP/npair_halffull_newtoff_trim_omp.h b/src/OPENMP/npair_halffull_newtoff_trim_omp.h new file mode 100644 index 0000000000..6d701b5cb6 --- /dev/null +++ b/src/OPENMP/npair_halffull_newtoff_trim_omp.h @@ -0,0 +1,44 @@ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#ifdef NPAIR_CLASS +// clang-format off +NPairStyle(halffull/newtoff/trim/omp, + NPairHalffullNewtoffTrimOmp, + NP_HALF_FULL | NP_NEWTOFF | NP_NSQ | NP_BIN | NP_MULTI | NP_MULTI_OLD | NP_HALF | + NP_ORTHO | NP_TRI | NP_TRIM | NP_OMP); + +NPairStyle(halffull/newtoff/skip/trim/omp, + NPairHalffullNewtoffTrimOmp, + NP_HALF_FULL | NP_NEWTOFF | NP_NSQ | NP_BIN | NP_MULTI | NP_MULTI_OLD | NP_HALF | + NP_ORTHO | NP_TRI | NP_SKIP | NP_TRIM | NP_OMP); +// clang-format on +#else + +#ifndef LMP_NPAIR_HALFFULL_NEWTOFF_TRIM_OMP_H +#define LMP_NPAIR_HALFFULL_NEWTOFF_TRIM_OMP_H + +#include "npair.h" + +namespace LAMMPS_NS { + +class NPairHalffullNewtoffTrimOmp : public NPair { + public: + NPairHalffullNewtoffTrimOmp(class LAMMPS *); + void build(class NeighList *) override; +}; + +} // namespace LAMMPS_NS + +#endif +#endif diff --git a/src/OPENMP/npair_halffull_newton_trim_omp.cpp b/src/OPENMP/npair_halffull_newton_trim_omp.cpp new file mode 100644 index 0000000000..cefca47a5e --- /dev/null +++ b/src/OPENMP/npair_halffull_newton_trim_omp.cpp @@ -0,0 +1,120 @@ +// clang-format off +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#include "npair_halffull_newton_trim_omp.h" + +#include "atom.h" +#include "error.h" +#include "my_page.h" +#include "neigh_list.h" +#include "npair_omp.h" + +#include "omp_compat.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NPairHalffullNewtonTrimOmp::NPairHalffullNewtonTrimOmp(LAMMPS *lmp) : NPair(lmp) {} + +/* ---------------------------------------------------------------------- + build half list from full list and trim to shorter cutoff + pair stored once if i,j are both owned and i < j + if j is ghost, only store if j coords are "above and to the right" of i + works if full list is a skip list +------------------------------------------------------------------------- */ + +void NPairHalffullNewtonTrimOmp::build(NeighList *list) +{ + const int inum_full = list->listfull->inum; + + NPAIR_OMP_INIT; +#if defined(_OPENMP) +#pragma omp parallel LMP_DEFAULT_NONE LMP_SHARED(list) +#endif + NPAIR_OMP_SETUP(inum_full); + + int i,j,ii,jj,n,jnum,joriginal; + int *neighptr,*jlist; + double xtmp,ytmp,ztmp; + double delx,dely,delz,rsq; + + double **x = atom->x; + int nlocal = atom->nlocal; + + int *ilist = list->ilist; + int *numneigh = list->numneigh; + int **firstneigh = list->firstneigh; + int *ilist_full = list->listfull->ilist; + int *numneigh_full = list->listfull->numneigh; + int **firstneigh_full = list->listfull->firstneigh; + + // each thread has its own page allocator + MyPage &ipage = list->ipage[tid]; + ipage.reset(); + + double cutsq_custom = cutoff_custom * cutoff_custom; + + // loop over parent full list + + for (ii = ifrom; ii < ito; ii++) { + + n = 0; + neighptr = ipage.vget(); + + i = ilist_full[ii]; + xtmp = x[i][0]; + ytmp = x[i][1]; + ztmp = x[i][2]; + + // loop over full neighbor list + + jlist = firstneigh_full[i]; + jnum = numneigh_full[i]; + + for (jj = 0; jj < jnum; jj++) { + joriginal = jlist[jj]; + j = joriginal & NEIGHMASK; + if (j < nlocal) { + if (i > j) continue; + } else { + if (x[j][2] < ztmp) continue; + if (x[j][2] == ztmp) { + if (x[j][1] < ytmp) continue; + if (x[j][1] == ytmp && x[j][0] < xtmp) continue; + } + } + + // trim to shorter cutoff + + delx = xtmp - x[j][0]; + dely = ytmp - x[j][1]; + delz = ztmp - x[j][2]; + rsq = delx * delx + dely * dely + delz * delz; + + if (rsq > cutsq_custom) continue; + + neighptr[n++] = joriginal; + } + + ilist[ii] = i; + firstneigh[i] = neighptr; + numneigh[i] = n; + ipage.vgot(n); + if (ipage.status()) + error->one(FLERR,"Neighbor list overflow, boost neigh_modify one"); + } + NPAIR_OMP_CLOSE; + list->inum = inum_full; +} diff --git a/src/OPENMP/npair_halffull_newton_trim_omp.h b/src/OPENMP/npair_halffull_newton_trim_omp.h new file mode 100644 index 0000000000..a774a68c2b --- /dev/null +++ b/src/OPENMP/npair_halffull_newton_trim_omp.h @@ -0,0 +1,44 @@ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#ifdef NPAIR_CLASS +// clang-format off +NPairStyle(halffull/newton/trim/omp, + NPairHalffullNewtonTrimOmp, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_MULTI_OLD | + NP_ORTHO | NP_TRI| NP_TRIM | NP_OMP); + +NPairStyle(halffull/newton/skip/trim/omp, + NPairHalffullNewtonTrimOmp, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_MULTI_OLD | + NP_ORTHO | NP_TRI | NP_SKIP | NP_TRIM | NP_OMP); +// clang-format on +#else + +#ifndef LMP_NPAIR_HALFFULL_NEWTON_TRIM_OMP_H +#define LMP_NPAIR_HALFFULL_NEWTON_TRIM_OMP_H + +#include "npair.h" + +namespace LAMMPS_NS { + +class NPairHalffullNewtonTrimOmp : public NPair { + public: + NPairHalffullNewtonTrimOmp(class LAMMPS *); + void build(class NeighList *) override; +}; + +} // namespace LAMMPS_NS + +#endif +#endif diff --git a/src/OPENMP/npair_trim_omp.cpp b/src/OPENMP/npair_trim_omp.cpp new file mode 100644 index 0000000000..4b3b835f8c --- /dev/null +++ b/src/OPENMP/npair_trim_omp.cpp @@ -0,0 +1,108 @@ +// clang-format off +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#include "npair_trim_omp.h" + +#include "atom.h" +#include "error.h" +#include "my_page.h" +#include "neigh_list.h" +#include "npair_omp.h" + +#include "omp_compat.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NPairTrimOmp::NPairTrimOmp(LAMMPS *lmp) : NPair(lmp) {} + +/* ---------------------------------------------------------------------- + trim from copy list to shorter cutoff +------------------------------------------------------------------------- */ + +void NPairTrimOmp::build(NeighList *list) +{ + const int inum_copy = list->listcopy->inum; + + NPAIR_OMP_INIT; + +#if defined(_OPENMP) +#pragma omp parallel LMP_DEFAULT_NONE LMP_SHARED(list) +#endif + NPAIR_OMP_SETUP(inum_copy); + + int i,j,ii,jj,n,jnum,joriginal; + int *neighptr,*jlist; + double xtmp,ytmp,ztmp; + double delx,dely,delz,rsq; + + double **x = atom->x; + + int *ilist = list->ilist; + int *numneigh = list->numneigh; + int **firstneigh = list->firstneigh; + int *ilist_copy = list->listcopy->ilist; + int *numneigh_copy = list->listcopy->numneigh; + int **firstneigh_copy = list->listcopy->firstneigh; + + // each thread has its own page allocator + MyPage &ipage = list->ipage[tid]; + ipage.reset(); + + double cutsq_custom = cutoff_custom * cutoff_custom; + + // loop over atoms in copy list + + for (ii = ifrom; ii < ito; ii++) { + + n = 0; + neighptr = ipage.vget(); + + // loop over parent copy list + + i = ilist_copy[ii]; + xtmp = x[i][0]; + ytmp = x[i][1]; + ztmp = x[i][2]; + + jlist = firstneigh_copy[i]; + jnum = numneigh_copy[i]; + + for (jj = 0; jj < jnum; jj++) { + joriginal = jlist[jj]; + j = joriginal & NEIGHMASK; + + // trim to shorter cutoff + + delx = xtmp - x[j][0]; + dely = ytmp - x[j][1]; + delz = ztmp - x[j][2]; + rsq = delx * delx + dely * dely + delz * delz; + + if (rsq > cutsq_custom) continue; + + neighptr[n++] = joriginal; + } + + ilist[ii] = i; + firstneigh[i] = neighptr; + numneigh[i] = n; + ipage.vgot(n); + if (ipage.status()) + error->one(FLERR,"Neighbor list overflow, boost neigh_modify one"); + } + NPAIR_OMP_CLOSE; + list->inum = inum_copy; +} diff --git a/src/OPENMP/npair_trim_omp.h b/src/OPENMP/npair_trim_omp.h new file mode 100644 index 0000000000..48cd9eb6a0 --- /dev/null +++ b/src/OPENMP/npair_trim_omp.h @@ -0,0 +1,39 @@ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#ifdef NPAIR_CLASS +// clang-format off +NPairStyle(trim/omp, + NPairTrimOmp, + NP_COPY | NP_TRIM | NP_OMP) + +// clang-format on +#else + +#ifndef LMP_NPAIR_TRIM_OMP_H +#define LMP_NPAIR_TRIM_OMP_H + +#include "npair.h" + +namespace LAMMPS_NS { + +class NPairTrimOmp : public NPair { + public: + NPairTrimOmp(class LAMMPS *); + void build(class NeighList *) override; +}; + +} // namespace LAMMPS_NS + +#endif +#endif diff --git a/src/neigh_list.cpp b/src/neigh_list.cpp index dc76b242e9..d8fe7a8c27 100644 --- a/src/neigh_list.cpp +++ b/src/neigh_list.cpp @@ -48,6 +48,7 @@ NeighList::NeighList(LAMMPS *lmp) : Pointers(lmp) respamiddle = 0; respainner = 0; copy = 0; + trim = 0; copymode = 0; // ptrs @@ -97,7 +98,8 @@ NeighList::NeighList(LAMMPS *lmp) : Pointers(lmp) NeighList::~NeighList() { if (copymode) return; - if (!copy) { + + if (!copy || trim || kk2cpu) { memory->destroy(ilist); memory->destroy(numneigh); memory->sfree(firstneigh); @@ -144,6 +146,7 @@ void NeighList::post_constructor(NeighRequest *nq) respamiddle = nq->respamiddle; respainner = nq->respainner; copy = nq->copy; + trim = nq->trim; id = nq->id; if (nq->copy) { @@ -286,6 +289,8 @@ void NeighList::print_attributes() printf(" %d = skip flag\n",rq->skip); printf(" %d = off2on\n",rq->off2on); printf(" %d = copy flag\n",rq->copy); + printf(" %d = trim flag\n",rq->trim); + printf(" %d = kk2cpu flag\n",kk2cpu); printf(" %d = half/full\n",rq->halffull); printf("\n"); } diff --git a/src/neigh_list.h b/src/neigh_list.h index abd9b8205e..9b1bc238e6 100644 --- a/src/neigh_list.h +++ b/src/neigh_list.h @@ -42,6 +42,7 @@ class NeighList : protected Pointers { int respamiddle; // 1 if there is also a rRespa middle list int respainner; // 1 if there is also a rRespa inner list int copy; // 1 if this list is copied from another list + int trim; // 1 if this list is trimmed from another list int kk2cpu; // 1 if this list is copied from Kokkos to CPU int copymode; // 1 if this is a Kokkos on-device copy int id; // copied from neighbor list request diff --git a/src/neigh_request.cpp b/src/neigh_request.cpp index f859ca0fd4..37da50eaea 100644 --- a/src/neigh_request.cpp +++ b/src/neigh_request.cpp @@ -86,6 +86,7 @@ NeighRequest::NeighRequest(LAMMPS *_lmp) : Pointers(_lmp) skiplist = -1; off2on = 0; copy = 0; + trim = 0; copylist = -1; halffull = 0; halffulllist = -1; diff --git a/src/neigh_request.h b/src/neigh_request.h index 71c5987dbe..13b5b11466 100644 --- a/src/neigh_request.h +++ b/src/neigh_request.h @@ -102,6 +102,7 @@ class NeighRequest : protected Pointers { int off2on; // 1 if this is newton on list, but skips from off list int copy; // 1 if this list copied from another list + int trim; // 1 if this list trimmed from another list int copylist; // index of list to copy from int halffull; // 1 if half list computed from another full list diff --git a/src/neighbor.cpp b/src/neighbor.cpp index 996f1f4594..03e22f46c5 100644 --- a/src/neighbor.cpp +++ b/src/neighbor.cpp @@ -162,6 +162,7 @@ pairclass(nullptr), pairnames(nullptr), pairmasks(nullptr) nrequest = maxrequest = 0; requests = nullptr; + j_sorted = nullptr; old_nrequest = 0; old_requests = nullptr; @@ -254,6 +255,8 @@ Neighbor::~Neighbor() if (old_requests[i]) delete old_requests[i]; memory->sfree(old_requests); + delete[] j_sorted; + delete[] binclass; delete[] binnames; delete[] binmasks; @@ -842,15 +845,22 @@ int Neighbor::init_pair() // (3) after (2), b/c it adjusts lists created by (2) // (4) after (2) and (3), // b/c (2) may create new full lists, (3) may change them - // (5) last, after all lists are finalized, so all possible copies found + // (5) last, after all lists are finalized, so all possible copies/trims found int nrequest_original = nrequest; morph_unique(); morph_skip(); morph_granular(); // this method can change flags set by requestor + + // sort requests by cutoff distance for trimming, used by + // morph_halffull and morph_copy_trim. Must come after + // morph_skip() which change the number of requests + + sort_requests(); + morph_halffull(); - morph_copy(); + morph_copy_trim(); // create new lists, one per request including added requests // wait to allocate initial pages until copy lists are detected @@ -1009,8 +1019,8 @@ int Neighbor::init_pair() // allocate initial pages for each list, except if copy flag set for (i = 0; i < nlist; i++) { - if (lists[i]->copy && !lists[i]->kk2cpu) - continue; + if (lists[i]->copy && !lists[i]->trim && !lists[i]->kk2cpu) + continue; lists[i]->setup_pages(pgsize,oneatom); } @@ -1022,7 +1032,7 @@ int Neighbor::init_pair() int maxatom = atom->nmax; for (i = 0; i < nlist; i++) { - if (neigh_pair[i] && (!lists[i]->copy || lists[i]->kk2cpu)) + if (neigh_pair[i] && (!lists[i]->copy || lists[i]->trim || lists[i]->kk2cpu)) lists[i]->grow(maxatom,maxatom); } @@ -1094,6 +1104,43 @@ int Neighbor::init_pair() return same; } +/* ---------------------------------------------------------------------- + sort NeighRequests by cutoff distance + to find smallest list for trimming +------------------------------------------------------------------------- */ + +void Neighbor::sort_requests() +{ + NeighRequest *jrq; + int i,j,jmin; + double jcut; + + delete[] j_sorted; + j_sorted = new int[nrequest]; + + for (i = 0; i < nrequest; i++) + j_sorted[i] = i; + + for (i = 0; i < nrequest; i++) { + double cutoff_min = cutneighmax; + jmin = i; + + for (j = i; j < nrequest-1; j++) { + jrq = requests[j_sorted[j]]; + if (jrq->cut) jcut = jrq->cutoff; + else jcut = cutneighmax; + + if (jcut <= cutoff_min) { + cutoff_min = jcut; + jmin = j; + } + } + int tmp = j_sorted[i]; + j_sorted[i] = j_sorted[jmin]; + j_sorted[jmin] = tmp; + } +} + /* ---------------------------------------------------------------------- scan NeighRequests to set additional flags: custom cutoff lists and accelerator lists @@ -1106,10 +1153,23 @@ void Neighbor::morph_unique() for (int i = 0; i < nrequest; i++) { irq = requests[i]; - // if cut flag set by requestor, set unique flag + // if cut flag set by requestor and cutoff is different than default, + // set unique flag, otherwise unset cut flag // this forces Pair,Stencil,Bin styles to be instantiated separately + // also add skin to cutoff of perpetual lists - if (irq->cut) irq->unique = 1; + if (irq->cut) { + if (!irq->occasional) + irq->cutoff += skin; + + if (irq->cutoff != cutneighmax) { + irq->unique = 1; + + } else { + irq->cut = 0; + irq->cutoff = 0.0; + } + } // avoid flagging a neighbor list as both INTEL and OPENMP @@ -1294,11 +1354,13 @@ void Neighbor::morph_granular() void Neighbor::morph_halffull() { - int i,j; + int i,j,jj,jmin; NeighRequest *irq,*jrq; + double icut,jcut; for (i = 0; i < nrequest; i++) { irq = requests[i]; + int trim_flag = irq->trim; // only processing half lists @@ -1311,8 +1373,10 @@ void Neighbor::morph_halffull() // check all other lists - for (j = 0; j < nrequest; j++) { - if (i == j) continue; + for (jj = 0; jj < nrequest; jj++) { + if (irq->cut) j = j_sorted[jj]; + else j = jj; + jrq = requests[j]; // can only derive from a perpetual full list @@ -1321,10 +1385,20 @@ void Neighbor::morph_halffull() if (jrq->occasional) continue; if (!jrq->full) continue; + // trim a list with longer cutoff + + if (irq->cut) icut = irq->cutoff; + else icut = cutneighmax; + + if (jrq->cut) jcut = jrq->cutoff; + else jcut = cutneighmax; + + if (icut > jcut) continue; + else if (icut != jcut) trim_flag = 1; + // these flags must be same, // else 2 lists do not store same pairs // or their data structures are different - // this includes custom cutoff set by requestor if (irq->ghost != jrq->ghost) continue; if (irq->size != jrq->size) continue; @@ -1335,8 +1409,6 @@ void Neighbor::morph_halffull() if (irq->kokkos_host != jrq->kokkos_host) continue; if (irq->kokkos_device != jrq->kokkos_device) continue; if (irq->ssa != jrq->ssa) continue; - if (irq->cut != jrq->cut) continue; - if (irq->cutoff != jrq->cutoff) continue; // skip flag must be same // if both are skip lists, skip info must match @@ -1351,25 +1423,28 @@ void Neighbor::morph_halffull() // if matching list exists, point to it - if (j < nrequest) { + if (jj < nrequest) { irq->halffull = 1; irq->halffulllist = j; + irq->trim = trim_flag; } } } /* ---------------------------------------------------------------------- - scan NeighRequests for possible copies - if 2 requests match, turn one into a copy of the other + scan NeighRequests for possible copies or trims + if 2 requests match, turn one into a copy or trim of the other ------------------------------------------------------------------------- */ -void Neighbor::morph_copy() +void Neighbor::morph_copy_trim() { - int i,j,inewton,jnewton; + int i,j,jj,jmin,inewton,jnewton; NeighRequest *irq,*jrq; + double icut,jcut; for (i = 0; i < nrequest; i++) { irq = requests[i]; + int trim_flag = irq->trim; // this list is already a copy list due to another morph method @@ -1377,7 +1452,10 @@ void Neighbor::morph_copy() // check all other lists - for (j = 0; j < nrequest; j++) { + for (jj = 0; jj < nrequest; jj++) { + if (irq->cut) j = j_sorted[jj]; + else j = jj; + if (i == j) continue; jrq = requests[j]; @@ -1385,13 +1463,24 @@ void Neighbor::morph_copy() if (jrq->copy && jrq->copylist == i) continue; + // trim a list with longer cutoff + + if (irq->cut) icut = irq->cutoff; + else icut = cutneighmax; + + if (jrq->cut) jcut = jrq->cutoff; + else jcut = cutneighmax; + + if (icut > jcut) continue; + else if (icut != jcut) trim_flag = 1; + // other list (jrq) to copy from must be perpetual // list that becomes a copy list (irq) can be perpetual or occasional // if both lists are perpetual, require j < i // to prevent circular dependence with 3 or more copies of a list if (jrq->occasional) continue; - if (!irq->occasional && j > i) continue; + if (!irq->occasional && !irq->cut && j > i) continue; // both lists must be half, or both full @@ -1420,7 +1509,6 @@ void Neighbor::morph_copy() // these flags must be same, // else 2 lists do not store same pairs // or their data structures are different - // this includes custom cutoff set by requestor // no need to check omp b/c it stores same pairs // NOTE: need check for 2 Kokkos flags? @@ -1431,8 +1519,6 @@ void Neighbor::morph_copy() if (irq->kokkos_host && !jrq->kokkos_host) continue; if (irq->kokkos_device && !jrq->kokkos_device) continue; if (irq->ssa != jrq->ssa) continue; - if (irq->cut != jrq->cut) continue; - if (irq->cutoff != jrq->cutoff) continue; // skip flag must be same // if both are skip lists, skip info must match @@ -1448,10 +1534,13 @@ void Neighbor::morph_copy() // turn list I into a copy of list J // do not copy a list from another copy list, but from its parent list - if (j < nrequest) { + if (jj < nrequest) { irq->copy = 1; - if (jrq->copy) irq->copylist = jrq->copylist; - else irq->copylist = j; + irq->trim = trim_flag; + if (jrq->copy && irq->cutoff == requests[jrq->copylist]->cutoff) + irq->copylist = jrq->copylist; + else + irq->copylist = j; } } } @@ -1666,10 +1755,16 @@ void Neighbor::print_pairwise_info() // order these to get single output of most relevant - if (rq->copy) - out += fmt::format(", copy from ({})",rq->copylist+1); - else if (rq->halffull) - out += fmt::format(", half/full from ({})",rq->halffulllist+1); + if (rq->copy) { + if (rq->trim) + out += fmt::format(", trim from ({})",rq->copylist+1); + else + out += fmt::format(", copy from ({})",rq->copylist+1); + } else if (rq->halffull) + if (rq->trim) + out += fmt::format(", half/full trim from ({})",rq->halffulllist+1); + else + out += fmt::format(", half/full from ({})",rq->halffulllist+1); else if (rq->skip) out += fmt::format(", skip from ({})",rq->skiplist+1); out += "\n"; @@ -1970,10 +2065,16 @@ int Neighbor::choose_pair(NeighRequest *rq) // pairnames[i],pairmasks[i]); // if copy request, no further checks needed, just return or continue - // Kokkos device/host flags must also match in order to copy + // trim and Kokkos device/host flags must also match in order to copy + // intel and omp flags must match to trim if (rq->copy) { if (!(mask & NP_COPY)) continue; + if (rq->trim) { + if (!rq->trim != !(mask & NP_TRIM)) continue; + if (!rq->omp != !(mask & NP_OMP)) continue; + if (!rq->intel != !(mask & NP_INTEL)) continue; + } if (rq->kokkos_device || rq->kokkos_host) { if (!rq->kokkos_device != !(mask & NP_KOKKOS_DEVICE)) continue; if (!rq->kokkos_host != !(mask & NP_KOKKOS_HOST)) continue; @@ -2025,6 +2126,8 @@ int Neighbor::choose_pair(NeighRequest *rq) if (!rq->skip != !(mask & NP_SKIP)) continue; + if (!rq->trim != !(mask & NP_TRIM)) continue; + if (!rq->halffull != !(mask & NP_HALF_FULL)) continue; if (!rq->off2on != !(mask & NP_OFF2ON)) continue; @@ -2329,7 +2432,7 @@ void Neighbor::build(int topoflag) for (i = 0; i < npair_perpetual; i++) { m = plist[i]; - if (!lists[m]->copy || lists[m]->kk2cpu) + if (!lists[m]->copy || lists[m]->trim || lists[m]->kk2cpu) lists[m]->grow(nlocal,nall); neigh_pair[m]->build_setup(); neigh_pair[m]->build(lists[m]); @@ -2420,7 +2523,7 @@ void Neighbor::build_one(class NeighList *mylist, int preflag) // build the list - if (!mylist->copy || mylist->kk2cpu) + if (!mylist->copy || mylist->trim || mylist->kk2cpu) mylist->grow(atom->nlocal,atom->nlocal+atom->nghost); np->build_setup(); np->build(mylist); diff --git a/src/neighbor.h b/src/neighbor.h index 3492693766..7d1711c149 100644 --- a/src/neighbor.h +++ b/src/neighbor.h @@ -94,6 +94,7 @@ class Neighbor : protected Pointers { NeighList **lists; NeighRequest **requests; // from Pair,Fix,Compute,Command classes NeighRequest **old_requests; // copy of requests to compare to + int* j_sorted; // index of requests sorted by cutoff distance // data from topology neighbor lists @@ -245,11 +246,13 @@ class Neighbor : protected Pointers { int init_pair(); virtual void init_topology(); + void sort_requests(); + void morph_unique(); void morph_skip(); void morph_granular(); void morph_halffull(); - void morph_copy(); + void morph_copy_trim(); void print_pairwise_info(); void requests_new2old(); @@ -323,7 +326,8 @@ namespace NeighConst { NP_SKIP = 1 << 22, NP_HALF_FULL = 1 << 23, NP_OFF2ON = 1 << 24, - NP_MULTI_OLD = 1 << 25 + NP_MULTI_OLD = 1 << 25, + NP_TRIM = 1 << 26 }; enum { diff --git a/src/npair_halffull_newtoff_trim.cpp b/src/npair_halffull_newtoff_trim.cpp new file mode 100644 index 0000000000..28b285ddf3 --- /dev/null +++ b/src/npair_halffull_newtoff_trim.cpp @@ -0,0 +1,101 @@ +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#include "npair_halffull_newtoff_trim.h" + +#include "atom.h" +#include "error.h" +#include "my_page.h" +#include "neigh_list.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NPairHalffullNewtoffTrim::NPairHalffullNewtoffTrim(LAMMPS *lmp) : NPair(lmp) {} + +/* ---------------------------------------------------------------------- + build half list from full list + 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) + works if full list is a skip list + works for owned (non-ghost) list, also for ghost list + if ghost, also store neighbors of ghost atoms & set inum,gnum correctly +------------------------------------------------------------------------- */ + +void NPairHalffullNewtoffTrim::build(NeighList *list) +{ + int i, j, ii, jj, n, jnum, joriginal; + int *neighptr, *jlist; + double xtmp,ytmp,ztmp; + double delx,dely,delz,rsq; + + double **x = atom->x; + + int *ilist = list->ilist; + int *numneigh = list->numneigh; + int **firstneigh = list->firstneigh; + MyPage *ipage = list->ipage; + + int *ilist_full = list->listfull->ilist; + int *numneigh_full = list->listfull->numneigh; + int **firstneigh_full = list->listfull->firstneigh; + int inum_full = list->listfull->inum; + if (list->ghost) inum_full += list->listfull->gnum; + + int inum = 0; + ipage->reset(); + + double cutsq_custom = cutoff_custom * cutoff_custom; + + // loop over atoms in full list + + for (ii = 0; ii < inum_full; ii++) { + n = 0; + neighptr = ipage->vget(); + + // loop over parent full list + + i = ilist_full[ii]; + xtmp = x[i][0]; + ytmp = x[i][1]; + ztmp = x[i][2]; + + jlist = firstneigh_full[i]; + jnum = numneigh_full[i]; + + for (jj = 0; jj < jnum; jj++) { + joriginal = jlist[jj]; + j = joriginal & NEIGHMASK; + if (j > i) { + delx = xtmp - x[j][0]; + dely = ytmp - x[j][1]; + delz = ztmp - x[j][2]; + rsq = delx * delx + dely * dely + delz * delz; + + if (rsq > cutsq_custom) continue; + + neighptr[n++] = joriginal; + } + } + + ilist[inum++] = i; + firstneigh[i] = neighptr; + numneigh[i] = n; + ipage->vgot(n); + if (ipage->status()) error->one(FLERR, "Neighbor list overflow, boost neigh_modify one"); + } + + list->inum = inum; + if (list->ghost) list->gnum = list->listfull->gnum; +} diff --git a/src/npair_halffull_newtoff_trim.h b/src/npair_halffull_newtoff_trim.h new file mode 100644 index 0000000000..a5e39e4800 --- /dev/null +++ b/src/npair_halffull_newtoff_trim.h @@ -0,0 +1,54 @@ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#ifdef NPAIR_CLASS +// clang-format off +NPairStyle(halffull/newtoff/trim, + NPairHalffullNewtoffTrim, + NP_HALF_FULL | NP_NEWTOFF | NP_NSQ | NP_BIN | NP_MULTI | NP_MULTI_OLD | NP_HALF | + NP_ORTHO | NP_TRI | NP_TRIM); + +NPairStyle(halffull/newtoff/skip/trim, + NPairHalffullNewtoffTrim, + NP_HALF_FULL | NP_NEWTOFF | NP_NSQ | NP_BIN | NP_MULTI | NP_MULTI_OLD | NP_HALF | + NP_ORTHO | NP_TRI | NP_SKIP | NP_TRIM); + +NPairStyle(halffull/newtoff/ghost/trim, + NPairHalffullNewtoffTrim, + NP_HALF_FULL | NP_NEWTOFF | NP_NSQ | NP_BIN | NP_MULTI | NP_MULTI_OLD | NP_HALF | + NP_ORTHO | NP_TRI | NP_GHOST | NP_TRIM); + +NPairStyle(halffull/newtoff/skip/ghost/trim, + NPairHalffullNewtoffTrim, + NP_HALF_FULL | NP_NEWTOFF | NP_NSQ | NP_BIN | NP_MULTI | NP_MULTI_OLD | NP_HALF | + NP_ORTHO | NP_TRI | NP_SKIP | NP_GHOST | NP_TRIM); +// clang-format on +#else + +#ifndef LMP_NPAIR_HALFFULL_NEWTOFF_TRIM_H +#define LMP_NPAIR_HALFFULL_NEWTOFF_TRIM_H + +#include "npair.h" + +namespace LAMMPS_NS { + +class NPairHalffullNewtoffTrim : public NPair { + public: + NPairHalffullNewtoffTrim(class LAMMPS *); + void build(class NeighList *) override; +}; + +} // namespace LAMMPS_NS + +#endif +#endif diff --git a/src/npair_halffull_newton_trim.cpp b/src/npair_halffull_newton_trim.cpp new file mode 100644 index 0000000000..26e46a0525 --- /dev/null +++ b/src/npair_halffull_newton_trim.cpp @@ -0,0 +1,108 @@ +// clang-format off +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#include "npair_halffull_newton_trim.h" + +#include "atom.h" +#include "error.h" +#include "my_page.h" +#include "neigh_list.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NPairHalffullNewtonTrim::NPairHalffullNewtonTrim(LAMMPS *lmp) : NPair(lmp) {} + +/* ---------------------------------------------------------------------- + build half list from full list + pair stored once if i,j are both owned and i < j + if j is ghost, only store if j coords are "above and to the right" of i + works if full list is a skip list +------------------------------------------------------------------------- */ + +void NPairHalffullNewtonTrim::build(NeighList *list) +{ + int i,j,ii,jj,n,jnum,joriginal; + int *neighptr,*jlist; + double xtmp,ytmp,ztmp; + double delx,dely,delz,rsq; + + double **x = atom->x; + int nlocal = atom->nlocal; + + int *ilist = list->ilist; + int *numneigh = list->numneigh; + int **firstneigh = list->firstneigh; + MyPage *ipage = list->ipage; + + int *ilist_full = list->listfull->ilist; + int *numneigh_full = list->listfull->numneigh; + int **firstneigh_full = list->listfull->firstneigh; + int inum_full = list->listfull->inum; + + int inum = 0; + ipage->reset(); + + double cutsq_custom = cutoff_custom * cutoff_custom; + + // loop over parent full list + + for (ii = 0; ii < inum_full; ii++) { + n = 0; + neighptr = ipage->vget(); + + i = ilist_full[ii]; + xtmp = x[i][0]; + ytmp = x[i][1]; + ztmp = x[i][2]; + + // loop over full neighbor list + + jlist = firstneigh_full[i]; + jnum = numneigh_full[i]; + + for (jj = 0; jj < jnum; jj++) { + joriginal = jlist[jj]; + j = joriginal & NEIGHMASK; + if (j < nlocal) { + if (i > j) continue; + } else { + if (x[j][2] < ztmp) continue; + if (x[j][2] == ztmp) { + if (x[j][1] < ytmp) continue; + if (x[j][1] == ytmp && x[j][0] < xtmp) continue; + } + } + + delx = xtmp - x[j][0]; + dely = ytmp - x[j][1]; + delz = ztmp - x[j][2]; + rsq = delx * delx + dely * dely + delz * delz; + + if (rsq > cutsq_custom) continue; + + neighptr[n++] = joriginal; + } + + ilist[inum++] = i; + firstneigh[i] = neighptr; + numneigh[i] = n; + ipage->vgot(n); + if (ipage->status()) + error->one(FLERR,"Neighbor list overflow, boost neigh_modify one"); + } + + list->inum = inum; +} diff --git a/src/npair_halffull_newton_trim.h b/src/npair_halffull_newton_trim.h new file mode 100644 index 0000000000..344ed618a8 --- /dev/null +++ b/src/npair_halffull_newton_trim.h @@ -0,0 +1,44 @@ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#ifdef NPAIR_CLASS +// clang-format off +NPairStyle(halffull/newton/trim, + NPairHalffullNewtonTrim, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_MULTI_OLD | + NP_ORTHO | NP_TRI | NP_TRIM); + +NPairStyle(halffull/newton/skip/trim, + NPairHalffullNewtonTrim, + NP_HALF_FULL | NP_NEWTON | NP_HALF | NP_NSQ | NP_BIN | NP_MULTI | NP_MULTI_OLD | + NP_ORTHO | NP_TRI | NP_SKIP | NP_TRIM); +// clang-format on +#else + +#ifndef LMP_NPAIR_HALFFULL_NEWTON_TRIM_H +#define LMP_NPAIR_HALFFULL_NEWTON_TRIM_H + +#include "npair.h" + +namespace LAMMPS_NS { + +class NPairHalffullNewtonTrim : public NPair { + public: + NPairHalffullNewtonTrim(class LAMMPS *); + void build(class NeighList *) override; +}; + +} // namespace LAMMPS_NS + +#endif +#endif diff --git a/src/npair_trim.cpp b/src/npair_trim.cpp new file mode 100644 index 0000000000..5d854b745e --- /dev/null +++ b/src/npair_trim.cpp @@ -0,0 +1,92 @@ +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Trimright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#include "npair_trim.h" +#include "neigh_list.h" +#include "atom.h" +#include "error.h" +#include "my_page.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NPairTrim::NPairTrim(LAMMPS *lmp) : NPair(lmp) {} + +/* ---------------------------------------------------------------------- + create list which is a trimmed version of parent list +------------------------------------------------------------------------- */ + +void NPairTrim::build(NeighList *list) +{ + NeighList *listcopy = list->listcopy; + + double cutsq_custom = cutoff_custom * cutoff_custom; + + int i,j,ii,jj,n,jnum,joriginal; + int *neighptr,*jlist; + double xtmp,ytmp,ztmp; + double delx,dely,delz,rsq; + + double **x = atom->x; + + int *ilist = list->ilist; + int *numneigh = list->numneigh; + int **firstneigh = list->firstneigh; + MyPage *ipage = list->ipage; + ipage->reset(); + + int *ilist_copy = listcopy->ilist; + int *numneigh_copy = listcopy->numneigh; + int **firstneigh_copy = listcopy->firstneigh; + int inum = listcopy->inum; + + list->inum = inum; + list->gnum = listcopy->gnum; + + for (ii = 0; ii < inum; ii++) { + n = 0; + neighptr = ipage->vget(); + + const int i = ilist_copy[ii]; + ilist[i] = i; + xtmp = x[i][0]; + ytmp = x[i][1]; + ztmp = x[i][2]; + + // loop over copy list with larger cutoff and trim to shorter cutoff + + jlist = firstneigh_copy[i]; + jnum = numneigh_copy[i]; + + for (jj = 0; jj < jnum; jj++) { + joriginal = jlist[jj]; + j = joriginal & NEIGHMASK; + + delx = xtmp - x[j][0]; + dely = ytmp - x[j][1]; + delz = ztmp - x[j][2]; + rsq = delx * delx + dely * dely + delz * delz; + + if (rsq > cutsq_custom) continue; + + neighptr[n++] = joriginal; + } + + firstneigh[i] = neighptr; + numneigh[i] = n; + ipage->vgot(n); + if (ipage->status()) + error->one(FLERR,"Neighbor list overflow, boost neigh_modify one"); + } +} diff --git a/src/npair_trim.h b/src/npair_trim.h new file mode 100644 index 0000000000..d8510051c8 --- /dev/null +++ b/src/npair_trim.h @@ -0,0 +1,38 @@ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + Trimright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +#ifdef NPAIR_CLASS +// clang-format off +NPairStyle(trim, + NPairTrim, + NP_COPY | NP_TRIM); +// clang-format on +#else + +#ifndef LMP_NPAIR_TRIM_H +#define LMP_NPAIR_TRIM_H + +#include "npair.h" + +namespace LAMMPS_NS { + +class NPairTrim : public NPair { + public: + NPairTrim(class LAMMPS *); + void build(class NeighList *) override; +}; + +} // namespace LAMMPS_NS + +#endif +#endif diff --git a/src/pair.cpp b/src/pair.cpp index 5ce4dc2d32..4d42035fd1 100644 --- a/src/pair.cpp +++ b/src/pair.cpp @@ -97,6 +97,7 @@ Pair::Pair(LAMMPS *lmp) : Pointers(lmp) tabinner_disp = sqrt(2.0); ftable = nullptr; fdisptable = nullptr; + trim_flag = 1; allocated = 0; suffix_flag = Suffix::NONE; @@ -204,6 +205,10 @@ void Pair::modify_params(int narg, char **arg) } else if (strcmp(arg[iarg],"nofdotr") == 0) { no_virial_fdotr_compute = 1; ++iarg; + } else if (strcmp(arg[iarg],"neigh/trim") == 0) { + if (iarg+2 > narg) error->all(FLERR,"Illegal pair_modify command"); + trim_flag = utils::logical(FLERR,arg[iarg+1],false,lmp); + iarg += 2; } else error->all(FLERR,"Illegal pair_modify command"); } } diff --git a/src/pair.h b/src/pair.h index 7ef18a4fac..5f2cd00373 100644 --- a/src/pair.h +++ b/src/pair.h @@ -82,6 +82,7 @@ class Pair : protected Pointers { int tail_flag; // pair_modify flag for LJ tail correction double etail, ptail; // energy/pressure tail corrections double etail_ij, ptail_ij; + int trim_flag; // pair_modify flag for trimming neigh list int evflag; // energy,virial settings int eflag_either, eflag_global, eflag_atom; diff --git a/src/pair_hybrid.cpp b/src/pair_hybrid.cpp index dbac77a3d2..375423c321 100644 --- a/src/pair_hybrid.cpp +++ b/src/pair_hybrid.cpp @@ -35,7 +35,7 @@ using namespace LAMMPS_NS; PairHybrid::PairHybrid(LAMMPS *lmp) : Pair(lmp), styles(nullptr), keywords(nullptr), multiple(nullptr), nmap(nullptr), - map(nullptr), special_lj(nullptr), special_coul(nullptr), compute_tally(nullptr) + map(nullptr), special_lj(nullptr), special_coul(nullptr), compute_tally(nullptr), cutmax_style(nullptr) { nstyles = 0; @@ -56,6 +56,7 @@ PairHybrid::~PairHybrid() } } delete[] styles; + delete[] cutmax_style; delete[] keywords; delete[] multiple; @@ -304,6 +305,8 @@ void PairHybrid::settings(int narg, char **arg) // allocate list of sub-styles as big as possibly needed if no extra args styles = new Pair *[narg]; + cutmax_style = new double[narg]; + memset(cutmax_style, 0.0, narg*sizeof(double)); keywords = new char *[narg]; multiple = new int[narg]; @@ -722,6 +725,24 @@ double PairHybrid::init_one(int i, int j) ptail_ij += styles[map[i][j][k]]->ptail_ij; } cutmax = MAX(cutmax,cut); + + int istyle; + for (istyle = 0; istyle < nstyles; istyle++) + if (styles[istyle] == styles[map[i][j][k]]) break; + + if (styles[istyle]->trim_flag) { + + if (cut > cutmax_style[istyle]) { + cutmax_style[istyle] = cut; + + for (auto &request : neighbor->get_pair_requests()) { + if (styles[istyle] == request->get_requestor() && styles[istyle]->trim_flag) { + request->set_cutoff(cutmax_style[istyle]); + break; + } + } + } + } } return cutmax; @@ -784,6 +805,8 @@ void PairHybrid::read_restart(FILE *fp) delete[] compute_tally; styles = new Pair*[nstyles]; + cutmax_style = new double[nstyles]; + memset(cutmax_style, 0.0, nstyles*sizeof(double)); keywords = new char*[nstyles]; multiple = new int[nstyles]; diff --git a/src/pair_hybrid.h b/src/pair_hybrid.h index e66cc4d3f6..8c7782bc4f 100644 --- a/src/pair_hybrid.h +++ b/src/pair_hybrid.h @@ -69,10 +69,11 @@ class PairHybrid : public Pair { double radii2cut(double, double) override; protected: - int nstyles; // # of sub-styles - Pair **styles; // list of Pair style classes - char **keywords; // style name of each Pair style - int *multiple; // 0 if style used once, else Mth instance + int nstyles; // # of sub-styles + Pair **styles; // list of Pair style classes + double *cutmax_style; // max cutoff for each style + char **keywords; // style name of each Pair style + int *multiple; // 0 if style used once, else Mth instance int outerflag; // toggle compute() when invoked by outer() int respaflag; // 1 if different substyles are assigned to diff --git a/src/pair_hybrid_scaled.cpp b/src/pair_hybrid_scaled.cpp index 2c0c574f70..a059deb3fb 100644 --- a/src/pair_hybrid_scaled.cpp +++ b/src/pair_hybrid_scaled.cpp @@ -276,6 +276,8 @@ void PairHybridScaled::settings(int narg, char **arg) // allocate list of sub-styles as big as possibly needed if no extra args styles = new Pair *[narg]; + cutmax_style = new double[narg]; + memset(cutmax_style, 0.0, narg*sizeof(double)); keywords = new char *[narg]; multiple = new int[narg];