N2 MSD Calculation

Dear all,

I want to verify my MD program output using lammps.
Using my program I calculated equilibrium state of N2 (as a point)
The experimental self-diffusion coefficient is 14.2 (in our units, nm^2 / ps), so in the simulation it can be slightly bigger or less.
My MD program does such calculation, and gets 13.8 by MSD and 17.3 by Green-Kubo (VACF).
The self-diffusion coefficient, calculated using MSD and VACF in LAMMPS differs a lot from the the experimental one.
Using LAMMPS with the script below I get next values for MSD:

Step twopoint fitslope
...
  6999000 38250.766 36939.342
  7000000 38250.995 36939.716

Step diffvacf
...
6999000 39858.095
7000000 39858.254

They are not correct, (in lammps nano units the self-diffusion coefficient should be close to 14200, and we have 36939 by MSD and 39858 by VACF)

units nano
dimension 3
boundary p p p
atom_style atomic

read_data test.data # equilibrium state of N2 (T = 273.14K)

pair_style mie/cut 0.91225 #such potential was used in my program
# eps sig ga gr
pair_coeff 1 1 0.0135 0.3649 6 11.5

mass 1 0.000046517 # mass of N2 molecula

timestep 0.000002 # the same timestep I used in my program
run_style verlet

variable dt_3 equal 0.000002

fix 1 all nve

compute msd all msd com yes
variable twopointmsd equal c_msd[4]/(3*2)/(step*v_dt_3+2.0e-6)
# the divison is from the einshtein formula for msd method, i looked through the LAMMPS src and find out that this operation do no performed by lammps.

fix 9 all vector 10 c_msd[4]
variable fitslopemsd equal slope(f_9)/6/(10*v_dt_3)

thermo_style custom step v_twopointmsd v_fitslopemsd

#for VACF it is
#compute VACF all vacf
#variable VACF4 equal c_VACF[4]
#fix 5 all vector 1 c_VACF[4]
#variable diffvacf equal v_dt_3*trap(f_5)
#variable diffvacfd3 equal v_diffvacf/3
#thermo_style custom step v_diffvacfd3

thermo $d
run 7000000

Can anybody tell me, where is the problem ?

Are the dynamics of the 2 codes the same. E.g. do

they agree in energy and forces. Nano units in LAMMPS

is not used much (I’'m guessing). I would verify from

a dump file (first snapshot, final snapshot) that atoms

have displaced consistent (on average) with what you

are getting as MSD from LAMMPS. If the MSD is

the correct numerical value, then the difference is

in your dynamics between the 2 codes, not the

way MSD/VACF is calcualted.

Steve