[lammps-users] tip4p water and Green-Kubo methods?


I am trying to calculate thermal conductivity and viscosity of tip4p water using Green-Kubo formulae. I started with the example given on the heat/flux page (http://lammps.sandia.gov/doc/compute_heat_flux.html) and modified it for water. I was able to reproduce the results for example given on the heat/flux page for the LJ fluid. Unfortunately the numbers I get for water are an order of magnitude low for viscosity, and seven orders of magnitude to high for thermal conductivity. The integrals appear nicely converged. Is there something fundamental I am missing? (Yes, I triple checked the unit conversions) The one thing that comes to mind is that maybe the shake algorithm affects the results?

BTW: There are some errors in the input script given on the heat/flux webpage:

The line:

fix 		NPT all npt 70 70 10 xyz 0.0 0.0 100.0 drag 0.2

should read:

fix 		NPT all npt temp 70 70 10 x 0.0 0.0 100.0 y 0.0 0.0 100.0 z 0.0 0.0 100.0 drag 0.2


Dave Schall

Hello Dave,

The possible way to check your results is to calculate thermal conductivity by Muller-Plathe method. Usually these two methods (GK and MP) gives results of the same order of magnitude at least. I tried for SiC, Si with MEAM and SW potentials the difference between results was ~ 50% for temperature 500 K and below and few persents for the temperatures
above 1000 K.

All the best,

German Samolyuk

Changed the doc page.


This is the equation for calculating the thermal conductivity given on the command page:

kappa = V/(k_B*(T^2)) x "thermo" output frequency x timestep x Integral value

In what unit is the “thermo” output frequency? time or steps?

For example if my timestep size is 0.001 ps and my thermo command is set to 10, should the thermo output frequency in the above be equation be 10 steps or 0.01 ps?


German, I’m curious about the value of thermal conductivity of Si you get using GK method. Please also indicate your sample size. Do you observe any size effects? Do you think your calculated value is corresponding to experimental case, that is bulk diffusive value?

Best regards,


I'm sorry it takes me so long to answer your question. Somehow I lost
your e-mail.
Attached please find plot for thermal conductivity as a function of
upper limit of integration
for three sets of modeling cells. You can see size effect from the plot.

And, yes, I do think it corresponds to diffusive thermal conductivity
in bulk material.




Thanks for your reply. From this plot, I dont think you’ve reached the diffusive regime, in which the thermal conductivity is independent of the sample size. Why do you think that corresponds to the bulk case?

My experience of using GK method to calculate thermal conductivity is that the auto-correlation is sensitively dependent upon the correlation time you choose. If you plot the figure longer than 180 ps, the curves may fall down or go up. Well, then people argue that you can do many of the simulations and do ensemble average. Some do converge and some just dont. So I just dont get why so many people in this field use this kind of ill-programmed scheme. Probably because it does not require large sample size. But this is obviously not true, as you demonstrated in this plot.