|
|
|
|
@ -33,32 +33,31 @@ Examples
|
|
|
|
|
Description
|
|
|
|
|
"""""""""""
|
|
|
|
|
|
|
|
|
|
Define a computation that computes per-atom stress
|
|
|
|
|
tensor for each atom in a group. In case of compute *stress/atom*,
|
|
|
|
|
the tensor for each atom is symmetric with 6
|
|
|
|
|
components and is stored as a 6-element vector in the following order:
|
|
|
|
|
:math:`xx`, :math:`yy`, :math:`zz`, :math:`xy`, :math:`xz`, :math:`yz`.
|
|
|
|
|
In case of compute *centroid/stress/atom*,
|
|
|
|
|
the tensor for each atom is asymmetric with 9 components
|
|
|
|
|
and is stored as a 9-element vector in the following order:
|
|
|
|
|
:math:`xx`, :math:`yy`, :math:`zz`, :math:`xy`, :math:`xz`, :math:`yz`,
|
|
|
|
|
:math:`yx`, :math:`zx`, :math:`zy`.
|
|
|
|
|
See the :doc:`compute pressure <compute_pressure>` command if you want the stress tensor
|
|
|
|
|
Define a computation that computes per-atom stress tensor for each
|
|
|
|
|
atom in a group. In case of compute *stress/atom*, the tensor for
|
|
|
|
|
each atom is symmetric with 6 components and is stored as a 6-element
|
|
|
|
|
vector in the following order: :math:`xx`, :math:`yy`, :math:`zz`,
|
|
|
|
|
:math:`xy`, :math:`xz`, :math:`yz`. In case of compute
|
|
|
|
|
*centroid/stress/atom*, the tensor for each atom is asymmetric with 9
|
|
|
|
|
components and is stored as a 9-element vector in the following order:
|
|
|
|
|
:math:`xx`, :math:`yy`, :math:`zz`, :math:`xy`, :math:`xz`,
|
|
|
|
|
:math:`yz`, :math:`yx`, :math:`zx`, :math:`zy`. See the :doc:`compute
|
|
|
|
|
pressure <compute_pressure>` command if you want the stress tensor
|
|
|
|
|
(pressure) of the entire system.
|
|
|
|
|
|
|
|
|
|
The stress tensor for atom :math:`I` is given by the following formula,
|
|
|
|
|
where :math:`a` and :math:`b` take on values :math:`x`, :math:`y`, :math:`z`
|
|
|
|
|
to generate the components of the tensor:
|
|
|
|
|
The stress tensor for atom :math:`I` is given by the following
|
|
|
|
|
formula, where :math:`a` and :math:`b` take on values :math:`x`,
|
|
|
|
|
:math:`y`, :math:`z` to generate the components of the tensor:
|
|
|
|
|
|
|
|
|
|
.. math::
|
|
|
|
|
|
|
|
|
|
S_{ab} = - m v_a v_b - W_{ab}
|
|
|
|
|
|
|
|
|
|
The first term is a kinetic energy contribution for atom :math:`I`. See
|
|
|
|
|
details below on how the specified *temp-ID* can affect the velocities
|
|
|
|
|
used in this calculation. The second term is the virial
|
|
|
|
|
contribution due to intra and intermolecular interactions,
|
|
|
|
|
where the exact computation details are determined by the compute style.
|
|
|
|
|
The first term is a kinetic energy contribution for atom :math:`I`.
|
|
|
|
|
See details below on how the specified *temp-ID* can affect the
|
|
|
|
|
velocities used in this calculation. The second term is the virial
|
|
|
|
|
contribution due to intra and intermolecular interactions, where the
|
|
|
|
|
exact computation details are determined by the compute style.
|
|
|
|
|
|
|
|
|
|
In case of compute *stress/atom*, the virial contribution is:
|
|
|
|
|
|
|
|
|
|
@ -68,29 +67,26 @@ In case of compute *stress/atom*, the virial contribution is:
|
|
|
|
|
& + \frac{1}{3} \sum_{n = 1}^{N_a} (r_{1_a} F_{1_b} + r_{2_a} F_{2_b} + r_{3_a} F_{3_b}) + \frac{1}{4} \sum_{n = 1}^{N_d} (r_{1_a} F_{1_b} + r_{2_a} F_{2_b} + r_{3_a} F_{3_b} + r_{4_a} F_{4_b}) \\
|
|
|
|
|
& + \frac{1}{4} \sum_{n = 1}^{N_i} (r_{1_a} F_{1_b} + r_{2_a} F_{2_b} + r_{3_a} F_{3_b} + r_{4_a} F_{4_b}) + {\rm Kspace}(r_{i_a},F_{i_b}) + \sum_{n = 1}^{N_f} r_{i_a} F_{i_b}
|
|
|
|
|
|
|
|
|
|
The first term is a pairwise energy
|
|
|
|
|
contribution where :math:`n` loops over the :math:`N_p`
|
|
|
|
|
neighbors of atom :math:`I`, :math:`\mathbf{r}_1` and :math:`\mathbf{r}_2`
|
|
|
|
|
are the positions of the 2 atoms in the pairwise interaction,
|
|
|
|
|
and :math:`\mathbf{F}_1` and :math:`\mathbf{F}_2` are the forces
|
|
|
|
|
on the 2 atoms resulting from the pairwise interaction.
|
|
|
|
|
The second term is a bond contribution of
|
|
|
|
|
similar form for the :math:`N_b` bonds which atom :math:`I` is part of.
|
|
|
|
|
There are similar terms for the :math:`N_a` angle, :math:`N_d` dihedral,
|
|
|
|
|
and :math:`N_i` improper interactions atom :math:`I` is part of.
|
|
|
|
|
There is also a term for the KSpace
|
|
|
|
|
contribution from long-range Coulombic interactions, if defined.
|
|
|
|
|
Finally, there is a term for the :math:`N_f` :doc:`fixes <fix>` that apply
|
|
|
|
|
internal constraint forces to atom :math:`I`. Currently, only the
|
|
|
|
|
:doc:`fix shake <fix_shake>` and :doc:`fix rigid <fix_rigid>` commands
|
|
|
|
|
contribute to this term.
|
|
|
|
|
As the coefficients in the formula imply, a virial contribution
|
|
|
|
|
produced by a small set of atoms (e.g. 4 atoms in a dihedral or 3
|
|
|
|
|
atoms in a Tersoff 3-body interaction) is assigned in equal portions
|
|
|
|
|
to each atom in the set. E.g. 1/4 of the dihedral virial to each of
|
|
|
|
|
the 4 atoms, or 1/3 of the fix virial due to SHAKE constraints applied
|
|
|
|
|
to atoms in a water molecule via the :doc:`fix shake <fix_shake>`
|
|
|
|
|
command.
|
|
|
|
|
The first term is a pairwise energy contribution where :math:`n` loops
|
|
|
|
|
over the :math:`N_p` neighbors of atom :math:`I`, :math:`\mathbf{r}_1`
|
|
|
|
|
and :math:`\mathbf{r}_2` are the positions of the 2 atoms in the
|
|
|
|
|
pairwise interaction, and :math:`\mathbf{F}_1` and
|
|
|
|
|
:math:`\mathbf{F}_2` are the forces on the 2 atoms resulting from the
|
|
|
|
|
pairwise interaction. The second term is a bond contribution of
|
|
|
|
|
similar form for the :math:`N_b` bonds which atom :math:`I` is part
|
|
|
|
|
of. There are similar terms for the :math:`N_a` angle, :math:`N_d`
|
|
|
|
|
dihedral, and :math:`N_i` improper interactions atom :math:`I` is part
|
|
|
|
|
of. There is also a term for the KSpace contribution from long-range
|
|
|
|
|
Coulombic interactions, if defined. Finally, there is a term for the
|
|
|
|
|
:math:`N_f` :doc:`fixes <fix>` that apply internal constraint forces
|
|
|
|
|
to atom :math:`I`. Currently, only the :doc:`fix shake <fix_shake>`
|
|
|
|
|
and :doc:`fix rigid <fix_rigid>` commands contribute to this term. As
|
|
|
|
|
the coefficients in the formula imply, a virial contribution produced
|
|
|
|
|
by a small set of atoms (e.g. 4 atoms in a dihedral or 3 atoms in a
|
|
|
|
|
Tersoff 3-body interaction) is assigned in equal portions to each atom
|
|
|
|
|
in the set. E.g. 1/4 of the dihedral virial to each of the 4 atoms,
|
|
|
|
|
or 1/3 of the fix virial due to SHAKE constraints applied to atoms in
|
|
|
|
|
a water molecule via the :doc:`fix shake <fix_shake>` command.
|
|
|
|
|
|
|
|
|
|
In case of compute *centroid/stress/atom*, the virial contribution is:
|
|
|
|
|
|
|
|
|
|
@ -99,71 +95,69 @@ In case of compute *centroid/stress/atom*, the virial contribution is:
|
|
|
|
|
W_{ab} & = \sum_{n = 1}^{N_p} r_{I0_a} F_{I_b} + \sum_{n = 1}^{N_b} r_{I0_a} F_{I_b} + \sum_{n = 1}^{N_a} r_{I0_a} F_{I_b} + \sum_{n = 1}^{N_d} r_{I0_a} F_{I_b} + \sum_{n = 1}^{N_i} r_{I0_a} F_{I_b} \\
|
|
|
|
|
& + {\rm Kspace}(r_{i_a},F_{i_b}) + \sum_{n = 1}^{N_f} r_{i_a} F_{i_b}
|
|
|
|
|
|
|
|
|
|
As with compute *stress/atom*, the first, second, third, fourth and fifth terms
|
|
|
|
|
are pairwise, bond, angle, dihedral and improper contributions,
|
|
|
|
|
but instead of assigning the virial contribution equally to each atom,
|
|
|
|
|
only the force :math:`\mathbf{F}_I` acting on atom :math:`I`
|
|
|
|
|
due to the interaction and the relative
|
|
|
|
|
position :math:`\mathbf{r}_{I0}` of the atom :math:`I` to the geometric center
|
|
|
|
|
of the interacting atoms, i.e. centroid, is used.
|
|
|
|
|
As the geometric center is different
|
|
|
|
|
for each interaction, the :math:`\mathbf{r}_{I0}` also differs.
|
|
|
|
|
The sixth and seventh terms, Kspace and :doc:`fix <fix>` contribution
|
|
|
|
|
respectively, are computed identical to compute *stress/atom*.
|
|
|
|
|
Although the total system virial is the same as compute *stress/atom*,
|
|
|
|
|
compute *centroid/stress/atom* is know to result in more consistent
|
|
|
|
|
heat flux values for angle, dihedrals and improper contributions
|
|
|
|
|
when computed via :doc:`compute heat/flux <compute_heat_flux>`.
|
|
|
|
|
As with compute *stress/atom*, the first, second, third, fourth and
|
|
|
|
|
fifth terms are pairwise, bond, angle, dihedral and improper
|
|
|
|
|
contributions, but instead of assigning the virial contribution
|
|
|
|
|
equally to each atom, only the force :math:`\mathbf{F}_I` acting on
|
|
|
|
|
atom :math:`I` due to the interaction and the relative position
|
|
|
|
|
:math:`\mathbf{r}_{I0}` of the atom :math:`I` to the geometric center
|
|
|
|
|
of the interacting atoms, i.e. centroid, is used. As the geometric
|
|
|
|
|
center is different for each interaction, the :math:`\mathbf{r}_{I0}`
|
|
|
|
|
also differs. The sixth and seventh terms, Kspace and :doc:`fix
|
|
|
|
|
<fix>` contribution respectively, are computed identical to compute
|
|
|
|
|
*stress/atom*. Although the total system virial is the same as
|
|
|
|
|
compute *stress/atom*, compute *centroid/stress/atom* is know to
|
|
|
|
|
result in more consistent heat flux values for angle, dihedrals and
|
|
|
|
|
improper contributions when computed via :doc:`compute heat/flux
|
|
|
|
|
<compute_heat_flux>`.
|
|
|
|
|
|
|
|
|
|
If no extra keywords are listed, the kinetic contribution
|
|
|
|
|
all of the virial contribution terms are
|
|
|
|
|
included in the per-atom stress tensor. If any extra keywords are
|
|
|
|
|
listed, only those terms are summed to compute the tensor. The
|
|
|
|
|
*virial* keyword means include all terms except the kinetic energy
|
|
|
|
|
*ke*\ .
|
|
|
|
|
If no extra keywords are listed, the kinetic contribution all of the
|
|
|
|
|
virial contribution terms are included in the per-atom stress tensor.
|
|
|
|
|
If any extra keywords are listed, only those terms are summed to
|
|
|
|
|
compute the tensor. The *virial* keyword means include all terms
|
|
|
|
|
except the kinetic energy *ke*\ .
|
|
|
|
|
|
|
|
|
|
Note that the stress for each atom is due to its interaction with all
|
|
|
|
|
other atoms in the simulation, not just with other atoms in the group.
|
|
|
|
|
|
|
|
|
|
Details of how compute *stress/atom* obtains the virial for individual atoms for
|
|
|
|
|
either pairwise or many-body potentials, and including the effects of
|
|
|
|
|
periodic boundary conditions is discussed in :ref:`(Thompson) <Thompson2>`.
|
|
|
|
|
The basic idea for many-body potentials is to treat each component of
|
|
|
|
|
the force computation between a small cluster of atoms in the same
|
|
|
|
|
manner as in the formula above for bond, angle, dihedral, etc
|
|
|
|
|
interactions. Namely the quantity :math:`\mathbf{r} \cdot \mathbf{F}`
|
|
|
|
|
is summed over the atoms in
|
|
|
|
|
the interaction, with the :math:`r` vectors unwrapped by periodic boundaries
|
|
|
|
|
so that the cluster of atoms is close together. The total
|
|
|
|
|
Details of how compute *stress/atom* obtains the virial for individual
|
|
|
|
|
atoms for either pairwise or many-body potentials, and including the
|
|
|
|
|
effects of periodic boundary conditions is discussed in
|
|
|
|
|
:ref:`(Thompson) <Thompson2>`. The basic idea for many-body
|
|
|
|
|
potentials is to treat each component of the force computation between
|
|
|
|
|
a small cluster of atoms in the same manner as in the formula above
|
|
|
|
|
for bond, angle, dihedral, etc interactions. Namely the quantity
|
|
|
|
|
:math:`\mathbf{r} \cdot \mathbf{F}` is summed over the atoms in the
|
|
|
|
|
interaction, with the :math:`r` vectors unwrapped by periodic
|
|
|
|
|
boundaries so that the cluster of atoms is close together. The total
|
|
|
|
|
contribution for the cluster interaction is divided evenly among those
|
|
|
|
|
atoms. Details of how compute *centroid/stress/atom* obtains
|
|
|
|
|
the virial for individual atoms
|
|
|
|
|
is given in :ref:`(Surblys) <Surblys1>`,
|
|
|
|
|
where the idea is that the virial of the atom :math:`I`
|
|
|
|
|
is the result of only the force :math:`\mathbf{F}_I` on the atom due
|
|
|
|
|
to the interaction
|
|
|
|
|
and its positional vector :math:`\mathbf{r}_{I0}`,
|
|
|
|
|
relative to the geometric center of the
|
|
|
|
|
interacting atoms, regardless of the number of participating atoms.
|
|
|
|
|
The periodic boundary treatment is identical to
|
|
|
|
|
atoms.
|
|
|
|
|
|
|
|
|
|
Details of how compute *centroid/stress/atom* obtains the virial for
|
|
|
|
|
individual atoms is given in :ref:`(Surblys) <Surblys1>`, where the
|
|
|
|
|
idea is that the virial of the atom :math:`I` is the result of only
|
|
|
|
|
the force :math:`\mathbf{F}_I` on the atom due to the interaction and
|
|
|
|
|
its positional vector :math:`\mathbf{r}_{I0}`, relative to the
|
|
|
|
|
geometric center of the interacting atoms, regardless of the number of
|
|
|
|
|
participating atoms. The periodic boundary treatment is identical to
|
|
|
|
|
that of compute *stress/atom*, and both of them reduce to identical
|
|
|
|
|
expressions for two-body interactions,
|
|
|
|
|
i.e. computed values for contributions from bonds and two-body pair styles,
|
|
|
|
|
such as :doc:`Lennard-Jones <pair_lj>`, will be the same,
|
|
|
|
|
while contributions from angles, dihedrals and impropers will be different.
|
|
|
|
|
expressions for two-body interactions, i.e. computed values for
|
|
|
|
|
contributions from bonds and two-body pair styles, such as
|
|
|
|
|
:doc:`Lennard-Jones <pair_lj>`, will be the same, while contributions
|
|
|
|
|
from angles, dihedrals and impropers will be different.
|
|
|
|
|
|
|
|
|
|
The :doc:`dihedral_style charmm <dihedral_charmm>` style calculates
|
|
|
|
|
pairwise interactions between 1-4 atoms. The virial contribution of
|
|
|
|
|
these terms is included in the pair virial, not the dihedral virial.
|
|
|
|
|
|
|
|
|
|
The KSpace contribution is calculated using the method in
|
|
|
|
|
:ref:`(Heyes) <Heyes2>` for the Ewald method and by the methodology described
|
|
|
|
|
in :ref:`(Sirk) <Sirk1>` for PPPM. The choice of KSpace solver is specified
|
|
|
|
|
by the :doc:`kspace_style pppm <kspace_style>` command. Note that for
|
|
|
|
|
PPPM, the calculation requires 6 extra FFTs each timestep that
|
|
|
|
|
per-atom stress is calculated. Thus it can significantly increase the
|
|
|
|
|
cost of the PPPM calculation if it is needed on a large fraction of
|
|
|
|
|
the simulation timesteps.
|
|
|
|
|
:ref:`(Heyes) <Heyes2>` for the Ewald method and by the methodology
|
|
|
|
|
described in :ref:`(Sirk) <Sirk1>` for PPPM. The choice of KSpace
|
|
|
|
|
solver is specified by the :doc:`kspace_style pppm <kspace_style>`
|
|
|
|
|
command. Note that for PPPM, the calculation requires 6 extra FFTs
|
|
|
|
|
each timestep that per-atom stress is calculated. Thus it can
|
|
|
|
|
significantly increase the cost of the PPPM calculation if it is
|
|
|
|
|
needed on a large fraction of the simulation timesteps.
|
|
|
|
|
|
|
|
|
|
The *temp-ID* argument can be used to affect the per-atom velocities
|
|
|
|
|
used in the kinetic energy contribution to the total stress. If the
|
|
|
|
|
@ -189,10 +183,10 @@ See the :doc:`compute voronoi/atom <compute_voronoi_atom>` command for
|
|
|
|
|
one possible way to estimate a per-atom volume.
|
|
|
|
|
|
|
|
|
|
Thus, if the diagonal components of the per-atom stress tensor are
|
|
|
|
|
summed for all atoms in the system and the sum is divided by :math:`dV`, where
|
|
|
|
|
:math:`d` = dimension and :math:`V` is the volume of the system,
|
|
|
|
|
the result should be :math:`-P`, where :math:`P`
|
|
|
|
|
is the total pressure of the system.
|
|
|
|
|
summed for all atoms in the system and the sum is divided by
|
|
|
|
|
:math:`dV`, where :math:`d` = dimension and :math:`V` is the volume of
|
|
|
|
|
the system, the result should be :math:`-P`, where :math:`P` is the
|
|
|
|
|
total pressure of the system.
|
|
|
|
|
|
|
|
|
|
These lines in an input script for a 3d system should yield that
|
|
|
|
|
result. I.e. the last 2 columns of thermo output will be the same:
|
|
|
|
|
@ -207,33 +201,43 @@ result. I.e. the last 2 columns of thermo output will be the same:
|
|
|
|
|
.. note::
|
|
|
|
|
|
|
|
|
|
The per-atom stress does not include any Lennard-Jones tail
|
|
|
|
|
corrections to the pressure added by the :doc:`pair_modify tail yes <pair_modify>` command, since those are contributions to the
|
|
|
|
|
global system pressure.
|
|
|
|
|
corrections to the pressure added by the :doc:`pair_modify tail yes
|
|
|
|
|
<pair_modify>` command, since those are contributions to the global
|
|
|
|
|
system pressure.
|
|
|
|
|
|
|
|
|
|
Output info
|
|
|
|
|
"""""""""""
|
|
|
|
|
|
|
|
|
|
This compute *stress/atom* calculates a per-atom array with 6 columns, which can be
|
|
|
|
|
accessed by indices 1-6 by any command that uses per-atom values from
|
|
|
|
|
a compute as input.
|
|
|
|
|
The compute *centroid/stress/atom* produces a per-atom array with 9 columns,
|
|
|
|
|
but otherwise can be used in an identical manner to compute *stress/atom*.
|
|
|
|
|
See the :doc:`Howto output <Howto_output>` doc page
|
|
|
|
|
for an overview of LAMMPS output options.
|
|
|
|
|
Compute *stress/atom* calculates a per-atom array with 6 columns,
|
|
|
|
|
which can be accessed by indices 1-6 by any command that uses per-atom
|
|
|
|
|
values from a compute as input. Compute *centroid/stress/atom*
|
|
|
|
|
produces a per-atom array with 9 columns, but otherwise can be used in
|
|
|
|
|
an identical manner to compute *stress/atom*. See the :doc:`Howto
|
|
|
|
|
output <Howto_output>` doc page for an overview of LAMMPS output
|
|
|
|
|
options.
|
|
|
|
|
|
|
|
|
|
The per-atom array values will be in pressure\*volume
|
|
|
|
|
:doc:`units <units>` as discussed above.
|
|
|
|
|
The per-atom array values will be in pressure\*volume :doc:`units
|
|
|
|
|
<units>` as discussed above.
|
|
|
|
|
|
|
|
|
|
Restrictions
|
|
|
|
|
""""""""""""
|
|
|
|
|
|
|
|
|
|
Currently (Spring 2020), compute *centroid/stress/atom* does not support
|
|
|
|
|
pair styles with many-body interactions, such as :doc:`Tersoff
|
|
|
|
|
<pair_tersoff>`, or pair styles with long-range Coulomb interactions.
|
|
|
|
|
LAMMPS will generate an error in such cases. In principal, equivalent
|
|
|
|
|
formulation to that of angle, dihedral and improper contributions in the
|
|
|
|
|
virial :math:`W_{ab}` formula can also be applied to the many-body pair
|
|
|
|
|
styles, and is planned in the future.
|
|
|
|
|
Currently, compute *centroid/stress/atom* does not support pair styles
|
|
|
|
|
with many-body interactions (:doc:`EAM <pair_eam>` is an exception,
|
|
|
|
|
since its computations are performed pairwise), nor granular pair
|
|
|
|
|
styles with pairwise forces which are not aligned with the vector
|
|
|
|
|
between the pair of particles. All bond styles are supported. All
|
|
|
|
|
angle, dihedral, improper styles are supported with the exception of
|
|
|
|
|
USER-INTEL and KOKKOS variants of specific styles. It also does not
|
|
|
|
|
support models with long-range Coulombic or dispersion forces,
|
|
|
|
|
i.e. the kspace_style command in LAMMPS. It also does not support the
|
|
|
|
|
following fixes which add rigid-body constraints: :doc:`fix shake
|
|
|
|
|
<fix_shake>`, :doc:`fix rattle <fix_rattle>`, :doc:`fix rigid
|
|
|
|
|
<fix_rigid>`, :doc:`fix rigid/small <fix_rigid>`.
|
|
|
|
|
|
|
|
|
|
LAMMPS will generate an error if one of these options is included in
|
|
|
|
|
your model. Extension of centroid stress calculations to these force
|
|
|
|
|
and fix styles is planned for the futre.
|
|
|
|
|
|
|
|
|
|
Related commands
|
|
|
|
|
""""""""""""""""
|
|
|
|
|
|
Stan Moore