Dear LAMMPS Users

I hope you are doing well.

I’m trying to model the interactions between two 2D layers in the vertical direction using a pair potential. I obtained some data from dft calculations and I performed a fit to obtain the parameters for LJ and for morse. The issue happens with both potentials but I’m going to illustrate the LJ case. The epsilon and sigma parameters simply are the min. energy and the corresponding distance between atoms. For example, I set the sigma to be 3 and the epsilon to be 1.53. However, I find out when minimizing the structures, the distance becomes approximately half of sigma and the minimum energy does not reach 1.53 kcal/mol per atom it was about 29 in the sample code below. Here is a sample script to calculate the binding energy between both layers and reproduce the issue I’m talking about. I’m using LAMMPS 27 May 2021.

**************************** BEGIN ****************************#

units real

dimension 3

boundary p p p

atom_style charge

#**************************** GEOMETRY DEFINITION ****************************#

box tilt large

read_data test.dat

#**************************** POTENTIAL ****************************#

pair_style hybrid tersoff tersoff lj/cut 6.0 0

pair_coeff * * tersoff 1 c_optimized.tersoff C NULL

pair_coeff * * tersoff 2 c_optimized.tersoff NULL C

pair_coeff 1 2 lj/cut 1.536325589 3.1262

#**************************** MINIMIZE ****************************#

variable etol equal 1.0E-8

variable ftol equal 0.0

variable iter equal 1000000

variable lat2 equal “ly/(sqrt(3)/2)”

group graphene type 1

group graphene2 type 2

thermo 100

thermo_style custom step temp vol press pe etotal evdwl ecoul fmax lx v_lat2

thermo_modify norm no

min_style cg

minimize {etol} {ftol} {iter} {iter}

delete_atoms group graphene

run 0

Kind Regards,

Khaled

test.dat (31.1 KB)

c_optimized.tersoff (492 Bytes)