Hello,
We are currently trying to learn how to use LAMMPS and we are having a begginer problem with the viscosity case in.wall.2d.
When we try to calculate the pressure manually (by reducing the result of compute atom/stress), we do not obtain the same results as the compute pressure command when the imposed shear rate on the top plate is not 0. However, when this shear rate is 0, the computed pressure is identical to the one we manually obtain. Furthermore, Pyy is always equal to the one we calculate manually.
What is the reason for this difference and how can we manully recompute the pressure tensor correctly when the shear rate is not 0?
Best regards,
Bruno
Joined LAMMPS input script :
sample LAMMPS input script for viscosity of 2d LJ liquid
use shearing wall, thermostat via fix langevin
settings
variable x equal 20
variable y equal 20
variable ylo equal -2.5
variable yhi equal 23
variable rho equal 0.6
variable t equal 1.0
variable rc equal 2.5
variable srate equal 2.7
problem setup
units lj
dimension 2
atom_style atomic
neigh_modify delay 0 every 1
lattice sq ${rho}
region simbox block 0 x {ylo} ${yhi} -0.1 0.1
create_box 3 simbox
create_atoms 1 box
mass * 1.0
pair_style lj/cut ${rc}
pair_coeff * * 1 1
region lower block INF INF INF 0.0 INF INF
region upper block INF INF $y INF INF INF
group lower region lower
group upper region upper
set group lower type 2
set group upper type 3
group wall union lower upper
group flow subtract all wall
mass 2 2.0
mass 3 2.0
velocity flow create t 97287 velocity upper set {srate} 0.0 0.0 units box
compute thermal flow temp/partial 0 1 0
compute flow flow temp
fix 1 all nve
fix 2 flow langevin $t $t 0.1 498094
fix_modify 2 temp thermal
fix 3 wall setforce 0.0 0.0 0.0
fix 4 flow ave/spatial 20 500 10000 y center 0.05 vx &
units reduced file profile.wall.2d
fix 5 all enforce2d
equilibration run
variable ybox equal y*ylat compute tilt flow temp/ramp vx 0 {srate} y 0 ${ybox} units box
thermo 1000
thermo_style custom step temp c_tilt epair etotal press pxy
dump 1 all image 1000 image.*.jpg mass type zoom 1.6 adiam 1.2
run 3000
data gathering run
compute atomstress all stress/atom
compute newpxx all reduce sum c_atomstress[1]
compute newpyy all reduce sum c_atomstress[2]
compute newpxy all reduce sum c_atomstress[4]
variable visc equal -pxy/(v_srate/ly)
variable newpxy equal -c_newpxy/(vol)
variable visc2 equal -v_newpxy/(v_srate/ly)
variable newpyy equal -c_newpyy/vol
variable newpxx equal -c_newpxx/vol
variable normal equal v_newpyy-v_newpxx
fix vave all ave/time 1000 1 1000 v_visc v_visc2 v_normal file visc.out
thermo_style custom step temp pxy v_newpxy pyy v_newpyy pxx v_newpxx v_visc v_visc2 v_normal
thermo_modify temp tilt
thermo 1
only need to run for 5000 steps to make a good 100-frame movie
#dump 1 all custom 50 dump.wall.2d id type x y z vx
dump 2 all image 1000 image.*.jpg mass type zoom 1.6 adiam 1.2
#dump_modify 2 pad 5 amap 0.0 ${srate} ca 0.0 2 min blue max red
run 50