I am doing a test run of shock wave propagation with Nickel using the mirror momentum method introduced in Oscar’s beginer’s guide of shockwaves. The shock velocity I got matches existing papers. But I am still not quite sure if I got the correct pressure.

I was calculating the stress in the direction of shock propagationi, i.e. Pzz. I used Voronoi/atom, ave/chunk for average stress.

Here is the part of command I used:

velocity piston set 0 0 v_Up sum no units box
fix 2 piston setforce 0.0 0.0 0.0
compute atom_vol bulk voronoi/atom
compute dvd all chunk/atom bin/1d z lower 3.52 nchunk once units box
variable meanpress atom -(c_peratom[1]+c_peratom[2]+c_peratom[3])/3/c_atom_vol[1]

My problem is: the stress obtained with voronoi at a certain chunk location is different from the stress obtained from this equation:

Pzz = rho_0U_pU_s.

I chose a location just behind the shock front for comparison. With voronoi, the pressure output is 210GPa (U_p = 2.7 km/s). With the equation above, the result is about 260GPa. Is this normal?

Which one is better if I want to know the local stress?

Another quick question: when we talk about “shock pressure”, is it the Pzz at the beginning of shock generation? I am asking this because I read in a paper saying the shock pressure is about 8GPa, but i got 25GPa of Pzz with same material and same shock condition.

velocity piston set 0 0 v_Up sum no units box
fix 2 piston setforce 0.0 0.0 0.0
compute atom_vol bulk voronoi/atom
compute dvd all chunk/atom bin/1d z lower 3.52 nchunk once units box
variable meanpress atom -(c_peratom[1]+c_peratom[2]+c_peratom[3])/3/c_atom_vol[1]

The last variable doesn;t make sense to me. There is no

compute with ID = peratom. I also don’t see what

chunks have to do with per-atom stress in your input

I am sorry, I missed my fix commands for output. The full input for property calculation is here:

velocity piston set 0 0 v_Up sum no units box
fix 2 piston setforce 0.0 0.0 0.0

compute peratom bulk stress/atom NULL

compute atom_vol bulk voronoi/atom
compute dvd all chunk/atom bin/1d z lower 3.52 nchunk once units box
variable meanpress atom -(c_peratom[1]+c_peratom[2]+c_peratom[3])/3/c_atom_vol[1]
fix pressure_profile bulk ave/chunk 10 5 100 dvd v_meanpress file pressure.profile

That’s a reasonable way to get local pressure. Someone more familiar

with shock simualtions will have to comment on this Q:

My problem is: the stress obtained with voronoi at a certain chunk location is different from the stress obtained from this equation:

Pzz = rho_0U_pU_s.

I chose a location just behind the shock front for comparison. With voronoi, the pressure output is 210GPa (U_p = 2.7 km/s). With the equation above, the result is about 260GPa. Is this normal?

Which one is better if I want to know the local stress?

Yes, I have seen a lot of researchers calculate pressure using this equation. But the method using LAMMPS commands sounds also a reasonable one. That’s why I am confused with the big difference. Maybe the voronoi estimate of volume has too much error?
Hope any shockwave colleagues could help.

Let me repeat my question here to make it easier to read: What is the difference between the two methods to calculate local stress in a shockwave plane:

Use equation: Pzz = rho_0U_pU_s

Use LAMMPS commands “compute stress/atom” and “compute voronoi/atom”, and compute stress by (S11+S22+S33)/3/vol.

Another quick question: when we talk about “shock pressure”, is it the Pzz at the beginning of shock generation? I am asking this because I read in a paper saying the shock pressure is about 8GPa, but i got 25GPa of Pzz with same material and same shock condition.