Dear LAMMPS users and developers,
I am trying to understand how the pressure tensor from “compute pressure” (https://lammps.sandia.gov/doc/compute_pressure.html) is calculated under PBC. According to this page, “the position and force vectors of ghost atoms are thus included in the summation”. It’s not clear whether the position and force vectors of original atoms are also included (force vectors of original and ghost atoms should be the same?).
In theory, the trajectory can output position (x, y, z), unwrapped position (xu,yu,zu), image flags(ix,iy,iz), force (fx,fy,fz). So whether we are only including ghost atoms (we get that information from image flags) or both original and ghost atoms, we should be able to calculate the pressure tensors that should be the same to pxx, pyy, pzz, pxy, pxz, pyz from the thermo_style, right?
I am using the following toy problem for this verification.
units real
variable T equal 86.4956
variable V equal vol
variable dt equal 4.0
variable p equal 400 # correlation length
variable s equal 5 # sample interval
variable d equal $p*$s # dump interval
convert from LAMMPS real units to SI
variable kB equal 1.3806504e-23 # [J/K] Boltzmann
variable atm2Pa equal 101325.0
variable A2m equal 1.0e-10
variable fs2s equal 1.0e-15
variable convert equal {atm2Pa}*{atm2Pa}{fs2s}*{A2m}{A2m}*{A2m}
setup problem
dimension 3
boundary p p p
lattice fcc 5.376 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
region box block 0 4 0 4 0 4
create_box 1 box
create_atoms 1 box
mass 1 39.948
pair_style lj/cut 13.0
pair_coeff * * 0.2381 3.405
timestep ${dt}
thermo $d
equilibration and thermalization
velocity all create $T 102486 mom yes rot yes dist gaussian
fix NVT all nvt temp $T $T 10 drag 0.2
run 8000
reset_timestep 0
compute 1 all temp
thermo_style custom step temp pxx pyy pzz pxy pxz pyz
dump dumpTrj all custom 10 traj.trj id x y z ix iy iz xu yu zu vx vy vz fx fy fz
dump_modify dumpTrj sort id
run 0