# Tutorial: Thermal conductivity from EMD

Jump to navigation
Jump to search
### The

### The

## Contents

## Introduction

In this example, we use the EMD (Green-Kubo) method to calculate the lattice thermal conductivity of graphene at 300 K and zero pressure.

## Preparing the Inputs

### The `xyz.in`

file

- The
`xyz.in`

file used is the same as in Tutorial: Density of states.

- Note that the thickness of the graphene sheet is set to 3.35 A according to the convention in the literature. This thickness is needed to calculate an effective 3D thermal conductivity for a 2D material.

### The `run.in`

file

- The
`run.in`

input file is given below:

potential potentials/tersoff/Graphene_Lindsay_2010_modified.txt velocity 300 # equilibration ensemble npt_ber 300 300 0.01 0 0 0 0.0005 time_step 1 dump_thermo 1000 run 1000000 # production ensemble nve compute_hac 20 50000 10 run 10000000

- The first line of command with the potential keyword tells that the potential to be used is specified in the file Graphene_Lindsay_2010_modified.txt.

- The second line of the command with the velocity keyword tells that the velocities will be initialized with a temperature of 300 K.

- There are two runs.
- The first run serves as the equilibration stage, where the NPT ensemble (the Berendsen barostat) is used. The temperature is 300 K and the pressures are zero in all the directions. The coupling constants are 0.01 (dimensionless) and 0.0005 (in the natural unit system adopted by GPUMD) for the thermostat and the barostat, respectively. The time step for integration is 1 fs. There are [math]10^6[/math] steps (1 ns) for this run and the thermodynamic quantities will be output every 1000 steps.
- The second run is for production. Here, the NVE ensemble is used. The line with the compute_hac keyword means that heat currents will be recorded every 20 steps (20 fs), 50000 HAC data (the maximum correlation time is then about 1 ns) will be calculated, and the HAC are averaged for every 10 data points before written out. The production time is 10 ns ([math]10^7[/math] steps), which is 10 times as long as the maximum correlation time. This is a reasonable choice.

## Results and discussion

### Results

- The above figure shows the results from three independent simulations, which took about two hours in total using a Tesla K40 card.
- As the system is essentially isotropic in the planar directions, we can average over the two in-plane directions.
- From (a), we can see that the in-plane component and the out-of-plane component of the HAC have different time scales. The latter decays much more slowly.
- Panel (b) shows the individual and averaged RTCs for the in-plane component [math]\kappa^{\rm in}(t)[/math]. The averaged RTC converges to about 1000 W/mK at around 100 ps.
- Panel (c) shows the individual and averaged RTCs for the out-of-plane component [math]\kappa^{\rm out}(t)[/math], and the convergence property is not very clear here. This is because the out-of-plane component converges very slowly and three independent simulations (each with 10 ns) are not enough to give accurate results.
- Summing up [math]\kappa^{\rm in}(t)[/math] and [math]\kappa^{\rm out}(t)[/math], we get [math]\kappa^{\rm tot}(t)[/math], as shown in panel (d).

### Discussion

- Accurately calculating thermal conductivity of graphene using the EMD method can be a very time consuming task. The results we presented here are from three independent simulations with a total production time of 30 ns. It can been seen that the HAC data already become very noisy when the correlation time is 100 ps. To obtain accurate results, one needs to do many independent simulations. Much more accurate data were presented in Fig. 2 of this paper. Here are the simulation parameters used in the reference paper which differ from those used in this example:
- The simulation cell size used in the reference paper is larger, which is about 25 nm x 25 nm (24000 atoms), instead of 15 nm x 15 nm (8640 atoms) here.
- The maximum correlation time used in the reference paper is larger, which is 10 ns, instead of 1 ns here.
- The production time used in the reference paper for one independent simulation is larger, which is 50 ns, instead of 10 ns here.
- There were 100 independent simulations in the reference paper, instead of three here.
- Each independent simulation in the reference paper took about 10 GPU hours (using Tesla K40) and about 1000 GPU hours were used to obtain the results shown in Fig. 2 of the reference paper.