Lammps openkim edip potential for carbon

OpenKIM
1

Hi,
I installed the EDIP potential for carbon described in Marks paper with OpenKIM. My system is an amorphous carbon with 300000 atoms.
Everything runs. EXCEPT I noticed something really weird. When I plot the stress per atom outputted by LAMMPS, I can only positive ones. This is not normal. In addition, if I output the pressure with LAMMPS and take into consideration sign and volume, it doesn’t lead me to the same values Tham with the stress per atom. Has anyone tried to output the stress per atom given by LAMMPS when given this potential? When I follow the same procedure with a different potential, I get a distribution that is much more expected and has negative and positive values.
here is the beginning of my code:

metal units, pressure in bar

kim init EDIP_LAMMPS_Marks_2000_C__MO_374144505645_000 metal
atom_style atomic
read_restart ${file}

neighbor 1.0 bin
neigh_modify every 1 delay 5 check yes

kim interactions C

Define output parameters

variable nrelax equal 10000
variable nthermo2 equal 10000
variable ndump equal 50000 # every 0.5% def
variable nrestart equal 500000 # every 5% def

Define MD parameters

variable temp equal ${temperature} # temperature of initial sample
variable timestep equal 0.001 # timestep, in ps because of metal unit
variable erate equal 0.0001 # strain rate

Setup simulation control

compute 1 all stress/atom NULL
compute 2 all pe
variable g equal c_2
variable k equal pxx
variable l equal lx

thermo_style custom step time temp pe pxx pyy pzz pxy pxz pyz lx ly lz
thermo_modify norm no # default parameter for metal units anyway. means we don’t normalize by the number of atoms
timestep ${timestep}

reset velocities to temp with random seed

velocity all create ${temp} 5 rot yes dist gaussian

relax for a short time before deformation

fix 1 all npt temp {temp} {temp} (100.0*dt) aniso 0 0 (1000.0*dt) drag 1.0 # is drag necessary?

run ${nrelax}

2

Hi Raphaelle,

I am not an expert in atomistic stress so I may have to ask around or direct you somewhere else to ask this question. It’s highly unlikely that someone who just happened to run a similar situation with this specific potential will come across your post.

First I will try to help, and I have some follow up questions:

  1. What other potentials have you tried that DO produce the result you expected? Were they KIM potentials or built-in LAMMPS pair_styles?
  2. When you say "In addition, if I output the pressure with LAMMPS and take into consideration sign and volume, it doesn’t lead me to the same values Tham [sic] with the stress per atom. ", do you mean something like that the following lines (from compute stress/atom command — LAMMPS documentation) don’t produce identical columns as they should?
compute        peratom all stress/atom NULL
compute        p all reduce sum c_peratom[1] c_peratom[2] c_peratom[3]
variable       press equal -(c_p[1]+c_p[2]+c_p[3])/(3*vol)
thermo_style   custom step temp etotal press v_press

Cheers, ilia

3

Hi Ilia,

  1. I tried ML-PACE potential which is built-in and not KIM potentials. Although I also want to say that I used to have an older version of the carbon EDIP potential that wasn’t OpenKIM but gave the same wrong results.
  2. something like: “variable k equal pxx” so I know pxx, compared to “compute 1 all stress/atom NULL” where I can extract sigmaxx
4

I used to have an older version of the carbon EDIP potential that wasn’t OpenKIM but gave the same wrong results.

In that case I think it’s pretty clear that this is either a peculiarity with this parametrization of EDIP, or your expectations of the results are incorrect (I have no specific reason to believe this since I’m not an expert in atomistic stress calculations, but always a possibility!)

In either case, I think you’d be able to get more help in the general LAMMPS forum. You may also want to consider contacting Nigel Marks directly to see if this is a known behavior in this potential.

1 Like