The thermal conductivity of a molecular dynamics system is calculated using the Green-Kubo fluctuation Dissipation Theorem
For this particular case the Tersoff potential was used for Silicon and the thermal conductivity was calculated at 1000 K. The calculations using the double exponential fit are done using references: Chen, J., Zhang, G., & Li, B. (2010). How to improve the accuracy of equilibrium molecular dynamics for computation of thermal conductivity?. Physics Letters A, 374(23), 2392-2396.
Files:
-
Readme.md : just description
-
1000kheatflux.in Input file for lammps MD code
-
SiC.tersoff Tersoff potential used for lammps MD simulation
-
"Green Kubo Thermal Conductivity Calculator.ipynb" Jupyter Python Notebook for calculating Thermal conductivity using the heat current flux outputted by LAMMPS
DESCRIPTION:
This notebook is to calculate the thermal conductivity using the Green-Kubo Heat Flux autocorrelation method. The steps and techniques are followed from Chen et al: Chen, J., Zhang, G., & Li, B. (2010). How to improve the accuracy of equilibrium molecular dynamics for computation of thermal conductivity?. Physics Letters A, 374(23), 2392-2396.
We have used a 15x15x15 (lattice parameter) Si system with a total of 27,000 atoms at 1,000K.
LAMMPS is told to calculate and output the total heat flux using compute heat/flux in the folowing way:
compute myKE all ke/atom
compute myPE all pe/atom
compute myStress all stress/atom NULL virial
compute flux all heat/flux myKE myPE myStress
variable Jx equal c_flux[1]/vol
variable Jy equal c_flux[2]/vol
variable Jz equal c_flux[3]/vol
The heat flux outputs for Jx, Jy, Jz were converted to a .csv file where it is then imported into python via pandas
Note: Chen et al uses Stillinger-Weber potential for Si. In this notebook, tersoff potential is used.
To do Solve for the Autocorrelation function we use the algorithm formulated in the wikipedia article and stackoverflow: https://stackoverflow.com/questions/14297012/estimate-autocorrelation-using-python
https://en.wikipedia.org/wiki/Autocorrelation#Estimation
Basically the autocorrelation function for a certain correlation length
can be described as:
where is the mean of the data set and
is the variance of the data set X. The total time (or length) of the data set is n
We need to do this for all three directions of the heat flux
Notice that we use the np.correlate function which is a full convolution of the data as shown in the following links: https://numpy.org/doc/stable/reference/generated/numpy.correlate.html
https://numpy.org/doc/stable/reference/generated/numpy.convolve.html#numpy.convolve