1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
|
<HTML>
<CENTER><A HREF = "http://lammps.sandia.gov">LAMMPS WWW Site</A> - <A HREF = "Manual.html">LAMMPS Documentation</A> - <A HREF = "Section_commands.html#comm">LAMMPS Commands</A>
</CENTER>
<HR>
<H3>compute stress/atom command
</H3>
<P><B>Syntax:</B>
</P>
<PRE>compute ID group-ID stress/atom keyword ...
</PRE>
<UL><LI>ID, group-ID are documented in <A HREF = "compute.html">compute</A> command
<LI>stress/atom = style name of this compute command
<LI>zero or more keywords may be appended
<LI>keyword = <I>ke</I> or <I>pair</I> or <I>bond</I> or <I>angle</I> or <I>dihedral</I> or <I>improper</I> or <I>kspace</I> or <I>fix</I> or <I>virial</I>
</UL>
<P><B>Examples:</B>
</P>
<PRE>compute 1 mobile stress/atom
compute 1 all stress/atom pair bond
</PRE>
<P><B>Description:</B>
</P>
<P>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 and is stored as a 6-element vector in the following order:
xx, yy, zz, xy, xz, yz. See the <A HREF = "compute_pressure.html">compute
pressure</A> command if you want the stress tensor
(pressure) of the entire system.
</P>
<P>The stress tensor for atom <I>I</I> is given by the following formula,
where <I>a</I> and <I>b</I> take on values x,y,z to generate the 6 components of
the symmetric tensor:
</P>
<CENTER><IMG SRC = "Eqs/stress_tensor.jpg">
</CENTER>
<P>The first term is a kinetic energy contribution for atom <I>I</I>. The
second term is a pairwise energy contribution where <I>n</I> loops over the
<I>Np</I> neighbors of atom <I>I</I>, <I>r1</I> and <I>r2</I> are the positions of the 2
atoms in the pairwise interaction, and <I>F1</I> and <I>F2</I> are the forces on
the 2 atoms resulting from the pairwise interaction. The third term
is a bond contribution of similar form for the <I>Nb</I> bonds which atom
<I>I</I> is part of. There are similar terms for the <I>Na</I> angle, <I>Nd</I>
dihedral, and <I>Ni</I> improper interactions atom <I>I</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 <I>Nf</I>
<A HREF = "fix.html">fixes</A> that apply internal constraint forces to atom <I>I</I>.
Currently, only the <A HREF = "fix_shake.html">fix shake</A> and <A HREF = "fix_rigid.html">fix
rigid</A> commands contribute to this term.
</P>
<P>IMPORTANT NOTE: For granular systems, this formular neglects the
contribution of average velocity in the kinetic energy contribution.
This is corrected in the compute ave/euler command (currently no doc
available).
</P>
<P>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 a water molecule via the <A HREF = "fix_shake.html">fix shake</A>
command.
</P>
<P>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. The
<I>virial</I> keyword means include all terms except the kinetic energy
<I>ke</I>.
</P>
<P>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.
</P>
<P>The <A HREF = "dihedral_charmm.html">dihedral_style charmm</A> style calculates
pairwise interactions between 1-4 atoms. The virial contribution of
these terms is included in the pair virial, not the dihedral virial.
</P>
<P>The KSpace contribution is calculated using the method in
<A HREF = "#Heyes">(Heyes)</A> for the Ewald method and by the methodology described
in <A HREF = "#Sirk">(Sirk)</A> for PPPM. The choice of KSpace solver is specified
by the <A HREF = "kspace_style.html">kspace_style pppm</A> command. Note that for
PPPM, the calcluation 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.
</P>
<P>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, meaning the computed quantity is in units of
pressure*volume. It would need to be divided by a per-atom volume to
have units of stress (pressure), but an individual atom's volume is
not well defined or 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 dV, where
d = dimension and V is the volume of the system, the result should be
-P, where P is the total pressure of the system.
</P>
<P>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:
</P>
<PRE>compute peratom all stress/atom
compute p all reduce sum c_peratom[1] c_peratom[2] c_peratom[3]
variable press equal -(c_p[1]+c_p[2]+c_p[3])/(3*vol)
thermo_style custom step temp etotal press v_press
</PRE>
<P><B>Output info:</B>
</P>
<P>This compute 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. See <A HREF = "Section_howto.html#howto_15">Section_howto
15</A> for an overview of LAMMPS output
options.
</P>
<P>The per-atom array values will be in pressure*volume
<A HREF = "units.html">units</A> as discussed above.
</P>
<P><B>Restrictions:</B> none
</P>
<P><B>Related commands:</B>
</P>
<P><A HREF = "compute_pe.html">compute pe</A>, <A HREF = "compute_pressure.html">compute pressure</A>
</P>
<P><B>Default:</B> none
</P>
<HR>
<A NAME = "Heyes"></A>
<P><B>(Heyes)</B> Heyes, Phys Rev B 49, 755 (1994),
</P>
<A NAME = "Sirk"></A>
<P><B>(Sirk)</B> Sirk, Moore, Brown, J Chem Phys, 138, 064505 (2013).
</P>
</HTML>
|
