My Community

Dear all,

We have run the simulation test on silicon provided with VASP 6.3. Now, we want to calculate the thermal conductivity (Green-Kubo formula), but we obtain irrelevant values. It could be a problem of units, so the questions are: what is the unit of the heat flow q(t) calculated by VASP and stored in ML_HEAT file? Should we perform unit conversions to get the thermal conductivity in W/(m.K)?

Thank you
Best regards,
Pascal

The units of the ML_HEAT file are in eVA/fs.
Also you must run the heat flux calculations in the NVE ensemble. Relax your system sufficiently long to a given temperature and after that run the heat flux calculation within NVE.
Silicon is suprisingly a very difficult material to obtain the correct thermal conductivity.
You will possibly need huge supercells and very long MD runs. Also you need multiple runs with different trajectories over which you need to average.
Please read our paper for the heat flux of Zirconia for the methodology (mind that this material was much easier to calculate than silicon):
https://arxiv.org/abs/2109.15002

Hello,

Thank you for your reply. I indeed had the sense that the units are eV/(Angs.fs). Thanks for confirming my intuition.

I wrote a piece of Fortran code to calculate the thermal conductivity, but the values are inconsistent. Even though silicon is a difficult case, which I did not know, I am a bit surprised to get 0.197743160497E+12 eV/(Angs.K.fs), which is 0.526091109303E-06 W/(m.K)... So something must be wrong in my code. BTW, I ran the simulation at 300K, so I know that the experimental value is around 150 W/(m/K).

Is there a code around that calculates Kappa with which I can compare?

Eventually, I can send my piece of code for checking, but I understand that people are too busy to spend time on this.

Thank you
Best
Pascal

Yes if you are 6 orders of magnitude wrong there is definitely something wrong in your calculation.
Unfortunately I don't have a script for the Green-Kubo to share with you.

Anyway you also need a quite good sampling, because too bad calculational parameters will also lead to total garbage results (but probably not 6 orders of magnitude).

I asked Carla Verdi who did an extensive study for the heat flux of silicon using VASP and here is what she wrote about her best possible calculations:

- As a preparation, the system is first equilibrated in the NPT ensemble; then you run an NVT trajectory at the equilibrated volume. You can sample random initial configurations from this trajectory. For each of them, you run long NVE simulations and average the resulting thermal conductivities (or heat-flux autocorrelations).

- Time step 1.5 fs, supercell of 1728 atoms (I went up to 4096), 15 independent trajectories of 2.4 ns each; kappa from the integral of the heat-flux autocorrelation truncated around 300ps is 93 W/(mK).

Note here that we are almost a factor 2 smaller than in experiment.
So expect to be significantly below experiment.

Dear Ferenc,

Sorry for my, very, very,..., late reply.

Thank you for your hellp. We managed to run the simulations properly. We got some results that look reasonable.

Best regards,
Pascal