spatial profile of stress printing


I have been trying to preheat a simulation box containing BCC crystal upto 2000K and then perform some studies on the spatial profile of stress. I have been using the following script and the obtained output files are attached. First NPT is ran to bring system to 0 bar and temperature 2000K. Then NVE is ran for some time to hold this temperature. And it accurately does so, as can be seen from log.lammps.
units metal

boundary p p p

velocity all create 300 12345 rot yes dist gaussian
fix equilibration all npt temp 300 2100 0.1 iso 0 0 1 drag 1
thermo 500
thermo_style custom step pxx pyy pzz lx ly lz temp etotal vol
run 100000
unfix equilibration
fix nve all nve
run 100000
unfix nve

And then to check the spatial profile of stress , box has been binned in Z-direction and centre of mass velocity of each bin has been subtracted to calculate Kinetic energy contribution to stress plus interaction energy contribution as shown, (This has been done since later on I want to check the stress due to impact in z direction). For this following lammps commands has been added to the script.

compute mytemp all temp/partial 1 1 0
compute peratom all stress/atom mytemp ke pair
compute vorvol all voronoi/atom

compute zavg all chunk/atom bin/1d z lower 10 units lattice
variable meanpress atom -(c_peratom[1]+c_peratom[2]+c_peratom[3])/3/(c_vorvol[1]+1e-99)
variable meanpresszz atom -(c_peratom[3])/(c_vorvol[1]+1e-99)
variable meanpressxx atom -(c_peratom[1])/(c_vorvol[1]+1e-99)
fix pres_prof all ave/chunk 100 5 500 zavg v_meanpress file pressure.file
fix presz_prof all ave/chunk 100 5 500 zavg v_meanpresszz file pressurezz.file
fix presx_profile all ave/chunk 100 5 500 zavg v_meanpressxx file pressurexzchunk.file
run 1500

My question is even if from log.lammps file we can see that the pzz value comes down to around -109 bar which I consider approximately equal to 0 bar. Why is the value of stress so high in z direction when it has been printed binwise. i.e. around 12000 bar. Check pressurezz.file for the chunkwise printing of stress in each bin.

What is that I can be doing wrong or may be understanding wrong! I feel the problem lies in stress/atom calculation. But I can’t figure it out.

Thanks and Regards,

pressurezz.file (4.33 KB)

log.lammps (55.4 KB)

pressurexzchunk.file (4.17 KB)

pressure.file (4.34 KB)

Here is what I suggest you do. Start with a simpler
analysis. Just bin zzz from compute stress/atom
and use the bin-volume to normalize and get zzz. Don’t do voronoi,
don’t subtract the mean, don’t time ave, etc. You
can tell compute stress/atom to just do virial
or include KE. Once that works and you verify
the answer is what you expect, then try adding
those other options.