[lammps-users] compute stress/atom

Dear LAMMPS users,

    I am using the serial version of LAMMPS (Jan 23, 2009) on cygwin on Windows
XP with additional packages: opt, molecule, kspace, class2, manybody.

    I am simulating 125 Lennard Jones particles (epsilon = sigma = 1) in LJ
units at the number density of 0.85 at reduced temperature of 1.2 in NVE
ensemble. The volume of my simulation box is 147.059.

    I'd appreciate if you can help me resolve the the following unexpected
observations:

1. I used compute stress/atom command to calculate the pressure as described in
the documentation. However, the pressure calculated by LAMMPS is not same as the
pressure calculated using the 'compute' command, see column #3 and #4 in the
output. As mentioned in the documentation, these two values should match. Here
is the link to this documentation page:
http://lammps.sandia.gov/doc/compute_stress_atom.html

2. The 'sum' mode of 'compute reduce' outputs the sum divided by the number of
particles for the global stress components, see the column #4, #5 and #6 in the
output. This is not the expected behavior as explained the manual. The manual
states that LAMMPS should output the extensive value for 'sum' mode, i.e. #4 =
-(#5+#6+#7)/(3*vol), but this is not true in the output. Here is the link to
this documentation page:
http://lammps.sandia.gov/doc/compute_reduce.html

   Here is the thermo output:

See the script below.

(1) if you comment out the pair_modify tail yes command then pressure
will match stress. You are adding in a tail correction to the global pressure,
so it doesn't match.

(2) Comment in/out the thermo_modify norm no command. You will
see that the output of compute reduce sum is extensive, but it is normalized
when output by the thermo command. Intensive quantities like temp or press are
not normalized.

Steve

# 3d Lennard-Jones melt

units lj
atom_style atomic

lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 1 box
create_atoms 1 box
mass 1 1.0

velocity all create 1.44 87287 loop geom

pair_style lj/cut 2.5
pair_coeff 1 1 1.0 1.0 2.5
#pair_modify tail yes

neighbor 0.3 bin
neigh_modify delay 0 every 20 check no

fix 1 all nve

compute stpa all stress/atom
compute stgb all reduce sum c_stpa[1] c_stpa[2] c_stpa[3]
variable pr equal -(c_stgb[1]+c_stgb[2]+c_stgb[3])/(3*vol)

thermo_style custom step temp press v_pr c_stgb[1]
#thermo_modify norm no

thermo 10
run 100