Application of atomic stress to compute heat flux via molecular dynamics for systems with many-body interactions

Abstract

Although the computation of heat flux and thermal conductivity either via Fourier's law or the Green-Kubo relation has become a common task in molecular dynamics simulation, contributions of three-body and larger many-body interactions have always proved problematic to compute. In recent years, due to the success when applying to pressure tensor computation, atomic stress approximation has been widely used to calculate heat flux, where the LAMMPS molecular dynamics package is the most prominent propagator. We demonstrated that the atomic stress approximation, while adequate for obtaining pressure, produces erroneous results in the case of heat flux when applied to systems with many-body interactions, such as angle, torsion, or improper potentials. This also produces incorrect thermal conductivity values. To remedy this deficiency, by starting from a strict formulation of heat flux with many-body interactions, we reworked the atomic stress definition which resulted in only a simple modification. We modified the LAMMPS package accordingly to demonstrate that the new atomic stress approximation produces excellent results close to that of a rigid formulation.

Unification of nonequilibrium molecular dynamics and the mode-resolved phonon Boltzmann equation for thermal transport simulations

“…The simulations are then run for 20 ns, with the data taken within the last 10 ns being used to extract the heat transport properties, e.g., the temperature and heat flux. Although it has recently been shown that LAMMPS might incorrectly implement the heat flux formula in some cases [48,49], it is accurate for crystalline silicon [50]. The heat flux values in the sample are also close to those obtained by energy conservation [30] in our cases.…”

Thermal conductivity reduction in carbon nanotube by fullerene encapsulation: A molecular dynamics study

“…From the literature results [29,35], it is known that the optimized Tersoff potential gives higher thermal conductivity for SWCNTs than both generations of the Brenner potential. The second difference comes from the fact that some previous works [16,17,20] used the LAMMPS package [39], which has an incorrect implementation of the heat current for many-body potentials [35,[40][41][42]. This incorrect heat current also has an overall effect of underestimating the thermal conductivity of a low-dimensional material described by a many-body potential [35].…”

Thermal conductivity of aqueous solutions of reline, ethaline, and glyceline deep eutectic solvents; a molecular dynamics simulation study

“…It is important to note that both EMD methods require the heat flux vector J for computing thermal conductivities [73]. As shown in the recent studies by Surbyls et al [85] and Boone et al [86], the computation of J in molecular simulation packages such as LAMMPS can be erroneous for molecular systems having angles, torsions, dihedrals, and improper dihedrals, producing unphysical values. In fact, our preliminary thermal conductivity computations for reline using the Einstein method (OCTP plugin in LAMMPS [84]) yielded values which were almost an order of magnitude higher than the respective experiments.…”