I just wanted to mention the procedure I am following for the thermal conductivity computation in silicon - in case I am making any obvious errors.
For my thermal conductivity computations in silicon, I have attached the input files. In the first step, I do an NVT simulation at 300 K (NVTstep.input), then followed by NVE simulation (NVEstep.input) where I use the fix/thermal conductivity command. I monitor the energy exchange via f_2 and the value gets printed out in the log file at the end of the NVE simulation. My thermal conductivity computation is :
k = (Energy exchange in eV over duration of NVE simulation)/ (2AreaNVEsimulation time in picosecond* temperature gradient in Kelvin/Angstrom)
The cross-sectional area is 27.15 X 27.15 Angstrom squared. I follow appropriate conversion factors to get W/mK.
One concern I have is about the temperature gradient, which does not appear to be as linear as when I was doing Lennard-Jones systems following the same procedure. In fact, for the Lennard-Jones system mentioned in the Muller-Plathe paper, I was able to get a thermal conductivity value close to that mentioned in the paper as well as experimental data.
NVEstep.input (505 Bytes)
NVTstep.input (607 Bytes)