Energy not conserved in NVE

I have run a LAMMPS run for SiC bulk system using the Tersoff/ZBL as supplied by LAMMPS distribution.

I have done the NPT first for 50 ps with a timestep of 1 fs and then switched it off and started NVE for 20 ps.

I noticed that during NPT T is maintained at 300 K. But as soon as I switched off NPT and switched on NVE, T increased from 300 K to 310 K and total energy also increased by ~75 eV.

I am pasting the input file below. I can’t find any obvious mistake in my input file.
Do you have any suggestion ?

echo both
units metal
dimension 3
#Periodic boundary conditions
boundary p p p


variable temp equal 300.0
variable seed equal 3402219

#set up
atom_style atomic
#kspace_style ewald 1.0e-5
neighbor 1.0 bin
neigh_modify delay 0 every 1 check yes

#Create a lattice of zincblende SiC ( 3C - SiC )

mass 1 28.1 # type 1 = Si
mass 2 12.0 # type 2 = C

#Choose the potentials to be used
pair_style tersoff/zbl
pair_coeff * * SiC.tersoff.zbl Si C

velocity all create {temp} {seed} rot yes mom yes dist gaussian

fix 1 all npt temp {temp} {temp} 0.1 iso 0.0 0.0 1
thermo 500
thermo_style custom step temp pe ke etotal press
run 50000


unfix 1


fix 2 all nve
compute strf all stress/atom
dump 1 all custom 5 all.dump id type x y z vx vy vz c_strf[1] c_strf[2] c_strf[3] c_strf[4] c_strf[5] c_strf[6]
dump_modify 1 sort id
run 20000
undump 1

#write a restart file with atoms still in excited positions
write_restart final.RS


Have you tried the usual things - longer equilibration,
smaller timestep?


Thanks for your comments.
Reducing stepsize to 0.1 fs helps to get energy conservation in NVE for SiC.

Note that if total energy (KE + PE) is conserved

when running NVE (after equil), then the
fact that temperature rises probably just
means you are not well-equlibrated yet.

If total E is not conserved, then reducing the
timestep is a good thing to try.


Thanks again.
I have checked the thermodynamic output and both T and P are seen to be stable at target values.