Syntax:
compute ID group-ID stress/atom keyword ...
Examples:
compute 1 mobile stress/atom
compute 1 all stress/atom pair bond
Description:
Define a computation that computes the symmetric per-atom stress tensor for each atom in a group. The tensor for each atom has 6 components: xx, yy, zz, xy, xz, yz. See the compute pressure command if you want the stress tensor (pressure) of the entire system. The 6 components can be accessed by indices 1-6 by any command that uses per-atom computes, e.g. the dump custom command or fix ave/spatial command or fix ave/atom command. See this section for an overview.
The stress tensor for atom I is given by the following formula, where a and b take on values x,y,z to generate the 6 components of the symmetric tensor:
The first term is a kinetic energy contribution for atom I. The second term is a pairwise energy contribution where n loops over the Np neighbors of atom I, r1 and r2 are the positions of the 2 atoms in the pairwise interaction, and F1 and F2 are the forces on the 2 atoms resulting from the pairwise interaction. The third term is a bond contribution of similar form for the Nb bonds which atom I is part of. The remaining terms are for the Na angle, Nd dihedral, and Ni improper interactions atom I is part of.
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.
If no extra keywords are listed, all of the terms in this formula are included in the per-atom stress tensor. If any extra keywords are listed, only those terms are summed to compute the tensor.
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.
The dihedral_style 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.
Note that as defined in the formula, per-atom stress is the negative of the per-atom pressure tensor. It is also really a stress-volume formulation. It would need to be divided by a per-atom volume to have units of stress, but an individual atom's volume is not easy to compute in a deformed solid or a liquid. 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 3V, where V is the volume of the system, the result should be -P, where P is the total pressure of the system.
These lines in an input script should yield that result (assuming there is no fix or long-range contribution to the stress):
compute global all pressure thermo_temp compute peratom all stress/atom compute p all sum peratom variable p equal div(add(add(c_p1,c_p2),c_p3),mult(3.0,vol)) thermo_style custom step temp etotal c_global v_p
IMPORTANT NOTE: The per-atom stress does NOT include contributions due to fixes (e.g. SHAKE) or long-range Coulombic interactions (via the kspace_style command). The former needs to be added to LAMMPS. We're not sure if the latter is possible to compute. There are also a few pair styles for manybody potentials that are not yet instrumented to yield per-atom stress. See the Restrictions below.
Restrictions:
These pair styles do not yet tabulate per-atom virial information to allow them to work with this compute: airebo, meam, and TIP4P.
Related commands:
compute pe, compute stress/atom
Default: none