Hi all,

lately I am dealing with thermal conductivity of silicon, precisely it is 5x5x300nm system of crystal silicon.

Here is a part of my code to LAMMPS:

#npt equillibration

velocity all create 300 345941 dist gaussian

fix 1 all npt temp 300 300 .1 iso 1 1 8

run 2000

My question is: why during my calculations temperature decreasing roughly to 150K, where I can observe some oscillations of temperature, and after that (after 1000 steps) reaches 300K?

Maybe that is why my final value of thermal conductivity is twice bigger than it should be?

I am enclosing graph of my temperature of my sample in a function of time.

Thank you in advance for your answer.

lately I am dealing with thermal conductivity of silicon, precisely it is 5x5x300nm system of crystal silicon.

Here is a part of my code to LAMMPS:

#npt equillibration

velocity all create 300 345941 dist gaussian

fix 1 all npt temp 300 300 .1 iso 1 1 8

run 2000

My question is: why during my calculations temperature decreasing roughly to 150K, where I can observe some oscillations of temperature, and after that (after 1000 steps) reaches 300K?

Is your npt "time" (Tdamp in http://lammps.sandia.gov/doc/fix_nh.html) big

enough? A too small Tdamp can cause oscillations.

Maybe that is why my final value of thermal conductivity is twice bigger than it should be?

No idea.

you are over trusting on the initial velocity distribution. you still need to equilibrate the system, and although the velocity distribution is indeed gaussian at equilibrium, the inverse is (trivially) not true. in this case, your initial layout has less potential energy than the equilibrium and it exchanges kinetic energy (ie, decreases temperature) to gain the potential energy. this is (apparently) much faster than your npt thermostat, and only once this happens you can see the thermostat “pulling” the system back to 300 K.

if you are taking this equilibration points into account, you will get wrong results.