From b00a281708aa16e5db6d96b58df0a91a97828b17 Mon Sep 17 00:00:00 2001 From: Donatas Surblys Date: Fri, 16 Jul 2021 14:29:07 +0900 Subject: [PATCH 1/6] a tally compute to obtain virial heat flux into group 1 due to group 2 --- src/TALLY/compute_heat_flux_virial_tally.cpp | 245 +++++++++++++++++++ src/TALLY/compute_heat_flux_virial_tally.h | 64 +++++ 2 files changed, 309 insertions(+) create mode 100644 src/TALLY/compute_heat_flux_virial_tally.cpp create mode 100644 src/TALLY/compute_heat_flux_virial_tally.h diff --git a/src/TALLY/compute_heat_flux_virial_tally.cpp b/src/TALLY/compute_heat_flux_virial_tally.cpp new file mode 100644 index 0000000000..0ac3607da6 --- /dev/null +++ b/src/TALLY/compute_heat_flux_virial_tally.cpp @@ -0,0 +1,245 @@ +// 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 "compute_heat_flux_virial_tally.h" + +#include +#include "atom.h" +#include "group.h" +#include "pair.h" +#include "update.h" +#include "memory.h" +#include "error.h" +#include "force.h" +#include "comm.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +ComputeHeatFluxVirialTally::ComputeHeatFluxVirialTally(LAMMPS *lmp, int narg, char **arg) : + Compute(lmp, narg, arg) +{ + if (narg < 4) error->all(FLERR,"Illegal compute heat/flux/virial/tally command"); + + igroup2 = group->find(arg[3]); + if (igroup2 == -1) + error->all(FLERR,"Could not find compute heat/flux/virial/tally second group ID"); + groupbit2 = group->bitmask[igroup2]; + + scalar_flag = 1; + vector_flag = 0; + peratom_flag = 1; + timeflag = 1; + dynamic_group_allow = 0; + + comm_reverse = size_peratom_cols = 3; + extscalar = 1; + peflag = 1; // we need Pair::ev_tally() to be run + + did_setup = invoked_peratom = invoked_scalar = -1; + nmax = -1; + fatom = nullptr; +} + +/* ---------------------------------------------------------------------- */ + +ComputeHeatFluxVirialTally::~ComputeHeatFluxVirialTally() +{ + if (force && force->pair) force->pair->del_tally_callback(this); + memory->destroy(fatom); +} + +/* ---------------------------------------------------------------------- */ + +void ComputeHeatFluxVirialTally::init() +{ + if (force->pair == nullptr) + error->all(FLERR,"Trying to use compute heat/flux/virial/tally without pair style"); + else + force->pair->add_tally_callback(this); + + if (comm->me == 0) { + if (force->pair->single_enable == 0 || force->pair->manybody_flag) + error->warning(FLERR,"Compute heat/flux/virial/tally used with incompatible pair style"); + + if (force->bond || force->angle || force->dihedral + || force->improper || force->kspace) + error->warning(FLERR,"Compute heat/flux/virial/tally only called from pair style"); + } + + // error out if any atoms are in both groups + for (int i = 0; i < atom->nlocal; i++) + { + if ((atom->mask[i] & groupbit) && (atom->mask[i] & groupbit2)) { + if (atom->tag_enable) { + error->all(FLERR,"Atom {} belonging to both groups is not allowed" + " with compute heat/flux/virial/tally", atom->tag[i]); + } else { + error->all(FLERR,"Atom belonging to both groups is not allowed" + " with compute heat/flux/virial/tally"); + } + } + } + + did_setup = -1; +} + +/* ---------------------------------------------------------------------- */ + +void ComputeHeatFluxVirialTally::pair_setup_callback(int, int) +{ + // run setup only once per time step. + // we may be called from multiple pair styles + + if (did_setup == update->ntimestep) return; + + const int ntotal = atom->nlocal + atom->nghost; + + // grow per-atom storage, if needed + + if (atom->nmax > nmax) { + memory->destroy(fatom); + nmax = atom->nmax; + memory->create(fatom,nmax,size_peratom_cols,"heat/flux/virial/tally:fatom"); + array_atom = fatom; + } + + // clear storage + + for (int i=0; i < ntotal; ++i) + for (int j=0; j < size_peratom_cols; ++j) + fatom[i][j] = 0.0; + + did_setup = update->ntimestep; +} + +/* ---------------------------------------------------------------------- */ +void ComputeHeatFluxVirialTally::pair_tally_callback(int i, int j, int nlocal, int newton, + double, double, double fpair, + double dx, double dy, double dz) +{ + const int * const mask = atom->mask; + + const bool ingroup1 = (mask[i] & groupbit); + if ( (ingroup1 && (mask[j] & groupbit2)) + || ((mask[i] & groupbit2) && (mask[j] & groupbit)) ) { + + // signs set to calculate heat flux from group2 into group1 + const double fx = (ingroup1?0.5:-0.5)*dx*fpair; + const double fy = (ingroup1?0.5:-0.5)*dy*fpair; + const double fz = (ingroup1?0.5:-0.5)*dz*fpair; + + if (newton || i < nlocal) { + fatom[i][0] += fx; + fatom[i][1] += fy; + fatom[i][2] += fz; + } + if (newton || j < nlocal) { + fatom[j][0] += fx; + fatom[j][1] += fy; + fatom[j][2] += fz; + } + } +} + +/* ---------------------------------------------------------------------- */ + +int ComputeHeatFluxVirialTally::pack_reverse_comm(int n, int first, double *buf) +{ + int i,m,last; + + m = 0; + last = first + n; + for (i = first; i < last; i++) { + buf[m++] = fatom[i][0]; + buf[m++] = fatom[i][1]; + buf[m++] = fatom[i][2]; + } + return m; +} + +/* ---------------------------------------------------------------------- */ + +void ComputeHeatFluxVirialTally::unpack_reverse_comm(int n, int *list, double *buf) +{ + int i,j,m; + + m = 0; + for (i = 0; i < n; i++) { + j = list[i]; + fatom[j][0] += buf[m++]; + fatom[j][1] += buf[m++]; + fatom[j][2] += buf[m++]; + } +} + +/* ---------------------------------------------------------------------- */ + +double ComputeHeatFluxVirialTally::compute_scalar() +{ + // need to collect contributions from ghost atoms + // because we need to use local velocities to compute heat flux + if (invoked_peratom != update->ntimestep) + compute_peratom(); + + invoked_scalar = update->ntimestep; + if ((did_setup != invoked_scalar) + || (update->eflag_global != invoked_scalar)) + error->all(FLERR,"Energy was not tallied on needed timestep"); + + // sum heat flux across procs + double hflux = 0.0; + for (int i = 0; i < atom->nlocal; i++) + hflux += fatom[i][0]*atom->v[i][0] + + fatom[i][1]*atom->v[i][1] + + fatom[i][2]*atom->v[i][2]; + + MPI_Allreduce(&hflux,&scalar,1,MPI_DOUBLE,MPI_SUM,world); + + return scalar; +} + +/* ---------------------------------------------------------------------- */ + +void ComputeHeatFluxVirialTally::compute_peratom() +{ + invoked_peratom = update->ntimestep; + if ((did_setup != invoked_peratom) + || (update->eflag_global != invoked_peratom)) + error->all(FLERR,"Energy was not tallied on needed timestep"); + + // collect contributions from ghost atoms + + if (force->newton_pair) { + comm->reverse_comm_compute(this); + + // clear out ghost atom data after it has been collected to local atoms + const int nall = atom->nlocal + atom->nghost; + for (int i = atom->nlocal; i < nall; ++i) + for (int j = 0; j < size_peratom_cols; ++j) + fatom[i][j] = 0.0; + } +} + +/* ---------------------------------------------------------------------- + memory usage of local atom-based array +------------------------------------------------------------------------- */ + +double ComputeHeatFluxVirialTally::memory_usage() +{ + double bytes = (nmax < 0) ? 0 : nmax*size_peratom_cols * sizeof(double); + return bytes; +} + diff --git a/src/TALLY/compute_heat_flux_virial_tally.h b/src/TALLY/compute_heat_flux_virial_tally.h new file mode 100644 index 0000000000..d2a3d39a7d --- /dev/null +++ b/src/TALLY/compute_heat_flux_virial_tally.h @@ -0,0 +1,64 @@ +/* -*- 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 COMPUTE_CLASS +// clang-format off +ComputeStyle(heat/flux/virial/tally,ComputeHeatFluxVirialTally); +// clang-format on +#else + +#ifndef LMP_COMPUTE_HEAT_FLUX_VIRIAL_TALLY_H +#define LMP_COMPUTE_HEAT_FLUX_VIRIAL_TALLY_H + +#include "compute.h" + +namespace LAMMPS_NS { + +class ComputeHeatFluxVirialTally : public Compute { + + public: + ComputeHeatFluxVirialTally(class LAMMPS *, int, char **); + virtual ~ComputeHeatFluxVirialTally(); + + void init(); + + double compute_scalar(); + void compute_peratom(); + + int pack_reverse_comm(int, int, double *); + void unpack_reverse_comm(int, int *, double *); + double memory_usage(); + + void pair_setup_callback(int, int); + void pair_tally_callback(int, int, int, int, double, double, double, double, double, double); + + private: + bigint did_setup; + int nmax, igroup2, groupbit2; + double **fatom; +}; + +} // namespace LAMMPS_NS + +#endif +#endif + +/* ERROR/WARNING messages: + +E: Illegal ... command + +Self-explanatory. Check the input script syntax and compare to the +documentation for the command. You can use -echo screen as a +command-line option when running LAMMPS to see the offending line. + +*/ From cfd9e31d18f11cd10a9c67856189de2e55294979 Mon Sep 17 00:00:00 2001 From: Donatas Surblys Date: Tue, 20 Jul 2021 13:31:41 +0900 Subject: [PATCH 2/6] add links to heat/flux/virial/tally and add short descriptions for tally computes --- doc/src/Commands_compute.rst | 1 + doc/src/compute.rst | 11 ++++++----- 2 files changed, 7 insertions(+), 5 deletions(-) diff --git a/doc/src/Commands_compute.rst b/doc/src/Commands_compute.rst index 0169826751..9dfb28fa8b 100644 --- a/doc/src/Commands_compute.rst +++ b/doc/src/Commands_compute.rst @@ -72,6 +72,7 @@ KOKKOS, o = OPENMP, t = OPT. * :doc:`gyration/shape/chunk ` * :doc:`heat/flux ` * :doc:`heat/flux/tally ` + * :doc:`heat/flux/virial/tally ` * :doc:`hexorder/atom ` * :doc:`hma ` * :doc:`improper ` diff --git a/doc/src/compute.rst b/doc/src/compute.rst index c2b2b92aa1..93916d701d 100644 --- a/doc/src/compute.rst +++ b/doc/src/compute.rst @@ -208,7 +208,7 @@ The individual style names on the :doc:`Commands compute ` doc * :doc:`event/displace ` - detect event on atom displacement * :doc:`fabric ` - calculates fabric tensors from pair interactions * :doc:`fep ` - -* :doc:`force/tally ` - +* :doc:`force/tally ` - force between two groups of atoms via the tally callback mechanism * :doc:`fragment/atom ` - fragment ID for each atom * :doc:`global/atom ` - * :doc:`group/group ` - energy/force between two groups of atoms @@ -217,7 +217,8 @@ The individual style names on the :doc:`Commands compute ` doc * :doc:`gyration/shape ` - shape parameters from gyration tensor * :doc:`gyration/shape/chunk ` - shape parameters from gyration tensor for each chunk * :doc:`heat/flux ` - heat flux through a group of atoms -* :doc:`heat/flux/tally ` - +* :doc:`heat/flux/tally ` - heat flux through a group of atoms via the tally callback mechanism +* :doc:`heat/flux/virial/tally ` - virial heat flux between two groups via the tally callback mechanism * :doc:`hexorder/atom ` - bond orientational order parameter q6 * :doc:`hma ` - harmonically mapped averaging for atomic crystals * :doc:`improper ` - energy of each improper sub-style @@ -240,8 +241,8 @@ The individual style names on the :doc:`Commands compute ` doc * :doc:`pe ` - potential energy * :doc:`pe/atom ` - potential energy for each atom * :doc:`mesont ` - Nanotube bending,stretching, and intertube energies -* :doc:`pe/mol/tally ` - -* :doc:`pe/tally ` - +* :doc:`pe/mol/tally ` - potential energy between two groups of atoms seperated into intermolecular and intramolecular components via the tally callback mechanism +* :doc:`pe/tally ` - potential energy between two groups of atoms via the tally callback mechanism * :doc:`plasticity/atom ` - Peridynamic plasticity for each atom * :doc:`pressure ` - total pressure and pressure tensor * :doc:`pressure/cylinder ` - pressure tensor in cylindrical coordinates @@ -289,7 +290,7 @@ The individual style names on the :doc:`Commands compute ` doc * :doc:`stress/atom ` - stress tensor for each atom * :doc:`stress/mop ` - normal components of the local stress tensor using the method of planes * :doc:`stress/mop/profile ` - profile of the normal components of the local stress tensor using the method of planes -* :doc:`stress/tally ` - +* :doc:`stress/tally ` - stress between two groups of atoms via the tally callback mechanism * :doc:`tdpd/cc/atom ` - per-atom chemical concentration of a specified species for each tDPD particle * :doc:`temp ` - temperature of group of atoms * :doc:`temp/asphere ` - temperature of aspherical particles From 6fd1cda2a603407308e111e947a0c8b86458e4b5 Mon Sep 17 00:00:00 2001 From: Donatas Surblys Date: Tue, 20 Jul 2021 17:48:12 +0900 Subject: [PATCH 3/6] update documentaion for heat/flux/tally and heat/flux/virial/tally --- doc/src/compute_tally.rst | 96 +++++++++++++++++++++++++++++++++++++-- 1 file changed, 92 insertions(+), 4 deletions(-) diff --git a/doc/src/compute_tally.rst b/doc/src/compute_tally.rst index 32e3e31030..d2e4569820 100644 --- a/doc/src/compute_tally.rst +++ b/doc/src/compute_tally.rst @@ -1,5 +1,6 @@ .. index:: compute force/tally .. index:: compute heat/flux/tally +.. index:: compute heat/flux/virial/tally .. index:: compute pe/tally .. index:: compute pe/mol/tally .. index:: compute stress/tally @@ -10,6 +11,9 @@ compute force/tally command compute heat/flux/tally command =============================== +compute heat/flux/virial/tally command +====================================== + compute pe/tally command ======================== @@ -27,7 +31,7 @@ Syntax compute ID group-ID style group2-ID * ID, group-ID are documented in :doc:`compute ` command -* style = *force/tally* or *pe/tally* or *pe/mol/tally* or *stress/tally* +* style = *force/tally* or *heat/flux/tally* or *heat/flux/virial/tally* or * or *pe/tally* or *pe/mol/tally* or *stress/tally* * group2-ID = group ID of second (or same) group Examples @@ -38,13 +42,17 @@ Examples compute 1 lower force/tally upper compute 1 left pe/tally right compute 1 lower stress/tally lower + compute 1 subregion heat/flux/tally all + compute 1 liquid heat/flux/virial/tally solid Description """"""""""" Define a computation that calculates properties between two groups of -atoms by accumulating them from pairwise non-bonded computations. The -two groups can be the same. This is similar to :doc:`compute group/group ` only that the data is +atoms by accumulating them from pairwise non-bonded computations. +Except for *heat/flux/virial/tally*, the two groups can be the same. +This is similar to :doc:`compute group/group ` +only that the data is accumulated directly during the non-bonded force computation. The computes *force/tally*\ , *pe/tally*\ , *stress/tally*\ , and *heat/flux/tally* are primarily provided as example how to program @@ -57,6 +65,76 @@ the based classes of LAMMPS. ---------- +Compute *heat/flux/tally* obtains the heat flux +(strictly speaking, heat flow) inside the first group, +which is the sum of the convective contribution +due to atoms in the first group and the virial contribution +due to interaction between the first and second groups: + +.. math:: + + \mathbf{Q}= \sum_{i \in \text{group 1}} e_i \mathbf{v}_i + \frac{1}{2} \sum_{i \in \text{group 1}} \sum_{\substack{j \in \text{group 2} \\ j \neq i } } \left( \mathbf{F}_{ij} \cdot \mathbf{v}_j \right) \mathbf{r}_{ij} + +When the second group in *heat/flux/tally* is set to "all", +the resulting values will be identical +to that obtained by :doc:`compute heat/flux `, +provided only pairwise interactions exist. + +Compute *heat/flux/virial/tally* obtains the total virial heat flux +(strictly speaking, heat flow) into the first group due to interaction +with the second group, and is defined as: + +.. math:: + + Q = \frac{1}{2} \sum_{i \in \text{group 1}} \sum_{j \in \text{group 2}} \mathbf{F}_{ij} \cdot \left(\mathbf{v}_i + \mathbf{v}_j \right) + +Although, the *heat/flux/virial/tally* compute +does not include the convective term, +it can be used to obtain the total heat flux over control surfaces, +when there are no particles crossing over, +such as is often in solid-solid and solid-liquid interfaces. +This would be identical to the method of planes method. +Note that the *heat/flux/virial/tally* compute is distinctly different +from the *heat/flux* and *heat/flux/tally* computes, +that are essentially volume averaging methods. +The following example demonstrates the difference: + +.. code-block:: LAMMPS + + # System with only pairwise interactions. + # Non-periodic boundaries in the x direction. + # Has LeftLiquid and RightWall groups along x direction. + + # Heat flux over the solid-liquid interface + compute hflow_hfvt LeftLiquid heat/flux/virial/tally RightWall + variable hflux_hfvt equal c_hflow_hfvt/(ly*lz) + + # x component of approximate heat flux vector inside the liquid region, + # two approaches. + # + compute myKE all ke/atom + compute myPE all pe/atom + compute myStress all stress/atom NULL virial + compute hflow_hf LeftLiquid heat/flux myKE myPE myStress + variable hflux_hf equal c_hflow_hf[1]/${volLiq} + # + compute hflow_hft LeftLiquid heat/flux/tally all + variable hflux_hft equal c_hflow_hft[1]/${volLiq} + + # Pressure over the solid-liquid interface, three approaches. + # + compute force_gg RightWall group/group LeftLiquid + variable press_gg equal c_force_gg[1]/(ly*lz) + # + compute force_ft RightWall force/tally LeftLiquid + compute rforce_ft RightWall reduce sum c_force_ft[1] + variable press_ft equal c_rforce_ft/(ly*lz) + # + compute rforce_hfvt all reduce sum c_hflow_hfvt[1] + variable press_hfvt equal -c_rforce_hfvt/(ly*lz) + +---------- + The pairwise contributions are computing via a callback that the compute registers with the non-bonded pairwise force computation. This limits the use to systems that have no bonds, no Kspace, and no @@ -83,7 +161,17 @@ magnitude) and a per atom 3-element vector (force contribution from each atom). Compute *stress/tally* calculates a global scalar (average of the diagonal elements of the stress tensor) and a per atom vector (the 6 elements of stress tensor contributions from the -individual atom). +individual atom). As in :doc:`compute heat/flux `, +compute *heat/flux/tally* calculates a global vector of length 6, +where the first 3 components are the :math:`x`, :math:`y`, :math:`z` +components of the full heat flux vector, +and the next 3 components are the corresponding components +of just the convective portion of the flux, i.e. the +first term in the equation for :math:`\mathbf{Q}`. +Compute *heat/flux/virial/tally* calculates a global scalar (heat flow) +and a per atom 3-element vector +(contribution to the force acting over atoms in the first group +from individual atoms in both groups). Both the scalar and vector values calculated by this compute are "extensive". From 2a4b60d5970c1b11339d0c70e70544efdbf435f2 Mon Sep 17 00:00:00 2001 From: Donatas Surblys Date: Tue, 20 Jul 2021 18:03:11 +0900 Subject: [PATCH 4/6] in the output info section: flux -> flow --- doc/src/compute_tally.rst | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/doc/src/compute_tally.rst b/doc/src/compute_tally.rst index d2e4569820..81756dda7c 100644 --- a/doc/src/compute_tally.rst +++ b/doc/src/compute_tally.rst @@ -164,9 +164,9 @@ vector (the 6 elements of stress tensor contributions from the individual atom). As in :doc:`compute heat/flux `, compute *heat/flux/tally* calculates a global vector of length 6, where the first 3 components are the :math:`x`, :math:`y`, :math:`z` -components of the full heat flux vector, +components of the full heat flow vector, and the next 3 components are the corresponding components -of just the convective portion of the flux, i.e. the +of just the convective portion of the flow, i.e. the first term in the equation for :math:`\mathbf{Q}`. Compute *heat/flux/virial/tally* calculates a global scalar (heat flow) and a per atom 3-element vector From b6e749f7fc9192baa916976241ddb90531511a3f Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Tue, 20 Jul 2021 14:27:29 -0400 Subject: [PATCH 5/6] reformat with clang-format --- src/TALLY/compute_heat_flux_virial_tally.cpp | 99 +++++++++----------- 1 file changed, 46 insertions(+), 53 deletions(-) diff --git a/src/TALLY/compute_heat_flux_virial_tally.cpp b/src/TALLY/compute_heat_flux_virial_tally.cpp index 0ac3607da6..1a594c1b36 100644 --- a/src/TALLY/compute_heat_flux_virial_tally.cpp +++ b/src/TALLY/compute_heat_flux_virial_tally.cpp @@ -1,4 +1,3 @@ -// clang-format off /* ---------------------------------------------------------------------- LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator https://www.lammps.org/, Sandia National Laboratories @@ -14,28 +13,29 @@ #include "compute_heat_flux_virial_tally.h" -#include #include "atom.h" -#include "group.h" -#include "pair.h" -#include "update.h" -#include "memory.h" +#include "comm.h" #include "error.h" #include "force.h" -#include "comm.h" +#include "group.h" +#include "memory.h" +#include "pair.h" +#include "update.h" + +#include using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ ComputeHeatFluxVirialTally::ComputeHeatFluxVirialTally(LAMMPS *lmp, int narg, char **arg) : - Compute(lmp, narg, arg) + Compute(lmp, narg, arg) { - if (narg < 4) error->all(FLERR,"Illegal compute heat/flux/virial/tally command"); + if (narg < 4) error->all(FLERR, "Illegal compute heat/flux/virial/tally command"); igroup2 = group->find(arg[3]); if (igroup2 == -1) - error->all(FLERR,"Could not find compute heat/flux/virial/tally second group ID"); + error->all(FLERR, "Could not find compute heat/flux/virial/tally second group ID"); groupbit2 = group->bitmask[igroup2]; scalar_flag = 1; @@ -46,7 +46,7 @@ ComputeHeatFluxVirialTally::ComputeHeatFluxVirialTally(LAMMPS *lmp, int narg, ch comm_reverse = size_peratom_cols = 3; extscalar = 1; - peflag = 1; // we need Pair::ev_tally() to be run + peflag = 1; // we need Pair::ev_tally() to be run did_setup = invoked_peratom = invoked_scalar = -1; nmax = -1; @@ -66,29 +66,30 @@ ComputeHeatFluxVirialTally::~ComputeHeatFluxVirialTally() void ComputeHeatFluxVirialTally::init() { if (force->pair == nullptr) - error->all(FLERR,"Trying to use compute heat/flux/virial/tally without pair style"); + error->all(FLERR, "Trying to use compute heat/flux/virial/tally without pair style"); else force->pair->add_tally_callback(this); if (comm->me == 0) { if (force->pair->single_enable == 0 || force->pair->manybody_flag) - error->warning(FLERR,"Compute heat/flux/virial/tally used with incompatible pair style"); + error->warning(FLERR, "Compute heat/flux/virial/tally used with incompatible pair style"); - if (force->bond || force->angle || force->dihedral - || force->improper || force->kspace) - error->warning(FLERR,"Compute heat/flux/virial/tally only called from pair style"); + if (force->bond || force->angle || force->dihedral || force->improper || force->kspace) + error->warning(FLERR, "Compute heat/flux/virial/tally only called from pair style"); } // error out if any atoms are in both groups - for (int i = 0; i < atom->nlocal; i++) - { + for (int i = 0; i < atom->nlocal; i++) { if ((atom->mask[i] & groupbit) && (atom->mask[i] & groupbit2)) { if (atom->tag_enable) { - error->all(FLERR,"Atom {} belonging to both groups is not allowed" - " with compute heat/flux/virial/tally", atom->tag[i]); + error->all(FLERR, + "Atom {} belonging to both groups is not allowed" + " with compute heat/flux/virial/tally", + atom->tag[i]); } else { - error->all(FLERR,"Atom belonging to both groups is not allowed" - " with compute heat/flux/virial/tally"); + error->all(FLERR, + "Atom belonging to both groups is not allowed" + " with compute heat/flux/virial/tally"); } } } @@ -112,34 +113,32 @@ void ComputeHeatFluxVirialTally::pair_setup_callback(int, int) if (atom->nmax > nmax) { memory->destroy(fatom); nmax = atom->nmax; - memory->create(fatom,nmax,size_peratom_cols,"heat/flux/virial/tally:fatom"); + memory->create(fatom, nmax, size_peratom_cols, "heat/flux/virial/tally:fatom"); array_atom = fatom; } // clear storage - for (int i=0; i < ntotal; ++i) - for (int j=0; j < size_peratom_cols; ++j) - fatom[i][j] = 0.0; + for (int i = 0; i < ntotal; ++i) + for (int j = 0; j < size_peratom_cols; ++j) fatom[i][j] = 0.0; did_setup = update->ntimestep; } /* ---------------------------------------------------------------------- */ -void ComputeHeatFluxVirialTally::pair_tally_callback(int i, int j, int nlocal, int newton, - double, double, double fpair, - double dx, double dy, double dz) +void ComputeHeatFluxVirialTally::pair_tally_callback(int i, int j, int nlocal, int newton, double, + double, double fpair, double dx, double dy, + double dz) { - const int * const mask = atom->mask; + const int *const mask = atom->mask; const bool ingroup1 = (mask[i] & groupbit); - if ( (ingroup1 && (mask[j] & groupbit2)) - || ((mask[i] & groupbit2) && (mask[j] & groupbit)) ) { + if ((ingroup1 && (mask[j] & groupbit2)) || ((mask[i] & groupbit2) && (mask[j] & groupbit))) { // signs set to calculate heat flux from group2 into group1 - const double fx = (ingroup1?0.5:-0.5)*dx*fpair; - const double fy = (ingroup1?0.5:-0.5)*dy*fpair; - const double fz = (ingroup1?0.5:-0.5)*dz*fpair; + const double fx = (ingroup1 ? 0.5 : -0.5) * dx * fpair; + const double fy = (ingroup1 ? 0.5 : -0.5) * dy * fpair; + const double fz = (ingroup1 ? 0.5 : -0.5) * dz * fpair; if (newton || i < nlocal) { fatom[i][0] += fx; @@ -158,7 +157,7 @@ void ComputeHeatFluxVirialTally::pair_tally_callback(int i, int j, int nlocal, i int ComputeHeatFluxVirialTally::pack_reverse_comm(int n, int first, double *buf) { - int i,m,last; + int i, m, last; m = 0; last = first + n; @@ -174,7 +173,7 @@ int ComputeHeatFluxVirialTally::pack_reverse_comm(int n, int first, double *buf) void ComputeHeatFluxVirialTally::unpack_reverse_comm(int n, int *list, double *buf) { - int i,j,m; + int i, j, m; m = 0; for (i = 0; i < n; i++) { @@ -191,22 +190,19 @@ double ComputeHeatFluxVirialTally::compute_scalar() { // need to collect contributions from ghost atoms // because we need to use local velocities to compute heat flux - if (invoked_peratom != update->ntimestep) - compute_peratom(); + if (invoked_peratom != update->ntimestep) compute_peratom(); invoked_scalar = update->ntimestep; - if ((did_setup != invoked_scalar) - || (update->eflag_global != invoked_scalar)) - error->all(FLERR,"Energy was not tallied on needed timestep"); + if ((did_setup != invoked_scalar) || (update->eflag_global != invoked_scalar)) + error->all(FLERR, "Energy was not tallied on needed timestep"); // sum heat flux across procs double hflux = 0.0; for (int i = 0; i < atom->nlocal; i++) - hflux += fatom[i][0]*atom->v[i][0] - + fatom[i][1]*atom->v[i][1] - + fatom[i][2]*atom->v[i][2]; + hflux += + fatom[i][0] * atom->v[i][0] + fatom[i][1] * atom->v[i][1] + fatom[i][2] * atom->v[i][2]; - MPI_Allreduce(&hflux,&scalar,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&hflux, &scalar, 1, MPI_DOUBLE, MPI_SUM, world); return scalar; } @@ -216,9 +212,8 @@ double ComputeHeatFluxVirialTally::compute_scalar() void ComputeHeatFluxVirialTally::compute_peratom() { invoked_peratom = update->ntimestep; - if ((did_setup != invoked_peratom) - || (update->eflag_global != invoked_peratom)) - error->all(FLERR,"Energy was not tallied on needed timestep"); + if ((did_setup != invoked_peratom) || (update->eflag_global != invoked_peratom)) + error->all(FLERR, "Energy was not tallied on needed timestep"); // collect contributions from ghost atoms @@ -228,8 +223,7 @@ void ComputeHeatFluxVirialTally::compute_peratom() // clear out ghost atom data after it has been collected to local atoms const int nall = atom->nlocal + atom->nghost; for (int i = atom->nlocal; i < nall; ++i) - for (int j = 0; j < size_peratom_cols; ++j) - fatom[i][j] = 0.0; + for (int j = 0; j < size_peratom_cols; ++j) fatom[i][j] = 0.0; } } @@ -239,7 +233,6 @@ void ComputeHeatFluxVirialTally::compute_peratom() double ComputeHeatFluxVirialTally::memory_usage() { - double bytes = (nmax < 0) ? 0 : nmax*size_peratom_cols * sizeof(double); + double bytes = (nmax < 0) ? 0 : nmax * size_peratom_cols * sizeof(double); return bytes; } - From eeea5660937edef190509d45eb30822dc2e9eb9a Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Tue, 20 Jul 2021 14:32:30 -0400 Subject: [PATCH 6/6] correct typo --- doc/src/compute.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/src/compute.rst b/doc/src/compute.rst index 93916d701d..7a89c6cc87 100644 --- a/doc/src/compute.rst +++ b/doc/src/compute.rst @@ -241,7 +241,7 @@ The individual style names on the :doc:`Commands compute ` doc * :doc:`pe ` - potential energy * :doc:`pe/atom ` - potential energy for each atom * :doc:`mesont ` - Nanotube bending,stretching, and intertube energies -* :doc:`pe/mol/tally ` - potential energy between two groups of atoms seperated into intermolecular and intramolecular components via the tally callback mechanism +* :doc:`pe/mol/tally ` - potential energy between two groups of atoms separated into intermolecular and intramolecular components via the tally callback mechanism * :doc:`pe/tally ` - potential energy between two groups of atoms via the tally callback mechanism * :doc:`plasticity/atom ` - Peridynamic plasticity for each atom * :doc:`pressure ` - total pressure and pressure tensor