Inconsistent results for units lj vs real for Argon

Hi,

I am trying to develop an input script in lj units . As a validation of the fact that results would be independent of the unit system chosen, I set up a small script for Argon as follows:

Hi,

I am trying to develop an input script in lj units . As a validation of the fact that results would be independent of the unit system chosen, I set up a small script for Argon as follows:


REAL UNIT


log logdim.log
units real
dimension 3
boundary p p p
lattice fcc 5.376 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
region box block 0 1 0 1 0 1
create_box 1 box
create_atoms 1 box
mass 1 39.948
pair_style lj/cut 13.0
pair_coeff * * 0.2381 3.405
timestep 1
thermo 100
dump 1 all custom 100 dim.lammpstrj x y z
thermo_style custom step temp pe ke etotal
velocity all create 119.91 102486 mom yes rot yes dist gaussian
fix NVT all nvt temp 119.91 119.91 10 drag 0.2
run 1


Assuming reference values of 0.2381 Kcal/mol for energy (epsilon) , 3.405 A for distance (sigma), 39.948 gm/mol for mass, I obtained time reference= (epsilon/m/sigma^2)^(-0.5)=2.161e-12 s and temperature reference = epsilon/Kb= 119.91 K (considering Kb = .0019872 Kcal/mol) . I cross checked this value from online references. The only other input to be determined was the lattice “scale” which here needs to be the reduced lj density rho. For fcc, N=no. of basis atoms 4, V=volume of Argon unit cell= (5.376 A)^3, Thus, using (N/V)sigma^3, I find rho=1.016. With this Lammps computes the lattice spacing 1.57902 which correctly corresponds to 1.57902sigma=5.376 A in real units.
The problem is, the energy values (pe, ke, etotal etc. ) get from lammps in lj units does not convert correctly to real values. Since epsilon is .2381 Kcal/mol, a value of energy in lj units, say 1.23 should correspond to 1.23
.2381 Kcal/mol, right? However, that is not the case with my results.

Here is the script for lj units


log lognondim.log
units lj

setup problem

dimension 3
boundary p p p
lattice fcc 1.016 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1 #reduced density is used
region box block 0 1 0 1 0 1
create_box 1 box
create_atoms 1 box
mass 1 1 #mass scale =39.948 g/mol
pair_style lj/cut 2.5
pair_coeff * * 1 1 #energy scale=0.2381 Kcal/mol Distance scale = 3.405 Angstrom
timestep 0.00046082 #time scale= 2.16*1e-12 s , which means time step is 1 fs
thermo 100
dump 1 all custom 100 nondim.lammpstrj x y z
thermo_style custom step temp pe ke etotal
velocity all create 1 102486 mom yes rot yes dist gaussian
fix NVT all nvt temp 1 1 0.0046082 drag 0.2
run 1


The output from real simulation:

Step Temp PotEng KinEng TotEng
0 119.91 -7.9023794 1.0722868 -6.8300926
1 119.90047 -7.9022939 1.0722016 -6.8300924

The output from lj simulation

Step Temp PotEng KinEng TotEng
0 1 -8.0073844 1.125 -6.8823844
1 0.99992157 -8.007296 1.1249118 -6.8823842

I fail to understand how the three energies in lj units can be related to real units?
I need to understand this relation to validate that I have setup my lj units sytem correctly.

Any pointer to what am I missing still is greatly appreciated.

Thermo output with lj units is normalized by default, with real units it is not.

Axel

Thanks Dr. Kohlmeyer for your very prompt response.
Indeed the values of temperature, energy etc. are normalized with respect to the temperature and energy scales in lj unit. I do understand that part. Thus the value ‘1’ in lj output at timestep ‘0’ corresponds to to 119.91 in real units. This is fine.

However, for energy, lets say pe, the system pe at timestep 0 is -7.9023794 Kcal/mol. Thus, when normalized by energy scale 0.2381 Kcal/mol, this becomes 33.1893. However, the value I get from lj for pe at time step ‘0’ is -8.007296. Thus I dont understand, what normalization factor is lammps using for energies? It should normalize with respect to 0.2381 Kcal/mol, but it seems that lammps is using some other factor here.

Thanks and regards,
Souvik.

Thanks Dr. Kohlmeyer for your very prompt response.
Indeed the values of temperature, energy etc. are normalized with respect to the temperature and energy scales in lj unit. I do understand that part. Thus the value ‘1’ in lj output at timestep ‘0’ corresponds to to 119.91 in real units. This is fine.

No, you don’t understand. This is about the thermo_modify option “norm” being set to yes in reduced units and to no by default for all other units settings in lammps.

Thank you Dr. Kohlmeyer again!

Now everything makes sense.

Thanks and regards,
Souvik.