diff --git a/doc/src/Commands_compute.rst b/doc/src/Commands_compute.rst index 9066531459..19e56fcc83 100644 --- a/doc/src/Commands_compute.rst +++ b/doc/src/Commands_compute.rst @@ -91,6 +91,7 @@ KOKKOS, o = OPENMP, t = OPT. * :doc:`ke/atom/eff ` * :doc:`ke/eff ` * :doc:`ke/rigid ` + * :doc:`local/comp/atom (k) ` * :doc:`mliap ` * :doc:`momentum ` * :doc:`msd ` diff --git a/doc/src/compute.rst b/doc/src/compute.rst index 4e5d8e6cb4..317efa4f41 100644 --- a/doc/src/compute.rst +++ b/doc/src/compute.rst @@ -245,6 +245,7 @@ The individual style names on the :doc:`Commands compute ` pag * :doc:`ke/atom/eff ` - per-atom translational and radial kinetic energy in the electron force field model * :doc:`ke/eff ` - kinetic energy of a group of nuclei and electrons in the electron force field model * :doc:`ke/rigid ` - translational kinetic energy of rigid bodies +* :doc:`local/comp/atom ` - local composition for each atom * :doc:`mliap ` - gradients of energy and forces with respect to model parameters and related quantities for training machine learning interatomic potentials * :doc:`momentum ` - translational momentum * :doc:`msd ` - mean-squared displacement of group of atoms diff --git a/doc/src/compute_local_comp_atom.rst b/doc/src/compute_local_comp_atom.rst new file mode 100644 index 0000000000..2d0dd8c6f2 --- /dev/null +++ b/doc/src/compute_local_comp_atom.rst @@ -0,0 +1,118 @@ +.. index:: compute local/comp/atom +.. index:: compute local/comp/atom/kk + +compute local/comp/atom command +=============================== + +Accelerator Variants: *local/comp/atom/kk* + +Syntax +"""""" + +.. code-block:: LAMMPS + + compute ID group-ID local/comp/atom keyword values ... + +* ID, group-ID are documented in :doc:`compute ` command +* local/comp/atom = style name of this compute command +* one or more keyword/value pairs may be appended + + .. parsed-literal:: + + keyword = *cutoff* + *cutoff* value = distance cutoff + +Examples +"""""""" + +.. code-block:: LAMMPS + + compute 1 all local/comp/atom + + compute 1 all local/comp/atom cutoff 9.0 + comm_modify cutoff 9.0 + + +Description +""""""""""" + +.. versionadded:: TBD + +Define a computation that calculates a local composition vector for each +atom. For a central atom with :math:`M` neighbors within the neighbor cutoff sphere, +composition is defined as the number of atoms of a given type +(including the central atom) divided by (:math:`M+1`). For a given central atom, +the sum of all compositions equals one. + +.. note:: + + This compute uses the number of atom types, not chemical species, assigned in + :doc:`pair_coeff ` command. If an interatomic potential has two + species (i.e., Cu and Ni) assigned to four different atom types in + :doc:`pair_coeff ` (i.e., 'Cu Cu Ni Ni'), the compute will + output four fractional values. In those cases, the user may desire an extra + calculation step to consolidate per-type fractions into per-species fractions. + This calculation can be conducted within LAMMPS using another compute such as + :doc:`compute reduce `, an atom-style :doc:`variable`, or as a + post-processing step. + +---------- + +The optional keyword *cutoff* defines the distance cutoff used when +searching for neighbors. The default value is the cutoff specified by +the pair style. If no pair style is defined, then a cutoff must be +defined using this keyword. If the specified cutoff is larger than +that of the pair_style plus neighbor skin (or no pair style is +defined), the *comm_modify cutoff* option must also be set to match +that of the *cutoff* keyword. + +The neighbor list needed to compute this quantity is constructed each +time the calculation is performed (i.e. each time a snapshot of atoms +is dumped). Thus it can be inefficient to compute/dump this quantity +too frequently. + +.. note:: + + If you have a bonded system, then the settings of + :doc:`special_bonds ` command can remove pairwise + interactions between atoms in the same bond, angle, or dihedral. + This is the default setting for the :doc:`special_bonds + ` command, and means those pairwise interactions do + not appear in the neighbor list. Because this compute uses the + neighbor list, it also means those pairs will not be included in + the order parameter. This difficulty can be circumvented by + writing a dump file, and using the :doc:`rerun ` command to + compute the order parameter for snapshots in the dump file. The + rerun script can use a :doc:`special_bonds ` command + that includes all pairs in the neighbor list. + +---------- + +Output info +""""""""""" + +This compute calculates a per-atom array with :math:`1 + N` columns, where :math:`N` +is the number of atom types. The first column is a count of the number of atoms +used to calculate composition (including the central atom), and each subsequent +column indicates the fraction of that atom type within the cutoff sphere. + +These values can be accessed by any command that uses per-atom values +from a compute as input. See the :doc:`Howto output ` +doc page for an overview of LAMMPS output options. + +Restrictions +"""""""""""" + +This compute is part of the EXTRA-COMPUTE package. It is only enabled +if LAMMPS was built with that package. See the :doc:`Build package +` page for more info. + +Related commands +"""""""""""""""" + +:doc:`comm_modify ` + +Default +""""""" + +The option defaults are *cutoff* = pair style cutoff. diff --git a/doc/utils/sphinx-config/false_positives.txt b/doc/utils/sphinx-config/false_positives.txt index 7352a921ef..1643b51fec 100644 --- a/doc/utils/sphinx-config/false_positives.txt +++ b/doc/utils/sphinx-config/false_positives.txt @@ -3621,6 +3621,7 @@ Tk Tkin tloop tlsph +tm tmax Tmax tmd diff --git a/src/.gitignore b/src/.gitignore index 6e2a8b92c9..115291c39a 100644 --- a/src/.gitignore +++ b/src/.gitignore @@ -580,6 +580,8 @@ /compute_ke_eff.h /compute_ke_rigid.cpp /compute_ke_rigid.h +/compute_local_comp_atom.cpp +/compute_local_comp_atom.h /compute_meso_e_atom.cpp /compute_meso_e_atom.h /compute_meso_rho_atom.cpp diff --git a/src/EXTRA-COMPUTE/compute_local_comp_atom.cpp b/src/EXTRA-COMPUTE/compute_local_comp_atom.cpp new file mode 100644 index 0000000000..77c4993996 --- /dev/null +++ b/src/EXTRA-COMPUTE/compute_local_comp_atom.cpp @@ -0,0 +1,214 @@ +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + LAMMPS development team: developers@lammps.org + + 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: Megan McCarthy (SNL) +------------------------------------------------------------------------- */ + +#include "compute_local_comp_atom.h" + +#include "atom.h" +#include "comm.h" +#include "domain.h" +#include "error.h" +#include "force.h" +#include "math_const.h" +#include "memory.h" +#include "neigh_list.h" +#include "neigh_request.h" +#include "neighbor.h" +#include "pair.h" +#include "update.h" + +#include + +using namespace LAMMPS_NS; +using namespace MathConst; + +/* ---------------------------------------------------------------------- */ + +ComputeLocalCompAtom::ComputeLocalCompAtom(LAMMPS *lmp, int narg, char **arg) : + Compute(lmp, narg, arg), result(nullptr) +{ + if (narg < 3 || narg > 5) error->all(FLERR, "Illegal compute local/comp/atom command"); + + cutoff = 0.0; + + int iarg = 3; + while (iarg < narg) { + if (strcmp(arg[iarg], "cutoff") == 0) { + if (iarg + 2 > narg) error->all(FLERR, "Illegal compute local/comp/atom command"); + cutoff = utils::numeric(FLERR, arg[iarg + 1], false, lmp); + if (cutoff <= 0.0) error->all(FLERR, "Illegal compute local/comp/atom command"); + iarg += 2; + } else + error->all(FLERR, "Illegal compute local/comp/atom command"); + } + + peratom_flag = 1; + + ntypes = atom->ntypes; + size_peratom_cols = 1 + ntypes; + + nmax = 0; +} + +/* ---------------------------------------------------------------------- */ + +ComputeLocalCompAtom::~ComputeLocalCompAtom() +{ + if (copymode) return; + + memory->destroy(result); +} + +/* ---------------------------------------------------------------------- */ + +void ComputeLocalCompAtom::init() +{ + if (!force->pair && cutoff == 0.0) + error->all(FLERR, + "Compute local/comp/atom requires a cutoff be specified " + "or a pair style be defined"); + + double skin = neighbor->skin; + if (cutoff != 0.0) { + double cutghost; // as computed by Neighbor and Comm + if (force->pair) + cutghost = MAX(force->pair->cutforce + skin, comm->cutghostuser); + else + cutghost = comm->cutghostuser; + + if (cutoff > cutghost) + error->all(FLERR, + "Compute local/comp/atom cutoff exceeds ghost atom range - " + "use comm_modify cutoff command"); + } + + int cutflag = 1; + if (force->pair) { + if (cutoff == 0.0) { cutoff = force->pair->cutforce; } + if (cutoff <= force->pair->cutforce + skin) cutflag = 0; + } + + cutsq = cutoff * cutoff; + + // need an occasional full neighbor list + + auto req = neighbor->add_request(this, NeighConst::REQ_FULL | NeighConst::REQ_OCCASIONAL); + if (cutflag) req->set_cutoff(cutoff); +} + +/* ---------------------------------------------------------------------- */ + +void ComputeLocalCompAtom::init_list(int /*id*/, NeighList *ptr) +{ + list = ptr; +} + +/* ---------------------------------------------------------------------- */ + +void ComputeLocalCompAtom::compute_peratom() +{ + int i, j, ii, jj, inum, jnum; + double xtmp, ytmp, ztmp, delx, dely, delz, rsq; + int *ilist, *jlist, *numneigh, **firstneigh; + int count, itype, jtype; + + invoked_peratom = update->ntimestep; + + // grow result array if necessary + + if (atom->nmax > nmax) { + memory->destroy(result); + nmax = atom->nmax; + memory->create(result, nmax, size_peratom_cols, "local/comp/atom:result"); + array_atom = result; + } + + // invoke full neighbor list (will copy or build if necessary) + + neighbor->build_one(list); + + inum = list->inum; + ilist = list->ilist; + numneigh = list->numneigh; + firstneigh = list->firstneigh; + + // compute properties for each atom in group + // use full neighbor list to count atoms less than cutoff + + double **x = atom->x; + int *type = atom->type; + int *mask = atom->mask; + + // get per-atom local compositions + + for (ii = 0; ii < inum; ii++) { + + i = ilist[ii]; + + if (mask[i] & groupbit) { + + xtmp = x[i][0]; + ytmp = x[i][1]; + ztmp = x[i][2]; + jlist = firstneigh[i]; + jnum = numneigh[i]; + + // i atom contribution + + count = 1; + + itype = type[i]; + result[i][itype]++; + + for (jj = 0; jj < jnum; jj++) { + j = jlist[jj]; + j &= NEIGHMASK; + + jtype = type[j]; + + 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) { + count++; + result[i][jtype]++; + } + } + + // total count of atoms found in sampled radius range + + result[i][0] = count; + + // local comp fractions per element + + double lfac = 1.0 / count; + for (int n = 1; n < size_peratom_cols; n++) + result[i][n+1] *= lfac; + + } + } +} + +/* ---------------------------------------------------------------------- + memory usage of local atom-based array +------------------------------------------------------------------------- */ + +double ComputeLocalCompAtom::memory_usage() +{ + double bytes = (double) 2 * nmax * sizeof(double); + return bytes; +} diff --git a/src/EXTRA-COMPUTE/compute_local_comp_atom.h b/src/EXTRA-COMPUTE/compute_local_comp_atom.h new file mode 100644 index 0000000000..76d92de08d --- /dev/null +++ b/src/EXTRA-COMPUTE/compute_local_comp_atom.h @@ -0,0 +1,47 @@ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + LAMMPS development team: developers@lammps.org + + 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 COMPUTE_CLASS +// clang-format off +ComputeStyle(local/comp/atom,ComputeLocalCompAtom); +// clang-format on +#else + +#ifndef LMP_COMPUTE_LOCAL_COMP_ATOM_H +#define LMP_COMPUTE_LOCAL_COMP_ATOM_H + +#include "compute.h" + +namespace LAMMPS_NS { + +class ComputeLocalCompAtom : public Compute { + public: + ComputeLocalCompAtom(class LAMMPS *, int, char **); + ~ComputeLocalCompAtom() override; + void init() override; + void init_list(int, class NeighList *) override; + void compute_peratom() override; + double memory_usage() override; + + protected: + int nmax, ntypes; + double cutoff; // global cutoff distance + double cutsq; // cutoff**2 + class NeighList *list; // neighbor list + double **result; // peratom array of local compositions +}; + +} // namespace LAMMPS_NS + +#endif +#endif diff --git a/src/KOKKOS/Install.sh b/src/KOKKOS/Install.sh index ede766cbf8..a842a169c0 100755 --- a/src/KOKKOS/Install.sh +++ b/src/KOKKOS/Install.sh @@ -96,6 +96,8 @@ action compute_coord_atom_kokkos.cpp action compute_coord_atom_kokkos.h action compute_erotate_sphere_kokkos.cpp action compute_erotate_sphere_kokkos.h +action compute_local_comp_atom_kokkos.cpp compute_local_comp_atom.cpp +action compute_local_comp_atom_kokkos.h compute_local_comp_atom.h action compute_orientorder_atom_kokkos.cpp action compute_orientorder_atom_kokkos.h action compute_temp_deform_kokkos.cpp diff --git a/src/KOKKOS/compute_local_comp_atom_kokkos.cpp b/src/KOKKOS/compute_local_comp_atom_kokkos.cpp new file mode 100644 index 0000000000..356efea2f1 --- /dev/null +++ b/src/KOKKOS/compute_local_comp_atom_kokkos.cpp @@ -0,0 +1,180 @@ +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + LAMMPS development team: developers@lammps.org + + 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 Gdirectoryeneral Public License. + + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +/* ---------------------------------------------------------------------- + Contributing author: Megan McCarthy (SNL) +------------------------------------------------------------------------- */ + +#include "compute_local_comp_atom_kokkos.h" + +#include "atom_kokkos.h" +#include "atom_masks.h" +#include "comm.h" +#include "domain.h" +#include "error.h" +#include "force.h" +#include "memory_kokkos.h" +#include "modify.h" +#include "neigh_list.h" +#include "neigh_request.h" +#include "neighbor_kokkos.h" +#include "pair.h" +#include "update.h" + +#include +#include + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +template +ComputeLocalCompAtomKokkos::ComputeLocalCompAtomKokkos(LAMMPS *lmp, int narg, char **arg) : + ComputeLocalCompAtom(lmp, narg, arg) +{ + kokkosable = 1; + atomKK = (AtomKokkos *) atom; + execution_space = ExecutionSpaceFromDevice::space; + datamask_read = EMPTY_MASK; + datamask_modify = EMPTY_MASK; +} + +/* ---------------------------------------------------------------------- */ + +template +ComputeLocalCompAtomKokkos::~ComputeLocalCompAtomKokkos() +{ + if (copymode) return; + + memoryKK->destroy_kokkos(k_result,result); +} + +/* ---------------------------------------------------------------------- */ + +template +void ComputeLocalCompAtomKokkos::init() +{ + ComputeLocalCompAtom::init(); + + // adjust neighbor list request for KOKKOS + + auto request = neighbor->find_request(this); + request->set_kokkos_host(std::is_same::value && + !std::is_same::value); + request->set_kokkos_device(std::is_same::value); +} + +/* ---------------------------------------------------------------------- */ + +template +void ComputeLocalCompAtomKokkos::compute_peratom() +{ + invoked_peratom = update->ntimestep; + + // grow result array if necessary + + int size_peratom_cols = 1 + atom->ntypes; + if (atom->nmax > nmax) { + memoryKK->destroy_kokkos(k_result,result); + nmax = atom->nmax; + memoryKK->create_kokkos(k_result,result,nmax,size_peratom_cols,"local/comp/atom:result"); + d_result = k_result.view(); + array_atom = result; + } + + // invoke full neighbor list (will copy or build if necessary) + + neighbor->build_one(list); + int inum = list->inum; + + NeighListKokkos* k_list = static_cast*>(list); + d_numneigh = k_list->d_numneigh; + d_neighbors = k_list->d_neighbors; + d_ilist = k_list->d_ilist; + + // compute properties for each atom in group + // use full neighbor list to count atoms less than cutoff + + atomKK->sync(execution_space,X_MASK|TYPE_MASK|MASK_MASK); + x = atomKK->k_x.view(); + type = atomKK->k_type.view(); + mask = atomKK->k_mask.view(); + ntypes = atom->ntypes; + + Kokkos::deep_copy(d_result,0.0); + + copymode = 1; + typename Kokkos::RangePolicy policy(0,inum); + Kokkos::parallel_for("ComputeLocalComp",policy,*this); + copymode = 0; + + k_result.modify(); + k_result.sync_host(); +} + +template +KOKKOS_INLINE_FUNCTION +void ComputeLocalCompAtomKokkos::operator()(TagComputeLocalCompAtom, const int &ii) const +{ + const int i = d_ilist[ii]; + + if (mask[i] & groupbit) { + + const X_FLOAT xtmp = x(i,0); + const X_FLOAT ytmp = x(i,1); + const X_FLOAT ztmp = x(i,2); + const int jnum = d_numneigh[i]; + + // i atom contribution + + int count = 1.0; + int itype = type[i]; + d_result(i,itype)++; + + for (int jj = 0; jj < jnum; jj++) { + + int j = d_neighbors(i,jj); + j &= NEIGHMASK; + + int jtype = type[j]; + + const F_FLOAT delx = x(j,0) - xtmp; + const F_FLOAT dely = x(j,1) - ytmp; + const F_FLOAT delz = x(j,2) - ztmp; + const F_FLOAT rsq = delx*delx + dely*dely + delz*delz; + if (rsq < cutsq) { + count++; + d_result(i,jtype) += 1.0; + } + } + + // total count of atoms found in sampled radius range + + d_result(i,0) = count; + + // local comp fractions per atom type + + double lfac = 1.0 / count; + + for (int n = 1; n < size_peratom_cols; n++) { + d_result(i,n) *= lfac; + } + } +} + +namespace LAMMPS_NS { +template class ComputeLocalCompAtomKokkos; +#ifdef LMP_KOKKOS_GPU +template class ComputeLocalCompAtomKokkos; +#endif +} diff --git a/src/KOKKOS/compute_local_comp_atom_kokkos.h b/src/KOKKOS/compute_local_comp_atom_kokkos.h new file mode 100644 index 0000000000..5fb52474b2 --- /dev/null +++ b/src/KOKKOS/compute_local_comp_atom_kokkos.h @@ -0,0 +1,64 @@ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + LAMMPS development team: developers@lammps.org + + 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 COMPUTE_CLASS +// clang-format off +ComputeStyle(local/comp/atom/kk,ComputeLocalCompAtomKokkos); +ComputeStyle(local/comp/atom/kk/device,ComputeLocalCompAtomKokkos); +ComputeStyle(local/comp/atom/kk/host,ComputeLocalCompAtomKokkos); +// clang-format on + +#else + +#ifndef LMP_COMPUTE_LOCAL_COMP_ATOM_KOKKOS_H +#define LMP_COMPUTE_LOCAL_COMP_ATOM_KOKKOS_H + +#include "compute_local_comp_atom.h" +#include "kokkos_type.h" + +namespace LAMMPS_NS { + +// clang-format off +struct TagComputeLocalCompAtom {}; +// clang-format on + +template class ComputeLocalCompAtomKokkos : public ComputeLocalCompAtom { + public: + typedef DeviceType device_type; + typedef ArrayTypes AT; + + ComputeLocalCompAtomKokkos(class LAMMPS *, int, char **); + ~ComputeLocalCompAtomKokkos() override; + void init() override; + void compute_peratom() override; + + KOKKOS_INLINE_FUNCTION + void operator()(TagComputeLocalCompAtom, const int &) const; + + private: + + typename AT::t_x_array x; + typename ArrayTypes::t_int_1d type; + typename ArrayTypes::t_int_1d mask; + + typename AT::t_neighbors_2d d_neighbors; + typename AT::t_int_1d d_ilist; + typename AT::t_int_1d d_numneigh; + DAT::tdual_float_2d k_result; + typename AT::t_float_2d d_result; +}; + +} // namespace LAMMPS_NS + +#endif +#endif diff --git a/unittest/formats/CMakeLists.txt b/unittest/formats/CMakeLists.txt index 93ea2f3b32..58c797b6e6 100644 --- a/unittest/formats/CMakeLists.txt +++ b/unittest/formats/CMakeLists.txt @@ -41,6 +41,8 @@ set_tests_properties(TextFileReader PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${ add_executable(test_file_operations test_file_operations.cpp) target_link_libraries(test_file_operations PRIVATE lammps GTest::GMock) add_test(NAME FileOperations COMMAND test_file_operations) +# try to mitigate possible OpenMPI bug +set_tests_properties(TextFileReader PROPERTIES ENVIRONMENT "OMPI_MCA_sharedfp=\"^sm\"") add_executable(test_dump_atom test_dump_atom.cpp) target_link_libraries(test_dump_atom PRIVATE lammps GTest::GMock)